From: bangerth Date: Tue, 23 Sep 2008 16:52:39 +0000 (+0000) Subject: Add a comparison of solvers/preconditioners. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e2428c3e09bd1c04f1cb5de12da12ef153040d76;p=dealii-svn.git Add a comparison of solvers/preconditioners. git-svn-id: https://svn.dealii.org/trunk@16902 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-6/doc/results.dox b/deal.II/examples/step-6/doc/results.dox index a639c71e5b..8f442eb90c 100644 --- a/deal.II/examples/step-6/doc/results.dox +++ b/deal.II/examples/step-6/doc/results.dox @@ -161,3 +161,111 @@ through a few more functions from the destructor of the LaplaceProblem class, exactly where we have commented out the call to DoFHandler::clear(). + + + +

Possibilities for extensions

+ +One thing that is always worth playing around with if one solves +problems of appreciable size (much bigger than the one we have here) +is to try different solvers or preconditioners. In the current case, +the linear system is symmetric and positive definite, which makes the +CG algorithm pretty much the canonical choice for solving. However, +the SSOR preconditioner we use in the solve() function is +up for grabs. + +In deal.II, it is relatively simple to change the preconditioner. For +example, by changing the existing lines of code +@code + PreconditionSSOR<> preconditioner; + preconditioner.initialize(system_matrix, 1.2); +@endcode +into +@code + PreconditionSSOR<> preconditioner; + preconditioner.initialize(system_matrix, 1.0); +@endcode +we can try out different relaxation parameters for SSOR. By using +(you have to also add the header file lac/sparse_ilu.h to +the include list at the top of the file) +@code + PreconditionJacobi<> preconditioner; + preconditioner.initialize(system_matrix); +@endcode +we can use Jacobi as a preconditioner. And by using +@code + SparseILU preconditioner; + preconditioner.initialize(system_matrix); +@endcode +we can use a very simply incomplete LU decomposition without any +thresholding or strengthening of the diagonal. + +Using these various different preconditioners, we can compare the +number of CG iterations needed (available through the +solver_control.last_step() call, see @ref step_4 +"step-4") as well as CPU time needed (using the Timer class, +discussed, for example, in @ref step_12 "step-12") and get the +following results (left: iterations; right: CPU time): + + + + + + + +
+ @image html step-6.q2.dofs_vs_iterations.png + + @image html step-6.q2.dofs_vs_time.png +
+ +As we can see, all preconditioners behave pretty much the same on this +simple problem, with the number of iterations growing like ${\cal +O}(N^{1/2})$ and because each iteration requires around ${\cal +O}(N)$ operations the total CPU time grows like ${\cal +O}(N^{3/2})$ (for the few smallest meshes, the CPU time is so small +that it doesn't record). Note that even though it is the simplest +method, Jacobi is the fastest for this problem. + +The situation changes slightly when the finite element is not a +bi-quadratic one as set in the constructor of this program, but a +bi-linear one. If one makes this change, the results are as follows: + + + + + + + +
+ @image html step-6.q1.dofs_vs_iterations.png + + @image html step-6.q1.dofs_vs_time.png +
+ +In other words, while the increase in iterations and CPU time is as +before, Jacobi is now the method that requires the most iterations; it +is still the fastest one, however, owing to the simplicity of the +operations it has to perform. + +The message to take away from this is not that simplicity in +preconditioners is always best. While this may be true for the current +problem, it definitely is not once we move to more complicated +problems (elasticity or Stokes, for examples @ref step_8 "step-8" or +@ref step_22 "step-22"). Secondly, all of these preconditioners still +lead to an increase in the number of iterations as the number $N$ of +degrees of freedom grows, for example ${\cal O}(N^\alpha)$; this, in +turn, leads to a total growth in effort as ${\cal O}(N^{1+\alpha})$ +since each iteration takes ${\cal O}(N)$ work. This behavior is +undesirable, we would really like to solve linear systems with $N$ +unknowns in a total of ${\cal O}(N)$ work; there is a class +of preconditioners that can achieve this, namely geometric (@ref +step_16 "step-16") or algebraic (@ref step_31 "step-31") multigrid +preconditioners. They are, however, significantly more complex than +the preconditioners outlined above. Finally, the last message to take +home is that today (in 2008), linear systems with 100,000 unknowns are +easily solved on a desktop machine in well under 10 seconds, making +the solution of relatively simple 2d problems even to very high +accuracy not that big a task as it used to be even in the recent +past. Of course, the situation for 3d problems is entirely different. +