Ifpack2 Templated Preconditioning Package
Version 1.0
|
Construct a level filled graph for use in computing an ILU(k) incomplete factorization. More...
#include <Ifpack2_IlukGraph.hpp>
Public Types | |
typedef Tpetra::RowGraph < local_ordinal_type, global_ordinal_type, node_type > | row_graph_type |
Tpetra::RowGraph specialization used by this class. More... | |
typedef Tpetra::CrsGraph < local_ordinal_type, global_ordinal_type, node_type > | crs_graph_type |
Tpetra::CrsGraph specialization used by this class. More... | |
Public Member Functions | |
IlukGraph (const Teuchos::RCP< const GraphType > &G, const int levelFill, const int levelOverlap, const double overalloc=2.) | |
Constructor. More... | |
virtual | ~IlukGraph () |
IlukGraph Destructor. More... | |
void | setParameters (const Teuchos::ParameterList ¶meterlist) |
Set parameters. More... | |
void | initialize () |
Set up the graph structure of the L and U factors. More... | |
int | getLevelFill () const |
The level of fill used to construct this graph. More... | |
int | getLevelOverlap () const |
The level of overlap used to construct this graph. More... | |
Teuchos::RCP< const GraphType > | getA_Graph () const |
Returns the original graph given. More... | |
Teuchos::RCP< crs_graph_type > | getL_Graph () const |
Returns the graph of lower triangle of the ILU(k) graph as a Tpetra::CrsGraph. More... | |
Teuchos::RCP< crs_graph_type > | getU_Graph () const |
Returns the graph of upper triangle of the ILU(k) graph as a Tpetra::CrsGraph. More... | |
Teuchos::RCP< const crs_graph_type > | getOverlapGraph () const |
Returns the the overlapped graph. More... | |
size_t | getNumGlobalDiagonals () const |
Returns the global number of diagonals in the ILU(k) graph. More... | |
Construct a level filled graph for use in computing an ILU(k) incomplete factorization.
GraphType | A specialization of Tpetra::CrsGraph (tested) or Tpetra::RowGraph (not tested). |
Ifpack2::IlukGraph enables the construction of matrix graphs using level-fill algorithms. The only function required for construction is a getGlobalRowView capability, i.e., the graph that is passed in to the constructor must implement the Tpetra::RowGraph interface defined in Tpetra_RowGraph.hpp.
Constructing an IlukGraph is a two step process. First, call the constructor, passing in a CrsGraph object and an integer indicating the desired level of fill. Then, call the initialize() method to complete the process. This naturally matches the three-stage initialization of Preconditioner objects, and in particular of RILUK.
It is worth noting that an IlukGraph object creates two Tpetra::CrsGraph objects containing L and U, the graphs for the lower and upper triangular parts of the ILU(k) graph. Thus, it is possible to manually insert and delete graph entries in L and U via the Tpetra::CrsGraph InsertIndices and RemoveIndices functions. However, in this case FillComplete must be called before the graph is used for subsequent operations.
typedef Tpetra::RowGraph<local_ordinal_type, global_ordinal_type, node_type> Ifpack2::IlukGraph< GraphType, KKHandleType >::row_graph_type |
Tpetra::RowGraph specialization used by this class.
typedef Tpetra::CrsGraph<local_ordinal_type, global_ordinal_type, node_type> Ifpack2::IlukGraph< GraphType, KKHandleType >::crs_graph_type |
Tpetra::CrsGraph specialization used by this class.
Ifpack2::IlukGraph< GraphType, KKHandleType >::IlukGraph | ( | const Teuchos::RCP< const GraphType > & | G, |
const int | levelFill, | ||
const int | levelOverlap, | ||
const double | overalloc = 2. |
||
) |
Constructor.
Create a IlukGraph object using the input graph and specified level of fill.
G | [in] An existing graph. |
levelFill | [in] The level of fill to compute; the k of ILU(k). |
levelOverlap | [in] The level of overlap between subdomains. |
overalloc | [in] The estimated number of nonzeros per row in the resulting matrices is (maxNodeNumRowEntries of the input * overalloc^levelFill). Must be greater than 1. Smaller values are more conservative with memory, but may require recomputation if the estimate is too low. Default value is two. |
|
virtual |
IlukGraph Destructor.
void Ifpack2::IlukGraph< GraphType, KKHandleType >::setParameters | ( | const Teuchos::ParameterList & | parameterlist | ) |
Set parameters.
This method recognizes two parameter names: Level_fill and Level_overlap. Both are case insensitive, and in both cases the parameter must have type int.
void Ifpack2::IlukGraph< GraphType, KKHandleType >::initialize | ( | ) |
Set up the graph structure of the L and U factors.
This method is called "initialize" by analogy with Preconditioner, where initialize() computes the symbolic (incomplete) factorization, and compute() computes the corresponding numeric factorization. IlukGraph is just a graph, so it can only compute a symbolic factorization (i.e., the graph structure of the factorization). Hence, it implements initialize(), but not compute(). RILUK calls IlukGraph's initialize() method in its own initialize() method, as one would expect.
|
inline |
The level of fill used to construct this graph.
|
inline |
The level of overlap used to construct this graph.
|
inline |
Returns the original graph given.
|
inline |
Returns the graph of lower triangle of the ILU(k) graph as a Tpetra::CrsGraph.
|
inline |
Returns the graph of upper triangle of the ILU(k) graph as a Tpetra::CrsGraph.
|
inline |
Returns the the overlapped graph.
|
inline |
Returns the global number of diagonals in the ILU(k) graph.