]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add a comparison of solvers/preconditioners.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Sep 2008 16:52:39 +0000 (16:52 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Sep 2008 16:52:39 +0000 (16:52 +0000)
git-svn-id: https://svn.dealii.org/trunk@16902 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-6/doc/results.dox

index a639c71e5b39d2a53dc46f139cffe55054d42ea4..8f442eb90cf191ad0da01bdb4499a7c03bfe066d 100644 (file)
@@ -161,3 +161,111 @@ through a few more functions from the destructor of the
 <code>LaplaceProblem</code> class, exactly where we have commented out
 the call to <code>DoFHandler::clear()</code>.
 
+
+
+<a name="extensions"></a>
+<h3>Possibilities for extensions</h3>
+
+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 <code>solve()</code> 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 <code>lac/sparse_ilu.h</code> 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<double> 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
+<code>solver_control.last_step()</code> 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):
+
+<TABLE WIDTH="60%" ALIGN="center">
+  <tr>
+    <td ALIGN="center">
+      @image html step-6.q2.dofs_vs_iterations.png
+    </td>
+
+    <td ALIGN="center">
+      @image html step-6.q2.dofs_vs_time.png
+    </td>
+  </tr>
+</table>
+
+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:
+
+<TABLE WIDTH="60%" ALIGN="center">
+  <tr>
+    <td ALIGN="center">
+      @image html step-6.q1.dofs_vs_iterations.png
+    </td>
+
+    <td ALIGN="center">
+      @image html step-6.q1.dofs_vs_time.png
+    </td>
+  </tr>
+</table>
+
+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.
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.