A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
More...
#include <ParallelOverlappingILU0.hpp>
|
using | matrix_type = typename std::remove_const< Matrix >::type |
| The matrix type the preconditioner is for.
|
|
using | domain_type = Domain |
| The domain type of the preconditioner.
|
|
using | range_type = Range |
| The range type of the preconditioner.
|
|
using | field_type = typename Domain::field_type |
| The field type of the preconditioner.
|
|
using | block_type = typename matrix_type::block_type |
|
using | size_type = typename matrix_type::size_type |
|
|
Dune::SolverCategory::Category | category () const override |
|
| ParallelOverlappingILU0 (const Matrix &A, const int n, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true) |
| Constructor. More...
|
|
| ParallelOverlappingILU0 (const Matrix &A, const ParallelInfo &comm, const int n, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true) |
| Constructor gets all parameters to operate the prec. More...
|
|
| ParallelOverlappingILU0 (const Matrix &A, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true) |
| Constructor. More...
|
|
| ParallelOverlappingILU0 (const Matrix &A, const ParallelInfo &comm, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true) |
| Constructor. More...
|
|
| ParallelOverlappingILU0 (const Matrix &A, const ParallelInfo &comm, const field_type w, MILU_VARIANT milu, size_type interiorSize, bool redblack=false, bool reorder_sphere=true) |
| Constructor. More...
|
|
void | pre (Domain &, Range &) override |
| Prepare the preconditioner. More...
|
|
void | apply (Domain &v, const Range &d) override |
| Apply the preconditoner. More...
|
|
template<class V > |
void | copyOwnerToAll (V &v) const |
|
void | post (Range &) override |
| Clean up. More...
|
|
void | update () override |
|
|
Range & | reorderD (const Range &d) |
| Reorder D if needed and return a reference to it.
|
|
Domain & | reorderV (Domain &v) |
| Reorder V if needed and return a reference to it.
|
|
void | reorderBack (const Range &reorderedV, Range &v) |
|
|
CRS | lower_ |
| The ILU0 decomposition of the matrix.
|
|
CRS | upper_ |
|
std::vector< block_type > | inv_ |
|
std::vector< std::size_t > | ordering_ |
| the reordering of the unknowns
|
|
Range | reorderedD_ |
| The reordered right hand side.
|
|
Domain | reorderedV_ |
| The reordered left hand side.
|
|
const ParallelInfo * | comm_ |
|
const field_type | w_ |
| The relaxation factor to use.
|
|
const bool | relaxation_ |
|
size_type | interiorSize_ |
|
const Matrix * | A_ |
|
int | iluIteration_ |
|
MILU_VARIANT | milu_ |
|
bool | redBlack_ |
|
bool | reorderSphere_ |
|
template<class Matrix, class Domain, class Range, class ParallelInfoT>
class Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
This preconditioner differs from a ParallelRestrictedOverlappingSchwarz with Dune:SeqILU0 in the following way: During apply we make sure that the current residual is consistent (i.e. each process knows the same value for each index. Then we solve Ly = d for y and make y consistent again. Last we solve Ux = y and make sure that x is consistent. In contrast for ParallelRestrictedOverlappingSchwarz we solve (LU)x = d for x without forcing consistency between the two steps.
- Template Parameters
-
Matrix | The type of the Matrix. |
Domain | The type of the Vector representing the domain. |
Range | The type of the Vector representing the range. |
ParallelInfo | The type of the parallel information object used, e.g. Dune::OwnerOverlapCommunication |
◆ ParallelOverlappingILU0() [1/5]
template<class Matrix , class Domain , class Range , class ParallelInfoT >
Constructor.
Constructor gets all parameters to operate the prec.
- Parameters
-
A | The matrix to operate on. |
n | ILU fill in level (for testing). This does not work in parallel. |
w | The relaxation factor. |
milu | The modified ILU variant to use. 0 means traditional ILU. |
- See also
- MILU_VARIANT.
- Parameters
-
redblack | Whether to use a red-black ordering. |
reorder_sphere | If true, we start the reordering at a root node. The vertices on each layer aound it (same distance) are ordered consecutivly. If false, we preserver the order of the vertices with the same color. |
◆ ParallelOverlappingILU0() [2/5]
template<class Matrix , class Domain , class Range , class ParallelInfoT >
Constructor gets all parameters to operate the prec.
- Parameters
-
A | The matrix to operate on. |
comm | communication object, e.g. Dune::OwnerOverlapCopyCommunication |
n | ILU fill in level (for testing). This does not work in parallel. |
w | The relaxation factor. |
milu | The modified ILU variant to use. 0 means traditional ILU. |
- See also
- MILU_VARIANT.
- Parameters
-
redblack | Whether to use a red-black ordering. |
reorder_sphere | If true, we start the reordering at a root node. The vertices on each layer aound it (same distance) are ordered consecutivly. If false, we preserver the order of the vertices with the same color. |
◆ ParallelOverlappingILU0() [3/5]
template<class Matrix , class Domain , class Range , class ParallelInfoT >
Constructor.
Constructor gets all parameters to operate the prec.
- Parameters
-
A | The matrix to operate on. |
w | The relaxation factor. |
milu | The modified ILU variant to use. 0 means traditional ILU. |
- See also
- MILU_VARIANT.
- Parameters
-
redblack | Whether to use a red-black ordering. |
reorder_sphere | If true, we start the reordering at a root node. The vertices on each layer aound it (same distance) are ordered consecutivly. If false, we preserver the order of the vertices with the same color. |
◆ ParallelOverlappingILU0() [4/5]
template<class Matrix , class Domain , class Range , class ParallelInfoT >
Constructor.
Constructor gets all parameters to operate the prec.
- Parameters
-
A | The matrix to operate on. |
comm | communication object, e.g. Dune::OwnerOverlapCopyCommunication |
w | The relaxation factor. |
milu | The modified ILU variant to use. 0 means traditional ILU. |
- See also
- MILU_VARIANT.
- Parameters
-
redblack | Whether to use a red-black ordering. |
reorder_sphere | If true, we start the reordering at a root node. The vertices on each layer aound it (same distance) are ordered consecutivly. If false, we preserver the order of the vertices with the same color. |
◆ ParallelOverlappingILU0() [5/5]
template<class Matrix , class Domain , class Range , class ParallelInfoT >
Constructor.
Constructor gets all parameters to operate the prec.
- Parameters
-
A | The matrix to operate on. |
n | ILU fill in level (for testing). This does not work in parallel. |
w | The relaxation factor. |
milu | The modified ILU variant to use. 0 means traditional ILU. |
- See also
- MILU_VARIANT.
- Parameters
-
interiorSize | The number of interior/owner rows in the matrix. |
redblack | Whether to use a red-black ordering. |
reorder_sphere | If true, we start the reordering at a root node. The vertices on each layer aound it (same distance) are ordered consecutivly. If false, we preserver the order of the vertices with the same color. |
◆ apply()
template<class Matrix , class Domain , class Range , class ParallelInfoT >
◆ post()
template<class Matrix , class Domain , class Range , class ParallelInfoT >
◆ pre()
template<class Matrix , class Domain , class Range , class ParallelInfoT >
Prepare the preconditioner.
The documentation for this class was generated from the following files: