Ifpack2 Templated Preconditioning Package
Version 1.0
|
Additive Schwarz domain decomposition for Tpetra sparse matrices. More...
#include <Ifpack2_AdditiveSchwarz_decl.hpp>
Public Types | |
Typedefs | |
using | scalar_type = typename MatrixType::scalar_type |
The type of the entries of the input MatrixType. More... | |
using | local_ordinal_type = typename MatrixType::local_ordinal_type |
The type of local indices in the input MatrixType. More... | |
using | global_ordinal_type = typename MatrixType::global_ordinal_type |
The type of global indices in the input MatrixType. More... | |
using | node_type = typename MatrixType::node_type |
The Node type used by the input MatrixType. More... | |
using | magnitude_type = typename Teuchos::ScalarTraits< scalar_type >::magnitudeType |
The type of the magnitude (absolute value) of a matrix entry. More... | |
using | row_matrix_type = Tpetra::RowMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > |
The Tpetra::RowMatrix specialization matching MatrixType. More... | |
using | crs_matrix_type = Tpetra::CrsMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > |
The Tpetra::CrsMatrix specialization that is a subclass of MatrixType. More... | |
Public Types inherited from Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > | |
typedef Teuchos::ScalarTraits < MatrixType::scalar_type > ::magnitudeType | magnitude_type |
The type of the magnitude (absolute value) of a matrix entry. More... | |
Public Member Functions | |
virtual Teuchos::RCP< const row_matrix_type > | getMatrix () const |
The input matrix. More... | |
virtual void | setParameters (const Teuchos::ParameterList &plist) |
Set the preconditioner's parameters. More... | |
void | setParameterList (const Teuchos::RCP< Teuchos::ParameterList > &plist) |
Set the preconditioner's parameters. More... | |
void | setZeroStartingSolution (bool zeroStartingSolution) |
Set this preconditioner's parameters. More... | |
Teuchos::RCP< const Teuchos::ParameterList > | getValidParameters () const |
Get a list of the preconditioner's default parameters. More... | |
virtual void | initialize () |
Computes all (graph-related) data necessary to initialize the preconditioner. More... | |
virtual bool | isInitialized () const |
Returns true if the preconditioner has been successfully initialized, false otherwise. More... | |
virtual void | compute () |
Computes all (coefficient) data necessary to apply the preconditioner. More... | |
virtual bool | isComputed () const |
Returns true if the preconditioner has been successfully computed, false otherwise. More... | |
virtual int | getNumInitialize () const |
Returns the number of calls to initialize(). More... | |
virtual int | getNumCompute () const |
Returns the number of calls to compute(). More... | |
virtual int | getNumApply () const |
Returns the number of calls to apply(). More... | |
virtual double | getInitializeTime () const |
Returns the time spent in initialize(). More... | |
virtual double | getComputeTime () const |
Returns the time spent in compute(). More... | |
virtual double | getApplyTime () const |
Returns the time spent in apply(). More... | |
virtual std::ostream & | print (std::ostream &os) const |
Prints basic information on iostream. This function is used by operator<<. More... | |
virtual int | getOverlapLevel () const |
Returns the level of overlap. More... | |
AdditiveSchwarz (const Teuchos::RCP< const row_matrix_type > &A) | |
Constructor that takes a matrix. More... | |
AdditiveSchwarz (const Teuchos::RCP< const row_matrix_type > &A, const int overlapLevel) | |
Constructor that takes a matrix and the level of overlap. More... | |
virtual | ~AdditiveSchwarz ()=default |
Destructor. More... | |
Implementation of Tpetra::Operator | |
virtual Teuchos::RCP< const Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > > | getDomainMap () const |
The domain Map of this operator. More... | |
virtual Teuchos::RCP< const Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > > | getRangeMap () const |
The range Map of this operator. More... | |
virtual void | apply (const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const |
Apply the preconditioner to X, putting the result in Y. More... | |
Implementation of Ifpack2::Details::NestedPreconditioner | |
virtual void | setInnerPreconditioner (const Teuchos::RCP< Preconditioner< scalar_type, local_ordinal_type, global_ordinal_type, node_type > > &innerPrec) |
Set the inner preconditioner. More... | |
Implementation of Ifpack2::Details::CanChangeMatrix | |
virtual void | setMatrix (const Teuchos::RCP< const row_matrix_type > &A) |
Change the matrix to be preconditioned. More... | |
Implementation of Teuchos::Describable | |
std::string | description () const |
Return a simple one-line description of this object. More... | |
void | describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const |
Print the object with some verbosity level to an FancyOStream object. More... | |
Public Member Functions inherited from Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > | |
virtual | ~Preconditioner () |
Destructor. More... | |
virtual void | apply (const Tpetra::MultiVector< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > &X, Tpetra::MultiVector< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, MatrixType::scalar_typealpha=Teuchos::ScalarTraits< MatrixType::scalar_type >::one(), MatrixType::scalar_typebeta=Teuchos::ScalarTraits< MatrixType::scalar_type >::zero()) const =0 |
Apply the preconditioner to X, putting the result in Y. More... | |
Public Member Functions inherited from Ifpack2::Details::CanChangeMatrix< Tpetra::RowMatrix< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > | |
virtual void | setMatrix (const Teuchos::RCP< const Tpetra::RowMatrix< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > &A)=0 |
Set the new matrix. More... | |
virtual | ~CanChangeMatrix () |
Destructor. More... | |
Public Member Functions inherited from Ifpack2::Details::NestedPreconditioner< Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > | |
virtual void | setInnerPreconditioner (const Teuchos::RCP< Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > &innerPrec)=0 |
Set the inner preconditioner. More... | |
Additive Schwarz domain decomposition for Tpetra sparse matrices.
MatrixType | A specialization of Tpetra::RowMatrix |
LocalInverseType | The type of the solver for the local (subdomain) problem. DO NOT USE ANY TYPE HERE OTHER THAN Ifpack2::Preconditioner. The default is perfectly fine. This template parameter only exists for backwards compatibility. |
This class implements Additive Schwarz domain decomposition, with optional overlap. It operates on a given Tpetra::RowMatrix. Each subdomain corresponds to exactly one MPI process in the given matrix's MPI communicator. For nonoverlapping domain decomposition within an MPI process, please use BlockRelaxation.
This class is a subclass of Tpetra::Operator, like all other subclasses of Preconditioner. Thus, the apply() method applies the preconditioner to a multivector.
One-level overlapping domain decomposition preconditioners use local solvers of Dirichlet type. This means that the solver effectively applies the inverse of the local matrix (possibly with overlap) to the residual to be preconditioned.
The preconditioner can be written as:
\[ P_{AS}^{-1} = \sum_{i=1}^M P_i A_i^{-1} R_i, \]
where \(M\) is the number of subdomains (in this case, the number of (MPI) processes in the computation), \(R_i\) is an operator that restricts the global vector to the vector lying on subdomain \(i\), \(P_i\) is the prolongator operator, and
\[ A_i = R_i A P_i. \]
Constructing a Schwarz preconditioner takes two steps:
The definition of the restriction and prolongation operators \(R_i\) and \(R_i^T\) depends on the level of overlap. If the overlap level is zero (no overlap), their implementation is trivial; \(R_i\) will return all the local components. For nonzero overlap, Tpetra's data redistribution facilities (Tpetra::Import) will be used to bring in the required data. Users may control how these data are combined with existing data, by setting the combine mode parameter. Valid combine modes include "ADD", "INSERT", "REPLACE", "ABSMAX", and "ZERO". These correspond to the valid values of Tpetra::CombineMode.
To solve linear systems involving \(A_i\) on each subdomain, the user can adopt any subclass of Preconditioner. Users must set this at run time by specifying the inner solver in the ParameterList, or by setting the inner solver instance explicitly.
The local matrix \(A_i\) can be filtered, to eliminate singletons, and reordered. At the present time, the only available reordering algorithm is RCM (reverse Cuthill-Mckee). Other orderings will be supported by the Zoltan2 package in the future.
The default is Restricted Additive Schwarz (RAS), which uses CombineMode Zero, see discussion below. Note that RAS does not preserve symmetry, so is generally not suitable as a preconditioner for CG. Classical Additive Schwarz is supported by setting the CombineMode to Add.
AdditiveSchwarz supports different combine modes. These are rules for combining overlapping entries from the subdomain solves' results, to form the final additive Schwarz result. The two modes that users will likely find most useful are "Add" and "Zero". "Add" sums overlapping results, and "Zero" ignores them.
The Zero combine mode can be a bit confusing. It helps to compare the Add combine mode without overlap, to the Zero combine mode with overlap. Consider the following 4 x 4 linear system \(Ax = b\), where
\[ A = \left( \begin{array}{rrrr} 2 & -1 & 0 & 0 \\ -1 & 2 & -1 & 0 \\ 0 & -1 & 2 & -1 \\ 0 & 0 & -1 & 2 \end{array} \right) \]
and
\[ b = \left( \begin{array}{rrrr} 1 \\ 4 \\ 9 \\ 16 \end{array} \right) . \]
Suppose that we give the first two rows of A and b to Process 0, and the last two rows of A and b to Process 1.
If we use additive Schwarz without overlap, and use the Add combine mode, then each process must solve a linear system with the following 2 x 2 matrix:
\[ A_{local} = \left( \begin{array}{rr} 2 & -1 \\ -1 & 2 \end{array} \right). \]
The inverse of this matrix is
\[ A_{local}^{-1} = \frac{1}{3} \cdot \left( \begin{array}{rr} 2 & 1 \\ 1 & 2 \end{array} \right). \]
Process 0 gets the right-hand side \((1, 4)^T\) and thus the local solution \((2, 3)^T\). Process 1 gets the right-hand side \((9, 16)^T\) and thus the local solution \((34/3, 41/3)^T\). The Add combine mode sums entries corresponding to overlap, but in this case there is no overlap, so we just concatenate the local results. Thus, the result of applying additive Schwarz once is
\[ x_{Add} = \begin{array}{r} 2 \\ 3 \\ 34/3 \\ 41/3 \end{array}. \]
If we introduce one level of overlap on each of these processes, and use the Zero combine mode with additive Schwarz, then each process has to solve a linear system with the following 3 x 3 matrix:
\[ A_{local} = \left( \begin{array}{rrr} 2 & -1 & 0 \\ -1 & 2 & -1 \\ 0 & -1 & 2 \end{array} \right). \]
The inverse of this matrix is
\[ A_{local}^{-1} = \frac{1}{4} \cdot \left( \begin{array}{rrr} 3 & 2 & 1 \\ 2 & 4 & 2 \\ 1 & 2 & 3 \end{array} \right). \]
Process 0 gets the right-hand side \((1, 4, 9)^T\) and thus the local solution \((5, 9, 9)^T\). Process 1 gets the right-hand side \((4, 9, 16)^T\) and thus the local solution \((23/2, 19, 35/2)^T\). The Zero combine mode ignores "remote entries," so that the result of applying additive Schwarz once is
\[ x_{Zero} = \begin{array}{r} 5 \\ 9 \\ 19 \\ 35/2 \end{array}. \]
Even though both of the above examples combine the local results in the same way, Zero produces a different final result than Add, because the Zero example uses overlap 1, whereas the Add example uses overlap 0.
This class gives you two ways to specify the subdomain solver. First, you may set the "subdomain solver name" or "inner preconditioner name" parameters in the input to setParameters() or setParameterList(). Second, you may construct the subdomain solver yourself, as an Ifpack2::Preconditioner instance, and give it to setInnerPreconditioner().
Please refer to the documentation of setParameters for a complete discussion of subdomain solvers and their parameters.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::scalar_type = typename MatrixType::scalar_type |
The type of the entries of the input MatrixType.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::local_ordinal_type = typename MatrixType::local_ordinal_type |
The type of local indices in the input MatrixType.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::global_ordinal_type = typename MatrixType::global_ordinal_type |
The type of global indices in the input MatrixType.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::node_type = typename MatrixType::node_type |
The Node type used by the input MatrixType.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::magnitude_type = typename Teuchos::ScalarTraits<scalar_type>::magnitudeType |
The type of the magnitude (absolute value) of a matrix entry.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::row_matrix_type = Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> |
The Tpetra::RowMatrix specialization matching MatrixType.
using Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::crs_matrix_type = Tpetra::CrsMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> |
The Tpetra::CrsMatrix specialization that is a subclass of MatrixType.
Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::AdditiveSchwarz | ( | const Teuchos::RCP< const row_matrix_type > & | A | ) |
Constructor that takes a matrix.
A | [in] The matrix to be preconditioned. |
Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::AdditiveSchwarz | ( | const Teuchos::RCP< const row_matrix_type > & | A, |
const int | overlapLevel | ||
) |
Constructor that takes a matrix and the level of overlap.
A | [in] The matrix to be preconditioned. |
overlapLevel | [in] The level of overlap. Must be nonnegative. Zero means no overlap. |
|
virtualdefault |
Destructor.
|
virtual |
The domain Map of this operator.
|
virtual |
The range Map of this operator.
|
virtual |
Apply the preconditioner to X, putting the result in Y.
|
virtual |
Set the inner preconditioner.
Most users do not need to call this method. Factory will handle calling this method if necessary. One case in which users might need to call this method, is if they implement their own subdomain solver as a subclass of Preconditioner, and want to use that subdomain solver in this class.
innerPrec | [in/out] The inner preconditioner. Its matrix (if it has one) may be replaced by a matrix specified by the outer (this) preconditioner. |
*this
) to itself as an inner solver, or otherwise set up a cyclic dependency between preconditioner instances. You MAY use an inner solver of the same TYPE as this object (as long as this makes sense mathematically), but it must be a different instance of that type.It does not make sense to nest different instances of AdditiveSchwarz, since it currently only sets up one subdomain per MPI process in the input matrix's communicator. Thus, if you were to use another AdditiveSchwarz instance as the inner preconditioner for AdditiveSchwarz, it would not do anything, since the inner preconditioner's input matrix would only have one process in its communicator. (AdditiveSchwarz does nothing in that case.)
&*innerPrec != this
.This method has collective semantics, because it may call initialize() or compute() on innerPrec
, in order to synchronize the inner preconditioner's state with that of the AdditiveSchwarz instance.
|
virtual |
Change the matrix to be preconditioned.
[in] | A | The new matrix. |
! isInitialized ()
! isComputed ()
Calling this method resets the preconditioner's state. After calling this method with a nonnull input, you must first call initialize() and compute() (in that order) before you may call apply().
You may call this method with a null input. If A is null, then you may not call initialize() or compute() until you first call this method again with a nonnull input. This method invalidates any previous factorization whether or not A is null, so calling setMatrix() with a null input is one way to clear the preconditioner's state (and free any memory that it may be using).
The new matrix A need not necessarily have the same Maps or even the same communicator as the original matrix.
|
virtual |
|
virtual |
Set the preconditioner's parameters.
plist | [in] List of parameters. |
This version of the method takes a const list, as required by the Preconditioner interface. setParameterList() takes a nonconst pointer to a list, in order to match the Teuchos::ParameterListAcceptor interface.
Both this method and setParameterList() have "delta" behavior. That is, if called twice with two different lists, any unspecified parameters in the new list retain their values in the old list.
In many cases, calling this method may require calling initialize() and compute() to recompute the preconditioner.
Accepted parameters include the following:
std::string
): the name of the subdomain solver. See discussion below.std::string
): The rule for combining incoming data with existing data in overlap regions. Valid values include "ADD", "INSERT", "REPLACE", "ABSMAX", and "ZERO". These correspond to the valid values of Tpetra::CombineMode. The default combine mode is "ZERO", meaning that overlapping incoming entries from different processes are not combined. (See the class documentation for an explanation and example of the "ZERO" vs. "ADD" combine modes.)int
): The level of overlap. The default is zero, meaning no overlap.bool
): Whether to use Zoltan2 to do reordering. If true, then Trilinos must have been built with Zoltan2 and Xpetra enabled. Default is false.int
): This option does not currently work.bool
): If true, exclude rows with just a single entry on the calling process. Default is false.int
): Numter of iterations to perform. Default is 1.This class lets users specify any subdomain solver they want, by calling setInnerPreconditioner(). However, users may instead specify the subdomain solver by setting the "inner preconditioner name" parameter (or any of its aliases). If they choose to do so, they may use any Ifpack2 preconditioner. These include but are not necessarily limited to the following:
This name need not necessarily correspond with LocalInverseType
. If the user does not specify this parameter, the following procedure specifies the default:
LocalInverseType
is just Preconditioner, then this class uses a default, which is currently "ILUT". LocalInverseType
is a concrete Preconditioner subclass, and if that subclass is in the above supported list of subdomain solver types, then this class uses that subclass as the subdomain solver. LocalInverseType
is a concrete Preconditioner subclass, and if that subclass is not in the above supported list of subdomain solver types, then users have one of two options, both of which we discuss below. The subdomain solver names "INVALID" and "CUSTOM" are reserved for internal use.
AdditiveSchwarz only knows, on its own, how to create "non-nested" preconditioners as inner preconditioners (i.e., subdomain solvers). It can't create nested preconditioners (e.g., AdditiveSchwarz and SupportGraph) on its own as inner preconditioners, and it doesn't know how to create arbitrary subclasses of Ifpack2::Preconditioner unless Ifpack2::Factory knows how to create them.
This leaves users two options in order to have any preconditioner as AdditiveSchwarz's inner preconditioner:
prec
, then users who don't want to create the inner preconditioner themselves must create AdditiveSchwarz using Factory, not by invoking AdditiveSchwarz's constructor themselves. Factory will set up the inner preconditioner for them before it returns the AdditiveSchwarz instance. prec
(for example, if it is not even implemented in Ifpack2), then users must create the inner preconditioner instance themselves, and give it to AdditiveSchwarz using setInnerPreconditioner. In this case, AdditiveSchwarz's ParameterList must not specify the inner preconditioner's name. Users are responsible for ensuring that the parameters they provide to setParameters() are up to date. For example, if the users first set an inner preconditioner using setInnerPreconditioner(), and then call setParameters() with the "inner preconditioner name" parameter set, AdditiveSchwarz will get rid of the users' inner preconditioner and attempt to create a new inner preconditioner itself. Remember that AdditiveSchwarz's ParameterList has "delta" (relative) semantics! If you don't specify a parameter, the current state is not changed.
If you specify a sublist of parameters to give to the subdomain solver, setInnerPreconditioner() does not pass that sublist to its argument. This is because we presume that if you call setInnerPreconditioner(), the input subdomain solver probably has a type that AdditiveSchwarz does not know how to create by itself, so the existing parameter list cannot apply.
On the other hand, if, after calling setInnerPreconditioner(), you then call setParameters(), we do pass any provided sublist of subdomain solver parameters to the inner solver. If no such sublist was provided, we do not call setParameters() on the inner solver.
The reason the last sentence matters is because not every inner solver necessarily has "delta" semantics for setParameters(). "Delta" or "relative" semantics means that an empty input ParameterList doesn't change any existing parameters. "Non-delta" or "absolute" semantics means that an empty input ParameterList causes all parameters to be set to their default values. (The difference matters if the user has called setParameters() before on the subdomain solver, with nondefault values.) If the user didn't specify a sublist for the inner solver, we assume that the user doesn't want to change the inner solver's parameters.
void Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::setParameterList | ( | const Teuchos::RCP< Teuchos::ParameterList > & | plist | ) |
Set the preconditioner's parameters.
This version of the method takes a nonconst pointer to a list, in order to match the Teuchos::ParameterListAcceptor interface. setParameters() takes a const list, as required by the Preconditioner interface.
Both this method and setParameterList() have "delta" behavior. That is, if called twice with two different lists, any unspecified parameters in the new list retain their values in the old list.
In many cases, calling this method may require calling initialize() and compute() to recompute the preconditioner.
plist | [in/out] On input: List of parameters, or Teuchos::null (meaning "do not change the current parameter values"). On output: If nonnull, any missing parameters are filled in with their current values (or their default values, if they have not yet been set). |
See the documentation of setParameters() for a list of the parameters this method accepts, and their default values.
|
inlinevirtual |
Set this preconditioner's parameters.
Reimplemented from Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type >.
Teuchos::RCP< const Teuchos::ParameterList > Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::getValidParameters | ( | ) | const |
Get a list of the preconditioner's default parameters.
See the documentation of setParameters() for a list of the parameters that AdditiveSchwarz accepts.
|
virtual |
Computes all (graph-related) data necessary to initialize the preconditioner.
|
virtual |
Returns true if the preconditioner has been successfully initialized, false otherwise.
|
virtual |
Computes all (coefficient) data necessary to apply the preconditioner.
|
virtual |
Returns true if the preconditioner has been successfully computed, false otherwise.
|
virtual |
Returns the number of calls to initialize().
|
virtual |
Returns the number of calls to compute().
|
virtual |
Returns the number of calls to apply().
|
virtual |
Returns the time spent in initialize().
|
virtual |
Returns the time spent in compute().
|
virtual |
Returns the time spent in apply().
std::string Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::description | ( | ) | const |
Return a simple one-line description of this object.
void Ifpack2::AdditiveSchwarz< MatrixType, LocalInverseType >::describe | ( | Teuchos::FancyOStream & | out, |
const Teuchos::EVerbosityLevel | verbLevel = Teuchos::Describable::verbLevel_default |
||
) | const |
Print the object with some verbosity level to an FancyOStream object.
|
virtual |
Prints basic information on iostream. This function is used by operator<<.
|
virtual |
Returns the level of overlap.