From: Matthias Maier Date: Tue, 12 Feb 2019 00:42:10 +0000 (-0600) Subject: Address Wolfgang's comments X-Git-Tag: v9.1.0-rc1~347^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c09b8083217dbf2a06ba4af9b6223a45f3218bed;p=dealii.git Address Wolfgang's comments - fix wording - remove author statement (the generic copyright is enough). --- diff --git a/examples/step-20/doc/intro.dox b/examples/step-20/doc/intro.dox index 97bcb989cf..7a836f0edf 100644 --- a/examples/step-20/doc/intro.dox +++ b/examples/step-20/doc/intro.dox @@ -384,7 +384,8 @@ the system. We note that the resulting solver is not optimal -- there are much better ways to efficiently compute the system, for example those explained in the results section of step-22 or the one we use in step-43 for a problem similar to the current one. Here, our goal shall be to -introduce new solution techniques how they can be implemented in deal.II. +introduce new solution techniques and how they can be implemented in +deal.II.

Solving using the Schur complement

@@ -494,8 +495,8 @@ factor (here, $10^{-10}$). In contrast the SolverControl class only checks for absolute tolerances. We have to use ReductionControl in our case to work around a minor issue: The right hand sides that we will feed to op_M_inv are essentially formed by residuals that naturally -have vastly differing norms. This makes control by an absolute tolerance -very error prone. +decrease vastly in norm as the outer iterations progress. This makes +control by an absolute tolerance very error prone. We now have a LinearOperator op_M_inv that we can use to construct more complicated operators such as the Schur complement $S$. @@ -515,7 +516,11 @@ similar to the following code: solver_M(M, tmp2, tmp1, preconditioner_M); // multiply with M^-1 B.Tvmult (dst, tmp2); // multiply with the bottom left block: B^T @endcode -(tmp1 and tmp2 are two temporary vectors). +(tmp1 and tmp2 are two temporary vectors). The +key point behind this approach is the fact that we never actually create an +inner product of matrices. Instead, whenever we have to perform a matrix +vector multiplication with op_S we simply run all individual +vmult operations in above sequence. @note We could have achieved the same goal of creating a "matrix like" object by implementing a specialized class SchurComplement @@ -539,9 +544,11 @@ class SchurComplement @endcode Even though both approaches are exactly equivalent, the LinearOperator class has a big advantage over this manual approach. -It provides so-called syntactic sugar: Mathematically, we think -about $S$ as being the composite matrix $S=B^TM^{-1}B$ and the -LinearOperator class allows you to write this out more or less verbatim, +It provides so-called +syntactic sugar: +Mathematically, we think about $S$ as being the composite matrix +$S=B^TM^{-1}B$ and the LinearOperator class allows you to write this out +more or less verbatim, @code const auto op_M_inv = inverse_operator(op_M, solver_M, preconditioner_M); const auto op_S = transpose_operator(op_B) * op_M_inv * op_B; @@ -570,13 +577,14 @@ us. Similarly in spirit to LinearOperator, a PackagedOperation stores a std::function apply; std::function apply_add; @endcode -The class allows lazy evaluation of expressions involving vectors and -linear operators. This is done by storing the computational expression and -only performing the computation when either the object is converted to a -vector object, or PackagedOperation::apply() (or -PackagedOperation::apply_add()) is invoked by hand. Assuming that -F and G are the two vectors of the right hand -side we can simply write: +The class allows +lazy evaluation +of expressions involving vectors and linear operators. This is done by +storing the computational expression and only performing the computation +when either the object is converted to a vector object, or +PackagedOperation::apply() (or PackagedOperation::apply_add()) is invoked +by hand. Assuming that F and G are the two +vectors of the right hand side we can simply write: @code const auto schur_rhs = transpose_operator(op_B) * op_M_inv * F - G; @endcode diff --git a/examples/step-20/step-20.cc b/examples/step-20/step-20.cc index 2c4ef1733f..5879d3b6a1 100644 --- a/examples/step-20/step-20.cc +++ b/examples/step-20/step-20.cc @@ -1,6 +1,6 @@ /* --------------------------------------------------------------------- * - * Copyright (C) 2005 - 2018 by the deal.II authors + * Copyright (C) 2005 - 2019 by the deal.II authors * * This file is part of the deal.II library. * @@ -12,10 +12,6 @@ * the top level directory of deal.II. * * --------------------------------------------------------------------- - - * - * Authors: Wolfgang Bangerth, Texas A&M University, 2005, 2006; - * (port to LinearOperator:) Matthias Maier, 2019 */ @@ -607,16 +603,15 @@ namespace Step20 const auto op_M_inv = inverse_operator(op_M, solver_M, preconditioner_M); - // This puts us in the position to be able to declare the Schur - // complement op_S and the approximate Schur complement - // op_aS: + // This allows us to declare the Schur complement op_S and + // the approximate Schur complement op_aS: const auto op_S = transpose_operator(op_B) * op_M_inv * op_B; const auto op_aS = transpose_operator(op_B) * linear_operator(preconditioner_M) * op_B; // We now create a preconditioner out of op_aS that - // applies a few CG iterations (until a very modest relative reduction - // of $10^{-3}$ is reached): + // applies a small number of CG iterations (until a very modest + // relative reduction of $10^{-3}$ is reached): ReductionControl reduction_control_aS(2000, 1.e-18, 1.0e-3); SolverCG<> solver_aS(reduction_control_aS); PreconditionIdentity preconditioner_aS;