From: Uwe Köcher Date: Wed, 17 Sep 2014 07:33:13 +0000 (+0200) Subject: Updated Documentation. Again, in the Introduction section of the step-20 tutorial... X-Git-Tag: v8.2.0-rc1~143^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F151%2Fhead;p=dealii.git Updated Documentation. Again, in the Introduction section of the step-20 tutorial a class InverseMatrix appears. Now its changed to use IterativeInverse<..> in all places to match with the code. --- diff --git a/examples/step-20/doc/intro.dox b/examples/step-20/doc/intro.dox index 74b1ae3569..dbea94c603 100644 --- a/examples/step-20/doc/intro.dox +++ b/examples/step-20/doc/intro.dox @@ -483,7 +483,13 @@ matrix and the mass matrix, respectively: template void MixedLaplaceProblem::solve () { - const InverseMatrix m_inverse (system_matrix.block(0,0)); + PreconditionIdentity identity; + IterativeInverse > m_inverse; + m_inverse.initialize(system_matrix.block(0,0), identity); + m_inverse.solver.select("cg"); + static ReductionControl inner_control(1000, 0., 1.e-13); + m_inverse.solver.set_control(inner_control); + Vector tmp (solution.block(0).size()); { @@ -600,8 +606,11 @@ the approximate Schur complement, i.e. we need code like this: ApproximateSchurComplement approximate_schur_complement (system_matrix); - InverseMatrix - preconditioner (approximate_schur_complement) + IterativeInverse > + preconditioner; + preconditioner.initialize(approximate_schur_complement, identity); + preconditioner.solver.select("cg"); + preconditioner.solver.set_control(inner_control); @endcode That's all! @@ -622,12 +631,15 @@ look like this: ApproximateSchurComplement approximate_schur_complement (system_matrix); - InverseMatrix - preconditioner (approximate_schur_complement); + IterativeInverse > + preconditioner; + preconditioner.initialize(approximate_schur_complement, identity); + preconditioner.solver.select("cg"); + preconditioner.solver.set_control(inner_control); - SolverControl solver_control (system_matrix.block(0,0).m(), - 1e-6*schur_rhs.l2_norm()); - SolverCG<> cg (solver_control); + SolverControl solver_control (solution.block(1).size(), + 1e-12*schur_rhs.l2_norm()); + SolverCG<> cg (solver_control); cg.solve (schur_complement, solution.block(1), schur_rhs, preconditioner);