Trilinos/IFPACK: Object-oriented Algebraic Preconditioner Package

Table of Contents

• Introduction
• Download IFPACK
• Configuring IFPACK
• Examples of Usage
• List of Supported Parameters
• Brief Description IFPACK Classes
• Error Table
• Known Bugs and Future Developments
• Copyright
• Browse all of Ifpack as a single doxygen collection

Introduction

IFPACK, developed by Marzio Sala (ETHZ/D-INFK) and Micheal Heroux (SNL 9214), provides a suite of object-oriented algebraic preconditioners for the solution of preconditioned iterative solvers. IFPACK constructors expect an Epetra_RowMatrix object for construction. IFPACK is part of the Trilinos Solver Project and IFPACK object interacts well with other Trilinos classes. In particular, IFPACK can be used as a preconditioner for AztecOO.

Most IFPACK preconditioners can be rewritten as additive Schwarz methods, of overlapping domain decomposition type. The user can adopt a minimal-overlap (that is, zero-row overlap), or ask IFPACK to extend the overlap. The resulting preconditioner reads:

\[ P_{IFPACK}^{-1} = \sum_{i = 0}^{NumProcs-1} P_i A_i^{-1} R_i \]

where \(R_i\) is the restriction operator from the global vector, to the overlapping subdomain \(i\), and \(P_i\) is the prolongator operator. \(P_i\) is generally the transpose of \(R_i\) (in which case the resulting preconditioner is symmetic). It is assumed that each subdomain is assigned to a different processor.

A key component of the previous formula is a strategy to apply \(A_i^{-1}\). Using IFPACK, this can be defined by one of the following:

• one or more steps of a point/block method (like Jacobi and Gauss-Seidel);
• an incomplete factorization (for symmetric and non-symmetric matrices, with dropping based on graph or on values);
• an exact LU solve on each subdomain (using Amesos' factorizations).

More precisely, the choices for \(A_i^{-1}\) are:

• Point relaxation preconditioners (through class Ifpack_PointRelaxation). Point relaxation preconditioners are probably the simplest iterative methods, and are generally used as smoothers in multilevel methods (for instance, within ML). Available choices are:

  • Jacobi
  • Gauss-Seidel
  • symmetric Guass-Seidel

• Block relaxation preconditioners (through class Ifpack_BlockRelaxation). Block relaxation preconditioners represent an extension of point preconditioners. The local matrix \(A_i\) is divided into blocks, then a relaxation method is applied on the block structure of \(A_i\). Available choices are:

  • block Jacobi
  • block Gauss-Seidel
  • block symmetric Gauss-Seidel

  The code extracts the diagonal blocks (the application of whose inverses is required), and stored them either as dense or as sparse matrices. Dense matrices should be used if the dimension of each block is small, since LAPACK routines are used to factorize the block and to perform the backsolve. For large blocks, instead, it is suggested to store the blocks are sparse matrices (Epetra_CrsMatrix's). Any IFPACK preconditioner derived from Ifpack_Preconditioner can be used to apply the inverse of the diagonal blocks.

• Point incomplete factorizations (through classes Ifpack_IC, Ifpack_ICT, Ifpack_ILU, Ifpack_ILUT). See General description of incomplete factorizations .

• Exact factorizations (through class Ifpack_Amesos):

  • all Amesos classes can be used to compute the LU factorization of a given matrix.

• Chebyshev polynomials (through class Ifpack_Chebyshev).

Note

IFPACK preconditioners can be used as smoothers for multilevel solvers, like ML.

Download IFPACK

IFPACK can be downloaded from the web page http://trilinos.sandia.gov/download/

Configuring IFPACK

IFPACK is configured with autotools. IFPACK is enabled by default if configured at the Trilinos level (please refer to the Trilinos documentation for more details). This section briefly recalls some of the configure parameters that affect the IFPACK compilation. For a complete list of parameters, please type

% $TRILINOS_HOME/packages/ifpack/configure --help

An (incomplete) list of parameters recognized by the IFPACK configure is reported below.

• "enable-amesos" enables the support for Amesos (off by default).
• "enable-aztecoo" enables the support for AztecOO (on by default).
• "enable-teuchos" enables the support for Teuchos (off by default).
• "enable-thyra" enables the support for Thyra (off by default).
• "enable-triutils" enables the support for Triutils (on by default).
• "enable-ifpack-metis" enables the support for the METIS package. The location of the header files should be specified using "--with-incdirs=-I/include/location", the path of the METIS library with "--with-ldflags=-L/lib/location", and the library name using "--with-libs=-lmetis-name.".

Note

IFPACK cannot be build without Epetra. Most of IFPACK preconditioners cannot be used without the Teuchos library. Most of the examples requires triutils and AztecOO.

Examples of Usage

This section details how to use the IFPACK Factory class.

Probably, the easiest way to use IFPACK is through the function class Ifpack. This is a factory class, that contains only one method, Create(). A call to Create() returns a pointer to a newly-created IFPACK object. The user is responsible of deleting this object.

An example of usage is reported by the following code:

#include "Ifpack.h"
...
Ifpack Factory;
Epetra_RowMatrix* A; // A is already FillComplete()'d
                     // but its values can still change (not its structure)
string PrecType = "Amesos"; // exact solve on each subdomain 
int OverlapLevel = 1; // one row of overlap among the processes
Ifpack_Preconditioner* Prec = Factory.Create(PrecType, A, OverlapLevel);
assert (Prec != 0);

IFPACK_CHK_ERR(Prec->SetParameters(List));
IFPACK_CHK_ERR(Prec->Initialize());
IFPACK_CHK_ERR(Prec->Compute());

// we can query the status of the preconditioner
assert (Prec->IsInitialized() == true);
assert (Prec->IsComputed() == true);

// information about the preconditioner
// 1.- number of initialization phases
cout << Prec->NumInitialize() << endl; 
// 2.- number of computation phases
cout << Prec->NumCompute() << endl; 
// 3.- number of applications of the prconditioner
cout << Prec->NumApplyInverse() << endl; 
// we can query for flops and CPU time as well, consult the
// documentation of Ifpack_Preconditioner, or simply write
cout << *Prec;

// now destroy the preconditioner
delete Prec;

Although the Ifpack factory is appropriate for most users, it is possible to create IFPACK preconditioners by handling directly an Ifpack_AdditiveSchwarz object. The preconditioner of the previous example is equivalent to the one created in the following fragment of code:

#include "Ifpack_AdditiveSchwarz.h"
...
Epetra_RowMatrix* A; 
int OverlapLevel = 1;
Ifpack_AdditiveSchwarz Prec<Ifpack_Amesos> (A, OverlapLevel);

Another example is the following, where the local solve is a block Jacobi method, which uses dense containers (LAPACK).

#include "Ifpack_AdditiveSchwarz.h"
#include "Ifpack_DenseContainer.h"
#include "Ifpack_BlockJacobi.h"
...
Epetra_RowMatrix* A; 
int OverlapLevel = 1;
Ifpack_AdditiveSchwarz Prec<Ifpack_BlockJacobi<Ifpack_DenseContainer> > (A, OverlapLevel);

Using sparse containers and Amesos the code may look as follows.

#include "Ifpack_AdditiveSchwarz.h"
#include "Ifpack_SparseContainer.h"
#include "Ifpack_BlockJacobi.h"
#include "Ifpack_Amesos.h"
...
Epetra_RowMatrix* A; 
int OverlapLevel = 1;
Ifpack_AdditiveSchwarz Prec<Ifpack_BlockJacobi<Ifpack_SparseContainer<Ifpack_Amesos> > > (A, OverlapLevel);

Other examples are reported in the example subdirectory with several comments:

• Ex_Factory shows how to use Ifpack factory to create a variety of IFPACK preconditioners. This is the simplest way to use IFPACK. The example details all the phases of IFPACK. Difficulty level: beginners.
• Ex_MatrixMarket reads an MatrixMarket matrix (a simple format that contains i,j and a(i,j) of all nonzero elements), redistribute the matrix to all avaialble processes (using a linear distribution), and solve the problem with AztecOO and IFPACK. Matrices downloaded from math.nist.gov/MatrixMarket/ in .mtx format can be used with this example. A small example is reported. Difficulty level: beginners.
• Ex_AmesosPrec shows how to use exact LU solvers on each subdomain, using Amesos. This example is very similar to the previous one. Difficulty level: beginners.
• Ex_vIctPrec shows how to define incomplete Cholesky factorizations (with dropping based on values), and tune parameters. Although an vICT factorization can be constructed using the Factory (as done in file Ifpack_ex_Factory.cpp), this example tells the real story: how to define an Ifpack_AdditiveSchwar which uses vICT on each subdomain. Difficulty level: intermediate.
• Ex_BlockRelaxationPrec shows how to define block preconditioner in the more general case, using Ifpack_AdditiveSchwarz. Difficulty: advanced.

List of Supported Parameters

The complete list of supported parameters is reported below.

• "relaxation: type": specifies the type of point and block relaxation scheme. Valid values:
  • "Jacobi";
  • "Gauss-Seidel";
  • "symmetric Gauss-Seidel";
• "relaxation: sweeps": specifies the number of sweeps in the application of point relaxation schemes.
• "relaxation: damping factor": specifies the damping factor in the application of point relaxation schemes.
• "relaxation: min diagonal value": replace diagonal values below this value with this value (only for point relaxation).
• "relaxation: zero starting solution": if true, uses the values in the vector to be preconditioned as starting solution. Otherwise, use the zero vector as starting solution.
• "relaxation: backward mode": Does Gauss-Seidel in backward order (rather than forward)
• "relaxation: use l1": Use the "L1-Gauss-Seidel" method of Baker, Falgout, Kolev and Yang (SISC 2011). This can be especially advantageous for using Gauss-Seidel and symmetric Gauss-Seidel in parallel.
• "relaxation: l1 eta": Sets the "eta" parameter used by the L1 GS/SGS method.
• "relaxation: number of local smoothing indices": Sets the number of local unknowns over which to smooth. This is used for selective (local) smoothing.
• "relaxation: local smoothing indices": A list of local unknowns over which to smooth. This can be used for selective (local) smoothing or to force the method to smooth unknowns in a particular order (such as upwind smoothing).
• "partitioner: type": specifies the scheme to adopt to partition the graph (for block relaxation methods only). Valid values:
  • "linear" (uses a linear decomposition);
  • "greedy" (uses a simple greedy algorithm);
  • "metis" (calls METIS, this required –enable-ifpack-metis).
  • "equation" (groups together all the nodes belonging to the same equation. The number of equation is specified using "partitioner: local parts", and it is supposed that the equation for each grid node are ordered consecutevely).
• "partitioner: overlap": specifies the overlap among the blocks which can differ from the overlap among the processors. Note that only the Jacobi block relaxation scheme can take advantage of non-zero overlaps.
• "partitioner: local parts": specifies the number of local blocks.
• "partitioner: root node": specifies the root node for greedy algorithm.
• "partitioner: use symmetric graph": if true, METIS will partition the graph of A + A^T. If false, METIS will partition the graph of A. Note that METIS can core dump if the input graph is non-symmetric. Users should set this option to false only when the graph is symmetric. If the non-symmetry is determined by Dirichlet nodes, then the singleton filter should create a symmetric graph. Note also that dropping techniques applied to non-symmetric matrices can result in non-symmetric graph.
• "amesos: solver type": defines the Amesos solver to be used by class Ifpack_Amesos. Valid values:
  • "Amesos_Lapack" (use Amesos interface to LAPACK);
  • "Amesos_Klu" (use Amesos' internal solver KLU);
  • "Amesos_Umfpack" (use UMFPACK interface);
  • "Amesos_Superlu" (use serial SuperLU interface);
  • "Amesos_Mumps" (use MUMPS interface);
  • "Amesos_Dscpack" (use DSCPACK interface).
• "fact: level-of-fill": defines the level of fill for ILU factorizations. This is based on powers of the graph, so the value 0 means no-fill.
• "fact: ilut level-of-fill": defines the level of fill for ILUT factorizations. This is a ratio, so 1.0 means same number of nonzeros for the ILU factors as in the original matrix.
• "fact: ict level-of-fill": defines the level of fill for ICT factorizations, and is a fill-ratio as for ILUT. Currently this parameter is used in both Ifpack_IC and Ifpack_ICT.
• "fact: absolute threshold": defines the value \(\alpha\) to add to each diagonal element (times the signum of the actual diagonal element), for incomplete factorizations only.
• "fact: relative threshold": defines \(\rho\), so that the diagonal element of the matrix (without the contribution specified by "fact: absolute threshold") is multiplied by \(rho \).
• "fact: relax value": if different from zero, the elements dropped during the factorization process will be added to the diagonal term, multiplied by the specified value.
• "schwarz: combine mode": defines how values corresponding to overlapping nodes are handled (after the solution of the local problem). Any Epetra_CombineMode is valid. Users should set this value to Add if interested in a symmetric preconditioner. Otherwise, the default value of "Zero" usually results in better convergence. Other details on combine mode are reported. All Ifpack 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 Ifpack_IlukGraph object. The application of the inverse of the local matrix results in redundant approximations for any elements of y that correspond to rows that are part of more than one local ILU factor. The combine mode defines how these redundancies are handled using the Epetra_CombineMode enum. The default (Zero) means to zero out all values of y for rows that were not part of the original matrix row distribution.
• "schwarz: compute condest": it true, compute a cheap estimate of the condition number. Warning: this is not the REAL conditioner number, but only a very light-weight indication of how well the current preconditioner is expected to behave.
• "schwarz: filter singletons": it true, singletons in the local matrix are removed.
• "schwarz: reordering type": it none, no reordering is applied. If rcm or metis, a reverse Cuthill-McKee or METIS are used to reorder to local matrix.

The table below reports the C++ type and the default value of each parameter.

Parameter name	C++ type	Defaut value
"relaxation: type"	string	"Jacobi"
"relaxation: sweeps"	int	"1"
"relaxation: damping factor"	double	"1.0"
"relaxation: min diagonal value"	double	"0.0"
"relaxation: zero starting solution"	bool	"true"
"relaxation: backward mode"	bool	"false"
"relaxation: use l1"	bool	"false"
"relaxation: l1 eta"	double	"1.5"
"relaxation: number of local smoothing indices"	int	"A->NumMyRows()"
"relaxation: local smoothing indices"	int*	"0"
"partitioner: type"	string	"greedy"
"partitioner: overlap"	int	"0"
"partitioner: local parts"	int	"1"
"partitioner: root node"	int	"0"
"partitioner: use symmetric graph"	bool	"true"
"amesos: solver type"	string	"Amesos_Klu"
"fact: absolute threshold"	double	"0.0"
"fact: relative threshold"	double	"0.0"
"fact: relax value"	double	"0.0"
"schwarz: combine mode"	Epetra_CombineMode	"Zero"
"schwarz: compute condest"	bool	"true"
"schwarz: reordering type"	string	"none"
"schwarz: filter singletons"	string	"false"

Note

Whenever possible, the Teuchos::ParameterList given in input to SetParameters() is only parsed. There are cases, however, where the solver needs to store a copy of this list, because it will be used in a later stage. This happens, for example, in block relaxation schemes that use Amesos as local solver. As the the list given in input to SetParameters() is parsed or copied, this list can go out of scope after this function returns.

Brief Description IFPACK Classes

IFPACK classes are divided into classes for users and classes for developers.

Users' Tools

• Ifpack_PointRelaxation and Ifpack_BlockRelaxation should be considered to define relaxation preconditioners (or smoothers);
• Ifpack_IC, Ifpack_ICT are incomplete Cholesky factorizations;
• Ifpack_ILU, Ifpack_ILUT are incomplete LU factorizations;
• Ifpack_AdditiveSchwarz is the general handle to define overlapping Schwarz methods with overlap.
• The matrix to be preconditioned can be filtered before performing the actual factorization. Available filters are Ifpack_LocalFilter, Ifpack_DropFilter, Ifpack_SparsityFilter, Ifpack_DiagonalFilter, Ifpack_ReoderFilter, Ifpack_SingletonFilter;
• Function Ifpack_PrintSparsity() and Ifpack_PrintSparsity_Simple() can be used to print out the sparsity pattern of an Epetra_RowMatrix (as PS, or on screen);
• Function Ifpack_Analyze() can be used to analyze basic properties of a given Epetra_RowMatrix; see an example of output Usage of Ifpack_Analyze()..
• Function Ifpack_AnalyzeMatrixElements() can be used to evaluate the distribution of matrix elements; see an example of output Usage of Ifpack_AnalyzeMatrixElements().
• Function Ifpack_AnalyzeVectorElements() can be used to evaluate the distribution of vector elements; see an example of output Usage of Ifpack_AnalyzeVectorElements().

Developers' Tools

• Graphs are derived from the pure virtual class Ifpack_Graph (whose documentation details the concept of graphs for IFPACK). Concrete implementations are Ifpack_Graph_Epetra_CrsGraph and Ifpack_Graph_Epetra_RowMatrix;
• Containers are derived from the pure virtual class Ifpack_Container (whose documentation describes what a container is). Concrete implementations are Ifpack_DenseContainer and Ifpack_SparseContainer;
• Reorderings are derived from the pure virtual class Ifpack_Reordering. Concrete implementations are Ifpack_RCMReordering and Ifpack_METISReordering;
• Partitioners are derived from the pure virtual class Ifpack_Partitioner. Concrete implementations are Ifpack_LinearPartitioner, Ifpack_METISPartitioner, Ifpack_GreedyPartitioner, Ifpack_EquationPartitioner, and Ifpack_UserPartitioner. See Ifpack_OverlappingPartitioner;
• Ifpack_ReorderFilter and Ifpack_SingletonFilter are used in Ifpack_AdditiveSchwarz.

Error Table

The following table defines the class of IFPACK errors. Return code should be checked using the IFPACK_CHK_ERR() macro. In general terms, we follow this convention:

• a return value of 0 means that the called function or method successfully completed;
• a negative return function means that an error occurred;
• a positive value is a warning for the user.

Error Code	Meaning
-1	Generic Error (called method or function returned an error)
-2	Input data not valid (wrong parameter, out-of-bounds, wrong dimensions, matrix is not square,...)
-3	Data has not been correctly pre-processed
-4	Problem encountered during application of the algorithm (division by zero, out-of-bounds, ...)
-5	Memory allocation error
-98	Feature is not supported
-99	Feature is not implemented yet (check Known Bugs and Future Developments, or submit a bug)

Known Bugs and Future Developments

• IFPACK does not support Epetra matrices whose OperatorDomainMap() is different from the OperatorRangeMap() and from RowMatrixRowMap(). Handling of this particular case is under development for class Ifpack_AdditiveSchwarz only.
• not all IFPACK preconditioners work correctly when UseTranspose() == true. Ifpack_PointRelaxation and Ifpack_BlockRelaxation do not support UseTranspose, while Ifpack_ILU, Ifpack_ILUT and Ifpack_Amesos do.

Copyright

Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
license for use of this work by or on behalf of the U.S. Government.

This library is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as
published by the Free Software Foundation; either version 2.1 of the
License, or (at your option) any later version.

This library is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public
License along with this library; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
USA

Browse all of Ifpack as a single doxygen collection

You can browse all of Ifpack as a single doxygen collection. Warning, this is not the recommended way to learn about Ifpack software. However, this is a good way to browse the directory structure of ifpack, to locate files, etc.