|
Ifpack2 Templated Preconditioning Package
Version 1.0
|
|
Ifpack2 Templated Preconditioning Package
Version 1.0
|
A class for constructing and using an incomplete lower/upper (ILU) factorization of a given Tpetra::RowMatrix. More...
#include <Ifpack2_RILUK_decl.hpp>
Public Member Functions | |
| RILUK (const Teuchos::RCP< const MatrixType > &A_in) | |
| RILUK constuctor with variable number of indices per row. | |
| virtual | ~RILUK () |
| Ifpack2_RILUK Destructor. | |
| void | SetRelaxValue (magnitudeType RelaxValue) |
| Set RILU(k) relaxation parameter. | |
| void | SetAbsoluteThreshold (magnitudeType Athresh) |
| Set absolute threshold value. | |
| void | SetRelativeThreshold (magnitudeType Rthresh) |
| Set relative threshold value. | |
| void | SetOverlapMode (Tpetra::CombineMode OverlapMode) |
| Set overlap mode type. | |
| void | setParameters (const Teuchos::ParameterList ¶meterlist) |
| Set parameters using a Teuchos::ParameterList object. | |
| void | initialize () |
| Computes all (graph-related) data necessary to initialize the preconditioner. | |
| bool | isInitialized () const |
| Returns true if the preconditioner has been successfully initialized, false otherwise. | |
| int | getNumInitialize () const |
| Returns the number of calls to initialize(). | |
| void | compute () |
| Compute ILU factors L and U using the specified diagonal perturbation thresholds and relaxation parameters. | |
| bool | isComputed () const |
| If compute() is completed, this query returns true, otherwise it returns false. | |
| int | getNumCompute () const |
| Returns the number of calls to compute(). | |
| int | getNumApply () const |
| Returns the number of calls to Apply(). | |
| double | getInitializeTime () const |
| Returns the time spent in Initialize(). | |
| double | getComputeTime () const |
| Returns the time spent in Compute(). | |
| double | getApplyTime () const |
| Returns the time spent in Apply(). | |
| void | apply (const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const |
| Returns the result of a RILUK forward/back solve on a Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> X in Y. | |
| int | Multiply (const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS) const |
| Returns the result of multiplying U, D and L in that order on an Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> X in Y. | |
| magnitudeType | computeCondEst (Teuchos::ETransp mode) const |
| Returns the maximum over all the condition number estimate for each local ILU set of factors. | |
| magnitudeType | computeCondEst (CondestType CT=Ifpack2::Cheap, LocalOrdinal MaxIters=1550, magnitudeType Tol=1e-9, const Teuchos::Ptr< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &Matrix=Teuchos::null) |
| Computes the condition number estimate and returns its value. | |
| magnitudeType | getCondEst () const |
| Returns the computed condition number estimate, or -1.0 if not computed. | |
| Teuchos::RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > | getMatrix () const |
| Returns a pointer to the input matrix. | |
| magnitudeType | GetRelaxValue () const |
| Get RILU(k) relaxation parameter. | |
| magnitudeType | getAbsoluteThreshold () const |
| Get absolute threshold value. | |
| magnitudeType | getRelativeThreshold () const |
| Get relative threshold value. | |
| Tpetra::CombineMode | getOverlapMode () |
| Get overlap mode type. | |
| int | getGlobalNumEntries () const |
| Returns the number of nonzero entries in the global graph. | |
| const Teuchos::RCP < Ifpack2::IlukGraph < LocalOrdinal, GlobalOrdinal, Node > > & | getGraph () const |
| Returns the Ifpack2::IlukGraph associated with this factored matrix. | |
| const MatrixType & | getL () const |
| Returns the L factor associated with this factored matrix. | |
| const Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & | getD () const |
| Returns the D factor associated with this factored matrix. | |
| const MatrixType & | getU () const |
| Returns the U factor associated with this factored matrix. | |
| const Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > & | getDomainMap () const |
| Returns the Tpetra::Map object associated with the domain of this operator. | |
| const Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > & | getRangeMap () const |
| Returns the Tpetra::Map object associated with the range of this operator. | |
A class for constructing and using an incomplete lower/upper (ILU) factorization of a given Tpetra::RowMatrix.
Ifpack2::RILUK computes a "Relaxed" ILU factorization with level k fill of a given Tpetra::RowMatrix.
For a complete list of valid parameters, see Ifpack2::RILUK::setParameters.
The factorization that is produced is a function of several parameters:
The pattern of the matrix - All fill is derived from the original matrix nonzero structure. Level zero fill is defined as the original matrix pattern (nonzero structure), even if the matrix value at an entry is stored as a zero. (Thus it is possible to add entries to the ILU factors by adding zero entries the original matrix.)
Level of fill - Starting with the original matrix pattern as level fill of zero, the next level of fill is determined by analyzing the graph of the previous level and determining nonzero fill that is a result of combining entries that were from previous level only (not the current level). This rule limits fill to entries that are direct decendents from the previous level graph. Fill for level k is determined by applying this rule recursively. For sufficiently large values of k, the fill would eventually be complete and an exact LU factorization would be computed.
Level of overlap - All Ifpack2 preconditioners work on parallel distributed memory computers by using the row partitioning the user input matrix to determine the partitioning for local ILU factors. If the level of overlap is set to zero, the rows of the user matrix that are stored on a given processor are treated as a self-contained local matrix and all column entries that reach to off-processor entries are ignored. Setting the level of overlap to one tells Ifpack to increase the size of the local matrix by adding rows that are reached to by rows owned by this processor. Increasing levels of overlap are defined recursively in the same way. For sufficiently large levels of overlap, the entire matrix would be part of each processor's local ILU factorization process. Level of overlap is defined during the construction of the Ifpack2_IlukGraph object.
Once the factorization is computed, applying the factorization \(LUy = x\) results in redundant approximations for any elements of y that correspond to rows that are part of more than one local ILU factor. The OverlapMode (changed by calling SetOverlapMode()) defines how these redundancies are handled using the Tpetra::CombineMode enum. The default is to zero out all values of y for rows that were not part of the original matrix row distribution.
Fraction of relaxation - Ifpack2_RILUK computes the ILU factorization row-by-row. As entries at a given row are computed, some number of them will be dropped because they do match the prescribed sparsity pattern. The relaxation factor determines how these dropped values will be handled. If the RelaxValue (changed by calling setRelaxValue()) is zero, then these extra entries will by dropped. This is a classical ILU approach. If the RelaxValue is 1, then the sum of the extra entries will be added to the diagonal. This is a classical Modified ILU (MILU) approach. If RelaxValue is between 0 and 1, then RelaxValue times the sum of extra entries will be added to the diagonal.
For most situations, RelaxValue should be set to zero. For certain kinds of problems, e.g., reservoir modeling, there is a conservation principle involved such that any operator should obey a zero row-sum property. MILU was designed for these cases and you should set the RelaxValue to 1. For other situations, setting RelaxValue to some nonzero value may improve the stability of factorization, and can be used if the computed ILU factors are poorly conditioned.
Estimating Preconditioner Condition Numbers
For ill-conditioned matrices, we often have difficulty computing usable incomplete factorizations. The most common source of problems is that the factorization may encounter a small or zero pivot, in which case the factorization can fail, or even if the factorization succeeds, the factors may be so poorly conditioned that use of them in the iterative phase produces meaningless results. Before we can fix this problem, we must be able to detect it. To this end, we use a simple but effective condition number estimate for \((LU)^{-1}\).
The condition of a matrix \(B\), called \(cond_p(B)\), is defined as \(cond_p(B) = \|B\|_p\|B^{-1}\|_p\) in some appropriate norm \(p\). \(cond_p(B)\) gives some indication of how many accurate floating point digits can be expected from operations involving the matrix and its inverse. A condition number approaching the accuracy of a given floating point number system, about 15 decimal digits in IEEE double precision, means that any results involving \(B\) or \(B^{-1}\) may be meaningless.
The \(\infty\)-norm of a vector \(y\) is defined as the maximum of the absolute values of the vector entries, and the \(\infty\)-norm of a matrix C is defined as \(\|C\|_\infty = \max_{\|y\|_\infty = 1} \|Cy\|_\infty\). A crude lower bound for the \(cond_\infty(C)\) is \(\|C^{-1}e\|_\infty\) where \(e = (1, 1, \ldots, 1)^T\). It is a lower bound because \(cond_\infty(C) = \|C\|_\infty\|C^{-1}\|_\infty \ge \|C^{-1}\|_\infty \ge |C^{-1}e\|_\infty\).
For our purposes, we want to estimate \(cond_\infty(LU)\), where \(L\) and \(U\) are our incomplete factors. Edmond in his Ph.D. thesis demonstrates that \(\|(LU)^{-1}e\|_\infty\) provides an effective estimate for \(cond_\infty(LU)\). Furthermore, since finding \(z\) such that \(LUz = y\) is a basic kernel for applying the preconditioner, computing this estimate of \(cond_\infty(LU)\) is performed by setting \(y = e\), calling the solve kernel to compute \(z\) and then computing \(\|z\|_\infty\).
A priori Diagonal Perturbations
Given the above method to estimate the conditioning of the incomplete factors, if we detect that our factorization is too ill-conditioned we can improve the conditioning by perturbing the matrix diagonal and restarting the factorization using this more diagonally dominant matrix. In order to apply perturbation, prior to starting the factorization, we compute a diagonal perturbation of our matrix \(A\) and perform the factorization on this perturbed matrix. The overhead cost of perturbing the diagonal is minimal since the first step in computing the incomplete factors is to copy the matrix \(A\) into the memory space for the incomplete factors. We simply compute the perturbed diagonal at this point.
The actual perturbation values we use are the diagonal values \((d_1, d_2, \ldots, d_n)\) with \(d_i = sgn(d_i)\alpha + d_i\rho\), \(i=1, 2, \ldots, n\), where \(n\) is the matrix dimension and \(sgn(d_i)\) returns the sign of the diagonal entry. This has the effect of forcing the diagonal values to have minimal magnitude of \(\alpha\) and to increase each by an amount proportional to \(\rho\), and still keep the sign of the original diagonal entry.
Counting Floating Point Operations
Each Ifpack2::RILUK object keeps track of the number of serial floating point operations performed using the specified object as the this argument to the function. The Flops() function returns this number as a double precision number. Using this information, in conjunction with the Teuchos::Time class, one can get accurate parallel performance numbers. The ResetFlops() function resets the floating point counter.
| Ifpack2::RILUK< MatrixType >::RILUK | ( | const Teuchos::RCP< const MatrixType > & | A_in | ) |
| Ifpack2::RILUK< MatrixType >::~RILUK | ( | ) | [virtual] |
Ifpack2_RILUK Destructor.
| void Ifpack2::RILUK< MatrixType >::SetRelaxValue | ( | magnitudeType | RelaxValue | ) | [inline] |
Set RILU(k) relaxation parameter.
| void Ifpack2::RILUK< MatrixType >::SetAbsoluteThreshold | ( | magnitudeType | Athresh | ) | [inline] |
Set absolute threshold value.
| void Ifpack2::RILUK< MatrixType >::SetRelativeThreshold | ( | magnitudeType | Rthresh | ) | [inline] |
Set relative threshold value.
| void Ifpack2::RILUK< MatrixType >::SetOverlapMode | ( | Tpetra::CombineMode | OverlapMode | ) | [inline] |
Set overlap mode type.
| void Ifpack2::RILUK< MatrixType >::setParameters | ( | const Teuchos::ParameterList & | parameterlist | ) | [virtual] |
Set parameters using a Teuchos::ParameterList object.
| void Ifpack2::RILUK< MatrixType >::initialize | ( | ) | [virtual] |
Computes all (graph-related) data necessary to initialize the preconditioner.
| bool Ifpack2::RILUK< MatrixType >::isInitialized | ( | ) | const [inline, virtual] |
Returns true if the preconditioner has been successfully initialized, false otherwise.
| int Ifpack2::RILUK< MatrixType >::getNumInitialize | ( | ) | const [inline, virtual] |
Returns the number of calls to initialize().
| void Ifpack2::RILUK< MatrixType >::compute | ( | ) | [virtual] |
Compute ILU factors L and U using the specified diagonal perturbation thresholds and relaxation parameters.
This function computes the RILU(k) factors L and U using the current:
initialize() must be called before the factorization can proceed.
| bool Ifpack2::RILUK< MatrixType >::isComputed | ( | ) | const [inline, virtual] |
If compute() is completed, this query returns true, otherwise it returns false.
| int Ifpack2::RILUK< MatrixType >::getNumCompute | ( | ) | const [inline, virtual] |
Returns the number of calls to compute().
| int Ifpack2::RILUK< MatrixType >::getNumApply | ( | ) | const [inline, virtual] |
Returns the number of calls to Apply().
| double Ifpack2::RILUK< MatrixType >::getInitializeTime | ( | ) | const [inline, virtual] |
Returns the time spent in Initialize().
| double Ifpack2::RILUK< MatrixType >::getComputeTime | ( | ) | const [inline, virtual] |
Returns the time spent in Compute().
| double Ifpack2::RILUK< MatrixType >::getApplyTime | ( | ) | const [inline, virtual] |
Returns the time spent in Apply().
| void Ifpack2::RILUK< MatrixType >::apply | ( | const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & | X, |
| Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & | Y, | ||
| Teuchos::ETransp | mode = Teuchos::NO_TRANS, |
||
| Scalar | alpha = Teuchos::ScalarTraits<Scalar>::one(), |
||
| Scalar | beta = Teuchos::ScalarTraits<Scalar>::zero() |
||
| ) | const [virtual] |
Returns the result of a RILUK forward/back solve on a Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> X in Y.
| In | Trans -If true, solve transpose problem. |
| In | X - A Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> of dimension NumVectors to solve for. |
| Out | Y -A Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> of dimension NumVectorscontaining result. |
| int Ifpack2::RILUK< MatrixType >::Multiply | ( | const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & | X, |
| Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & | Y, | ||
| Teuchos::ETransp | mode = Teuchos::NO_TRANS |
||
| ) | const |
Returns the result of multiplying U, D and L in that order on an Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> X in Y.
| In | Trans -If true, multiply by L^T, D and U^T in that order. |
| In | X - A Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> of dimension NumVectors to solve for. |
| Out | Y -A Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> of dimension NumVectorscontaining result. |
| Teuchos::ScalarTraits< typename MatrixType::scalar_type >::magnitudeType Ifpack2::RILUK< MatrixType >::computeCondEst | ( | Teuchos::ETransp | mode | ) | const |
Returns the maximum over all the condition number estimate for each local ILU set of factors.
This functions computes a local condition number estimate on each processor and return the maximum over all processors of the estimate.
| In | Trans -If true, solve transpose problem. |
| Out | ConditionNumberEstimate - The maximum across all processors of the infinity-norm estimate of the condition number of the inverse of LDU. |
| magnitudeType Ifpack2::RILUK< MatrixType >::computeCondEst | ( | CondestType | CT = Ifpack2::Cheap, |
| LocalOrdinal | MaxIters = 1550, |
||
| magnitudeType | Tol = 1e-9, |
||
| const Teuchos::Ptr< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > & | Matrix = Teuchos::null |
||
| ) | [inline, virtual] |
Computes the condition number estimate and returns its value.
| magnitudeType Ifpack2::RILUK< MatrixType >::getCondEst | ( | ) | const [inline, virtual] |
Returns the computed condition number estimate, or -1.0 if not computed.
| Teuchos::RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Ifpack2::RILUK< MatrixType >::getMatrix | ( | ) | const [inline, virtual] |
Returns a pointer to the input matrix.
| magnitudeType Ifpack2::RILUK< MatrixType >::GetRelaxValue | ( | ) | const [inline] |
Get RILU(k) relaxation parameter.
| magnitudeType Ifpack2::RILUK< MatrixType >::getAbsoluteThreshold | ( | ) | const [inline] |
Get absolute threshold value.
| magnitudeType Ifpack2::RILUK< MatrixType >::getRelativeThreshold | ( | ) | const [inline] |
Get relative threshold value.
| Tpetra::CombineMode Ifpack2::RILUK< MatrixType >::getOverlapMode | ( | ) | [inline] |
Get overlap mode type.
| int Ifpack2::RILUK< MatrixType >::getGlobalNumEntries | ( | ) | const [inline] |
Returns the number of nonzero entries in the global graph.
| const Teuchos::RCP<Ifpack2::IlukGraph<LocalOrdinal,GlobalOrdinal,Node> >& Ifpack2::RILUK< MatrixType >::getGraph | ( | ) | const [inline] |
Returns the Ifpack2::IlukGraph associated with this factored matrix.
| const MatrixType& Ifpack2::RILUK< MatrixType >::getL | ( | ) | const [inline] |
Returns the L factor associated with this factored matrix.
| const Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Ifpack2::RILUK< MatrixType >::getD | ( | ) | const [inline] |
Returns the D factor associated with this factored matrix.
| const MatrixType& Ifpack2::RILUK< MatrixType >::getU | ( | ) | const [inline] |
Returns the U factor associated with this factored matrix.
| const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >& Ifpack2::RILUK< MatrixType >::getDomainMap | ( | ) | const [inline, virtual] |
Returns the Tpetra::Map object associated with the domain of this operator.
| const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >& Ifpack2::RILUK< MatrixType >::getRangeMap | ( | ) | const [inline, virtual] |
Returns the Tpetra::Map object associated with the range of this operator.
1.7.6.1
