 |
Reference documentation for deal.II version 8.5.1
|
|
Reference documentation for deal.II version 8.5.1
|
Creation of a BlockLinearOperator | |
| template<typename Range , typename Domain , typename BlockPayload , typename BlockMatrixType > | |
| BlockLinearOperator< Range, Domain, BlockPayload > | block_operator (const BlockMatrixType &block_matrix) |
| template<size_t m, size_t n, typename Range , typename Domain , typename BlockPayload > | |
| BlockLinearOperator< Range, Domain, BlockPayload > | block_operator (const std::array< std::array< LinearOperator< typename Range::BlockType, typename Domain::BlockType, typename BlockPayload::BlockType >, n >, m > &ops) |
| template<typename Range , typename Domain , typename BlockPayload , typename BlockMatrixType > | |
| BlockLinearOperator< Range, Domain, BlockPayload > | block_diagonal_operator (const BlockMatrixType &block_matrix) |
| template<size_t m, typename Range , typename Domain , typename BlockPayload > | |
| BlockLinearOperator< Range, Domain, BlockPayload > | block_diagonal_operator (const std::array< LinearOperator< typename Range::BlockType, typename Domain::BlockType, typename BlockPayload::BlockType >, m > &ops) |
| template<size_t m, typename Range , typename Domain , typename BlockPayload > | |
| BlockLinearOperator< Range, Domain, BlockPayload > | block_diagonal_operator (const LinearOperator< typename Range::BlockType, typename Domain::BlockType, typename BlockPayload::BlockType > &op) |
Manipulation of a BlockLinearOperator | |
| template<typename Range , typename Domain , typename BlockPayload > | |
| LinearOperator< Domain, Range, typename BlockPayload::BlockType > | block_forward_substitution (const BlockLinearOperator< Range, Domain, BlockPayload > &block_operator, const BlockLinearOperator< Domain, Range, BlockPayload > &diagonal_inverse) |
| template<typename Range , typename Domain , typename BlockPayload > | |
| LinearOperator< Domain, Range, typename BlockPayload::BlockType > | block_back_substitution (const BlockLinearOperator< Range, Domain, BlockPayload > &block_operator, const BlockLinearOperator< Domain, Range, BlockPayload > &diagonal_inverse) |
Vector space operations | |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | operator+ (const LinearOperator< Range, Domain, Payload > &first_op, const LinearOperator< Range, Domain, Payload > &second_op) |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | operator- (const LinearOperator< Range, Domain, Payload > &first_op, const LinearOperator< Range, Domain, Payload > &second_op) |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | operator* (typename Range::value_type number, const LinearOperator< Range, Domain, Payload > &op) |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | operator* (const LinearOperator< Range, Domain, Payload > &op, typename Domain::value_type number) |
Composition and manipulation of a LinearOperator | |
| template<typename Range , typename Intermediate , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | operator* (const LinearOperator< Range, Intermediate, Payload > &first_op, const LinearOperator< Intermediate, Domain, Payload > &second_op) |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Domain, Range, Payload > | transpose_operator (const LinearOperator< Range, Domain, Payload > &op) |
| template<typename Payload , typename Solver , typename Preconditioner , typename Range = typename Solver::vector_type, typename Domain = Range> | |
| LinearOperator< Domain, Range, Payload > | inverse_operator (const LinearOperator< Range, Domain, Payload > &op, Solver &solver, const Preconditioner &preconditioner) |
Creation of a LinearOperator | |
| template<typename Range , typename Payload = internal::LinearOperator::EmptyPayload> | |
| LinearOperator< Range, Range, Payload > | identity_operator (const std::function< void(Range &, bool)> &reinit_vector) |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | identity_operator (const LinearOperator< Range, Domain, Payload > &op) |
| template<typename Range , typename Domain , typename Payload > | |
| LinearOperator< Range, Domain, Payload > | null_operator (const LinearOperator< Range, Domain, Payload > &op) |
| template<typename Range , typename Domain , typename Payload , typename Matrix > | |
| LinearOperator< Range, Domain, Payload > | linear_operator (const Matrix &matrix) |
| template<typename Range , typename Domain , typename Payload , typename OperatorExemplar , typename Matrix > | |
| LinearOperator< Range, Domain, Payload > | linear_operator (const OperatorExemplar &operator_exemplar, const Matrix &matrix) |
Vector space operations | |
| template<typename Range > | |
| PackagedOperation< Range > | operator+ (const PackagedOperation< Range > &first_comp, const PackagedOperation< Range > &second_comp) |
| template<typename Range > | |
| PackagedOperation< Range > | operator- (const PackagedOperation< Range > &first_comp, const PackagedOperation< Range > &second_comp) |
| template<typename Range > | |
| PackagedOperation< Range > | operator* (const PackagedOperation< Range > &comp, typename Range::value_type number) |
| template<typename Range > | |
| PackagedOperation< Range > | operator* (typename Range::value_type number, const PackagedOperation< Range > &comp) |
| template<typename Range > | |
| PackagedOperation< Range > | operator+ (const PackagedOperation< Range > &comp, const Range &offset) |
| template<typename Range > | |
| PackagedOperation< Range > | operator+ (const Range &offset, const PackagedOperation< Range > &comp) |
| template<typename Range > | |
| PackagedOperation< Range > | operator- (const PackagedOperation< Range > &comp, const Range &offset) |
| template<typename Range > | |
| PackagedOperation< Range > | operator- (const Range &offset, const PackagedOperation< Range > &comp) |
Creation of a PackagedOperation object | |
| template<typename Range , typename = typename std::enable_if<has_vector_interface<Range>::type::value>::type> | |
| PackagedOperation< Range > | operator+ (const Range &u, const Range &v) |
| template<typename Range , typename = typename std::enable_if<has_vector_interface<Range>::type::value>::type> | |
| PackagedOperation< Range > | operator- (const Range &u, const Range &v) |
| template<typename Range , typename = typename std::enable_if<has_vector_interface<Range>::type::value>::type> | |
| PackagedOperation< Range > | operator* (const Range &u, typename Range::value_type number) |
| template<typename Range , typename = typename std::enable_if<has_vector_interface<Range>::type::value>::type> | |
| PackagedOperation< Range > | operator* (typename Range::value_type number, const Range &u) |
| template<typename Range , typename Domain , typename Payload > | |
| PackagedOperation< Range > | operator* (const LinearOperator< Range, Domain, Payload > &op, const Domain &u) |
| template<typename Range , typename Domain , typename Payload > | |
| PackagedOperation< Domain > | operator* (const Range &u, const LinearOperator< Range, Domain, Payload > &op) |
| template<typename Range , typename Domain , typename Payload > | |
| PackagedOperation< Range > | operator* (const LinearOperator< Range, Domain, Payload > &op, const PackagedOperation< Domain > &comp) |
| template<typename Range , typename Domain , typename Payload > | |
| PackagedOperation< Domain > | operator* (const PackagedOperation< Range > &comp, const LinearOperator< Range, Domain, Payload > &op) |
Creation of a LinearOperator related to the Schur Complement | |
| template<typename Range_1 , typename Domain_1 , typename Range_2 , typename Domain_2 , typename Payload > | |
| LinearOperator< Range_2, Domain_2, Payload > | schur_complement (const LinearOperator< Domain_1, Range_1, Payload > &A_inv, const LinearOperator< Range_1, Domain_2, Payload > &B, const LinearOperator< Range_2, Domain_1, Payload > &C, const LinearOperator< Range_2, Domain_2, Payload > &D) |
Creation of PackagedOperation objects related to the Schur Complement | |
| template<typename Range_1 , typename Domain_1 , typename Range_2 , typename Payload > | |
| PackagedOperation< Range_2 > | condense_schur_rhs (const LinearOperator< Range_1, Domain_1, Payload > &A_inv, const LinearOperator< Range_2, Domain_1, Payload > &C, const Range_1 &f, const Range_2 &g) |
| template<typename Range_1 , typename Domain_1 , typename Domain_2 , typename Payload > | |
| PackagedOperation< Domain_1 > | postprocess_schur_solution (const LinearOperator< Range_1, Domain_1, Payload > &A_inv, const LinearOperator< Range_1, Domain_2, Payload > &B, const Domain_2 &y, const Range_1 &f) |
If deal.II is configured with C++11 support (i.e., DEAL_II_WITH_CXX11=ON or DEAL_II_WITH_CXX14=ON during configuration) a versatile mechanism for storing the concept of a linear operator is available. (For questions about C++11, see deal.II and the C++11 standard .)
This is done with a LinearOperator class that, like the MatrixType concept, defines a minimal interface for applying a linear operation on a vector.
Thus, such an object can be used as a matrix object in all iterative solver classes, either as a matrix object, or as preconditioner.
The big advantage of the LinearOperator class is that it provides syntactic sugar for complex matrix-vector operations. As an example consider the operation \((A+k\,B)\,C\), where \(A\), \(B\) and \(C\) denote (possibly different) SparseMatrix objects. In order to construct a LinearOperator op that performs above computation when applied on a vector, one can write:
Now, op can be used as a matrix object for further computation.
The linear_operator() function can be used to wrap an ordinary matrix or preconditioner object into a LinearOperator. A linear operator can be transposed with transpose_operator(), or inverted by using the inverse_operator() together with an iterative solver.
For objects of type LinearOperator, all vector space operations, i.e., addition and subtraction, scalar multiplication and composition (of compatible linear operators) are implemented:
block_operator() and block_diagonal_operator() provide further encapsulation of individual linear operators into blocked linear operator variants.
res = op_a*x a PackagedOperation class instance is generated behind-the-curtains. Consequently, the user program has to include header files for both classes for compilation to be successful. In an attempt to make easier the decision of which headers to include in what circumstances and to prevent hidden templates-related compiler errors, all headers relevant to LinearOperator are grouped in <deal.ii/lac/linear_operator_tools.h>.An application of a LinearOperator object to a vector via operator* yields a PackagedOperation object that stores this computation.
The PackagedOperation 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 implicitly converted to a vector object, or PackagedOperation::apply() (or PackagedOperation::apply_add()) is invoked by hand. This avoids unnecessary temporary storage of intermediate results.
As an example consider the addition of multiple vectors:
Converting the PackagedOperation a + b - c + d to a vector results in code equivalent to the following code
that avoids any intermediate storage. As a second example (involving a LinearOperator object) consider the computation of a residual \(b-Ax\):
Here, the expression b - op_a * x results again in an object of type PackagedOperation that stores the sequence of operations that should be performed using the two vectors and the linear operator. Converting the expression to a vector (as happens here with the assignment to the vector residual) executes the computation (see the following note).
residual_expr object b and matrix A. It does not perform any computation at this point. In particular, if b or A are changed after the creation of residual_expr every subsequent evaluation of the expression is performed with the new values auto):
|
related |
A function that encapsulates a block_matrix into a BlockLinearOperator.
All changes made on the block structure and individual blocks of block_matrix after the creation of the BlockLinearOperator object are reflected by the operator object.
Definition at line 441 of file block_linear_operator.h.
|
related |
A variant of above function that encapsulates a given collection ops of LinearOperators into a block structure. Here, it is assumed that Range and Domain are blockvectors, i.e., derived from BlockVectorBase. The individual linear operators in ops must act on the underlying vector type of the block vectors, i.e., on Domain::BlockType yielding a result in Range::BlockType.
The list ops is best passed as an initializer list. Consider for example a linear operator block (acting on Vector<double>)
The corresponding block_operator invocation takes the form
Definition at line 504 of file block_linear_operator.h.
|
related |
This function extracts the diagonal blocks of block_matrix (either a block matrix type or a BlockLinearOperator) and creates a BlockLinearOperator with the diagonal. Off-diagonal elements are initialized as null_operator (with correct reinit_range_vector and reinit_domain_vector methods).
All changes made on the individual diagonal blocks of block_matrix after the creation of the BlockLinearOperator object are reflected by the operator object.
Definition at line 557 of file block_linear_operator.h.
|
related |
A variant of above function that builds up a block diagonal linear operator from an array ops of diagonal elements (off-diagonal blocks are assumed to be 0).
The list ops is best passed as an initializer list. Consider for example a linear operator block (acting on Vector<double>) diag(op_a0, op_a1, ..., op_am). The corresponding block_operator invocation takes the form
Definition at line 613 of file block_linear_operator.h.
|
related |
A variant of above function that only takes a single LinearOperator argument op and creates a blockdiagonal linear operator with m copies of it.
Definition at line 657 of file block_linear_operator.h.
|
related |
This function implements forward substitution to invert a lower block triangular matrix. As arguments, it takes a BlockLinearOperator block_operator representing a block lower triangular matrix, as well as a BlockLinearOperator diagonal_inverse representing inverses of diagonal blocks of block_operator.
Let us assume we have a linear system with the following block structure:
First of all, x0 = A00^-1 y0. Then, we can use x0 to recover x1:
and therefore:
block_operator is used as well as the diagonal of diagonal_inverse. Definition at line 715 of file block_linear_operator.h.
|
related |
This function implements back substitution to invert an upper block triangular matrix. As arguments, it takes a BlockLinearOperator block_operator representing an upper block triangular matrix, as well as a BlockLinearOperator diagonal_inverse representing inverses of diagonal blocks of block_operator.
Let us assume we have a linear system with the following block structure:
First of all, xn = Ann^-1 yn. Then, we can use xn to recover x(n-1):
and therefore:
block_operator is used as well as the diagonal of diagonal_inverse. Definition at line 827 of file block_linear_operator.h.
|
related |
Addition of two linear operators first_op and second_op given by \((\text{first\_op}+\text{second\_op})x := \text{first\_op}(x) + \text{second\_op}(x)\)
Definition at line 365 of file linear_operator.h.
|
related |
Subtraction of two linear operators first_op and second_op given by \((\text{first\_op}-\text{second\_op})x := \text{first\_op}(x) - \text{second\_op}(x)\)
Definition at line 428 of file linear_operator.h.
|
related |
Scalar multiplication of a ScalarOperator object op with number from the left.
The Domain and Range types must implement the following operator*= member functions accepting the appropriate scalar Number type for rescaling:
Definition at line 466 of file linear_operator.h.
|
related |
Scalar multiplication of a ScalarOperator object from the right.
The Domain and Range types must implement the following operator*= member functions for rescaling:
Definition at line 536 of file linear_operator.h.
|
related |
Composition of two linear operators first_op and second_op given by \((\text{first\_op}*\text{second\_op})x := \text{first\_op}(\text{second\_op}(x))\)
Definition at line 565 of file linear_operator.h.
|
related |
Return the transpose linear operations of op.
Definition at line 646 of file linear_operator.h.
|
related |
Return an object representing the inverse of the LinearOperator op.
The function takes references solver and preconditioner to an iterative solver and a preconditioner that are used in the vmult and Tvmult implementations of the LinearOperator object.
The LinearOperator object that is created stores a reference to solver and preconditioner. Thus, both objects must remain a valid reference for the whole lifetime of the LinearOperator object. Internal data structures of the solver object will be modified upon invocation of vmult or Tvmult.
Definition at line 688 of file linear_operator.h.
|
related |
Return a LinearOperator that is the identity of the vector space Range.
The function takes an std::function object reinit_vector as an argument to initialize the reinit_range_vector and reinit_domain_vector objects of the LinearOperator object.
Definition at line 760 of file linear_operator.h.
|
related |
Return a LinearOperator that is the identity of the vector space Range.
The function takes a LinearOperator op and uses its range initializer to create an identiy operator. In contrast to the function above, this function also ensures that the underlying Payload matches that of the input op.
Definition at line 805 of file linear_operator.h.
|
related |
Return a nulled variant of the LinearOperator op, i.e. with optimized LinearOperator::vmult, LinearOperator::vmult_add, etc. functions and with LinearOperator::is_null_operator set to true.
Definition at line 835 of file linear_operator.h.
|
related |
A function that encapsulates generic matrix objects that act on a compatible Vector type into a LinearOperator. The LinearOperator object that is created stores a reference to the matrix object. Thus, matrix must remain a valid reference for the whole lifetime of the LinearOperator object.
All changes made on matrix after the creation of the LinearOperator object are reflected by the operator object. For example, it is a valid procedure to first create a LinearOperator and resize, reassemble the matrix later.
The Matrix class in question must provide the following minimal interface:
The following (optional) interface is used if available:
If the Matrix does not provide vmult_add and Tvmult_add, they are implemented in terms of vmult and Tvmult (requiring intermediate storage).
Definition at line 1242 of file linear_operator.h.
|
related |
Variant of above function that takes an operator object operator_exemplar as an additional reference. This object is used to populate the reinit_domain_vector and reinit_range_vector function objects. The reference matrix is used to construct vmult, Tvmult, etc.
This variant can, for example, be used to encapsulate preconditioners (that typically do not expose any information about the underlying matrix).
Definition at line 1270 of file linear_operator.h.
|
related |
Addition of two PackagedOperation objects first_comp and second_comp given by vector space addition of the corresponding results.
Definition at line 289 of file packaged_operation.h.
|
related |
Subtraction of two PackagedOperation objects first_comp and second_comp given by vector space addition of the corresponding results.
Definition at line 324 of file packaged_operation.h.
|
related |
Scalar multiplication of a PackagedOperation objects comp with a scalar number given by a scaling PackagedOperation result with number.
Definition at line 362 of file packaged_operation.h.
|
related |
Scalar multiplication of a PackagedOperation objects comp with a scalar number given by a scaling PackagedOperation result with number.
Definition at line 410 of file packaged_operation.h.
|
related |
Add a constant offset (of the Range space) to the result of a PackagedOperation.
Definition at line 425 of file packaged_operation.h.
|
related |
Add a constant offset (of the Range space) to the result of a PackagedOperation.
Definition at line 440 of file packaged_operation.h.
|
related |
Subtract a constant offset (of the Range space) from the result of a PackagedOperation.
Definition at line 455 of file packaged_operation.h.
|
related |
Subtract a computational result from a constant offset (of the Range space). The result is a PackagedOperation object that applies this computation.
Definition at line 472 of file packaged_operation.h.
|
related |
Create a PackagedOperation object that stores the addition of two vectors.
The PackagedOperation object that is created stores a reference to u and v. Thus, the vectors must remain valid references for the whole lifetime of the PackagedOperation object. All changes made on u or v after the creation of the PackagedOperation object are reflected by the operator object.
Definition at line 528 of file packaged_operation.h.
|
related |
Create a PackagedOperation object that stores the subtraction of two vectors.
The PackagedOperation object that is created stores a reference to u and v. Thus, the vectors must remain valid references for the whole lifetime of the PackagedOperation object. All changes made on u or v after the creation of the PackagedOperation object are reflected by the operator object.
Definition at line 573 of file packaged_operation.h.
|
related |
Create a PackagedOperation object that stores the scaling of a vector with a number.
The PackagedOperation object that is created stores a reference to u. Thus, the vectors must remain valid references for the whole lifetime of the PackagedOperation object. All changes made on u or v after the creation of the PackagedOperation object are reflected by the operator object.
Definition at line 617 of file packaged_operation.h.
|
related |
Create a PackagedOperation object that stores the scaling of a vector with a number.
The PackagedOperation object that is created stores a reference to u. Thus, the vectors must remain valid references for the whole lifetime of the PackagedOperation object. All changes made on u or v after the creation of the PackagedOperation object are reflected by the operator object.
Definition at line 640 of file packaged_operation.h.
|
related |
Create a PackagedOperation object from a LinearOperator and a reference to a vector u of the Domain space. The object stores the PackagedOperation \(\text{op} \,u\) (in matrix notation). return (return_add) are implemented with vmult(__1,u) (vmult_add(__1,u)).
The PackagedOperation object that is created stores a reference to u. Thus, the vector must remain a valid reference for the whole lifetime of the PackagedOperation object. All changes made on u after the creation of the PackagedOperation object are reflected by the operator object.
Definition at line 665 of file packaged_operation.h.
|
related |
Create a PackagedOperation object from a LinearOperator and a reference to a vector u of the Range space. The object stores the PackagedOperation \(\text{op}^T \,u\) (in matrix notation). return (return_add) are implemented with Tvmult(__1,u) (Tvmult_add(__1,u)).
The PackagedOperation object that is created stores a reference to u. Thus, the vector must remain a valid reference for the whole lifetime of the PackagedOperation object. All changes made on u after the creation of the PackagedOperation object are reflected by the operator object.
Definition at line 707 of file packaged_operation.h.
|
related |
Composition of a PackagedOperation object with a LinearOperator. The object stores the computation \(\text{op} \,comp\) (in matrix notation).
Definition at line 741 of file packaged_operation.h.
|
related |
Composition of a PackagedOperation object with a LinearOperator. The object stores the computation \(\text{op}^T \,comp\) (in matrix notation).
Definition at line 791 of file packaged_operation.h.
|
related |
Return a LinearOperator that performs the operations associated with the Schur complement. There are two additional helper functions, condense_schur_rhs() and postprocess_schur_solution(), that are likely necessary to be used in order to perform any useful tasks in linear algebra with this operator.
We construct the definition of the Schur complement in the following way:
Consider a general system of linear equations that can be decomposed into two major sets of equations:
\begin{eqnarray*} \mathbf{K}\mathbf{d} = \mathbf{f} \quad \Rightarrow\quad \left(\begin{array}{cc} A & B \\ C & D \end{array}\right) \left(\begin{array}{cc} x \\ y \end{array}\right) = \left(\begin{array}{cc} f \\ g \end{array}\right), \end{eqnarray*}
where \( A,B,C,D \) represent general subblocks of the matrix \( \mathbf{K} \) and, similarly, general subvectors of \( \mathbf{d},\mathbf{f} \) are given by \( x,y,f,g \) .
This is equivalent to the following two statements:
\begin{eqnarray*} (1) \quad Ax + By &=& f \\ (2) \quad Cx + Dy &=& g \quad . \end{eqnarray*}
Assuming that \( A,D \) are both square and invertible, we could then perform one of two possible substitutions,
\begin{eqnarray*} (3) \quad x &=& A^{-1}(f - By) \quad \text{from} \quad (1) \\ (4) \quad y &=& D^{-1}(g - Cx) \quad \text{from} \quad (2) , \end{eqnarray*}
which amount to performing block Gaussian elimination on this system of equations.
For the purpose of the current implementation, we choose to substitute (3) into (2)
\begin{eqnarray*} C \: A^{-1}(f - By) + Dy &=& g \\ -C \: A^{-1} \: By + Dy &=& g - C \: A^{-1} \: f \quad . \end{eqnarray*}
This leads to the result
\[ (5) \quad (D - C\: A^{-1} \:B)y = g - C \: A^{-1} f \quad \Rightarrow \quad Sy = g' \]
with \( S = (D - C\: A^{-1} \:B) \) being the Schur complement and the modified right-hand side vector \( g' = g - C \: A^{-1} f \) arising from the condensation step. Note that for this choice of \( S \), submatrix \( D \) need not be invertible and may thus be the null matrix. Ideally \( A \) should be well-conditioned.
So for any arbitrary vector \( a \), the Schur complement performs the following operation:
\[ (6) \quad Sa = (D - C \: A^{-1} \: B)a \]
A typical set of steps needed the solve a linear system (1),(2) would be:
A_inv (using inverse_operator()).\[ g' = g - C \: A^{-1} \: f \]
\[ y = S^{-1} g' \]
\[ x = A^{-1} (f - By) \]
An illustration of typical usage of this operator for a fully coupled system is given below.
In the above example, the preconditioner for \( S \) was defined as the preconditioner for \( D \), which is valid since they operate on the same space. However, if \( D \) and \( S \) are too dissimilar, then this may lead to a large number of solver iterations as \( \text{prec}(D) \) is not a good approximation for \( S^{-1} \).
A better preconditioner in such a case would be one that provides a more representative approximation for \( S^{-1} \). One approach is shown in step-22, where \( D \) is the null matrix and the preconditioner for \( S^{-1} \) is derived from the mass matrix over this space.
From another viewpoint, a similar result can be achieved by first constructing an object that represents an approximation for \( S \) wherein expensive operation, namely \( A^{-1} \), is approximated. Thereafter we construct the approximate inverse operator \( \tilde{S}^{-1} \) which is then used as the preconditioner for computing \( S^{-1} \).
Note that due to the construction of S_inv_approx and subsequently S_inv, there are a pair of nested iterative solvers which could collectively consume a lot of resources. Therefore care should be taken in the choices leading to the construction of the iterative inverse_operators. One might consider the use of a IterationNumberControl (or a similar mechanism) to limit the number of inner solver iterations. This controls the accuracy of the approximate inverse operation \( \tilde{S}^{-1} \) which acts only as the preconditioner for \( S^{-1} \). Furthermore, the preconditioner to \( \tilde{S}^{-1} \), which in this example is \( \text{prec}(D) \), should ideally be computationally inexpensive.
However, if an iterative solver based on IterationNumberControl is used as a preconditioner then the preconditioning operation is not a linear operation. Here a flexible solver like SolverFGMRES (flexible GMRES) is best employed as an outer solver in order to deal with the variable behaviour of the preconditioner. Otherwise the iterative solver can stagnate somewhere near the tolerance of the preconditioner or generally behave erratically. Alternatively, using a ReductionControl would ensure that the preconditioner always solves to the same tolerance, thereby rendering its behaviour constant.
Further examples of this functionality can be found in the test-suite, such as tests/lac/schur_complement_01.cc . The solution of a multi- component problem (namely step-22) using the schur_complement can be found in tests/lac/schur_complement_03.cc .
Definition at line 230 of file schur_complement.h.
|
related |
For the system of equations
\begin{eqnarray*} Ax + By &=& f \\ Cx + Dy &=& g \quad , \end{eqnarray*}
this operation performs the pre-processing (condensation) step on the RHS subvector g so that the Schur complement can be used to solve this system of equations. More specifically, it produces an object that represents the condensed form of the subvector g, namely
\[ g' = g - C \: A^{-1} \: f \]
Definition at line 279 of file schur_complement.h.
|
related |
For the system of equations
\begin{eqnarray*} Ax + By &=& f \\ Cx + Dy &=& g \quad , \end{eqnarray*}
this operation performs the post-processing step of the Schur complement to solve for the second subvector x once subvector y is known, with the result that
\[ x = A^{-1}(f - By) \]
Definition at line 316 of file schur_complement.h.
1.8.14