My Project
Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT > Class Template Reference

A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as. More...

#include <ParallelOverlappingILU0.hpp>

Inheritance diagram for Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >:
Dune::PreconditionerWithUpdate< Domain, Range >

Classes

struct  CRS
 

Public Types

typedef std::remove_const< Matrix >::type matrix_type
 The matrix type the preconditioner is for.
 
typedef Domain domain_type
 The domain type of the preconditioner.
 
typedef Range range_type
 The range type of the preconditioner.
 
typedef Domain::field_type field_type
 The field type of the preconditioner.
 
typedef matrix_type::block_type block_type
 
typedef matrix_type::size_type size_type
 

Public Member Functions

Dune::SolverCategory::Category category () const override
 
template<class BlockType , class Alloc >
 ParallelOverlappingILU0 (const Dune::BCRSMatrix< BlockType, Alloc > &A, const int n, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true)
 Constructor. More...
 
template<class BlockType , class Alloc >
 ParallelOverlappingILU0 (const Dune::BCRSMatrix< BlockType, Alloc > &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...
 
template<class BlockType , class Alloc >
 ParallelOverlappingILU0 (const Dune::BCRSMatrix< BlockType, Alloc > &A, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true)
 Constructor. More...
 
template<class BlockType , class Alloc >
 ParallelOverlappingILU0 (const Dune::BCRSMatrix< BlockType, Alloc > &A, const ParallelInfo &comm, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true)
 Constructor. More...
 
template<class BlockType , class Alloc >
 ParallelOverlappingILU0 (const Dune::BCRSMatrix< BlockType, Alloc > &A, const ParallelInfo &comm, const field_type w, MILU_VARIANT milu, size_type interiorSize, bool redblack=false, bool reorder_sphere=true)
 Constructor. More...
 
virtual void pre (Domain &x, Range &b) override
 Prepare the preconditioner. More...
 
virtual void apply (Domain &v, const Range &d) override
 Apply the preconditoner. More...
 
template<class V >
void copyOwnerToAll (V &v) const
 
virtual void post (Range &x) override
 Clean up. More...
 
virtual void update () override
 

Protected Member Functions

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)
 

Protected Attributes

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_
 

Detailed Description

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 follwing way: During apply we make sure that the current residual is consistent (i.e. each process knows the same value for each index. The 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
MatrixThe type of the Matrix.
DomainThe type of the Vector representing the domain.
RangeThe type of the Vector representing the range.
ParallelInfoThe type of the parallel information object used, e.g. Dune::OwnerOverlapCommunication

Constructor & Destructor Documentation

◆ ParallelOverlappingILU0() [1/5]

template<class Matrix , class Domain , class Range , class ParallelInfoT >
template<class BlockType , class Alloc >
Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::ParallelOverlappingILU0 ( const Dune::BCRSMatrix< BlockType, Alloc > &  A,
const int  n,
const field_type  w,
MILU_VARIANT  milu,
bool  redblack = false,
bool  reorder_sphere = true 
)
inline

Constructor.

Constructor gets all parameters to operate the prec.

Parameters
AThe matrix to operate on.
nILU fill in level (for testing). This does not work in parallel.
wThe relaxation factor.
miluThe modified ILU variant to use. 0 means traditional ILU.
See also
MILU_VARIANT.
Parameters
redblackWhether to use a red-black ordering.
reorder_sphereIf 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 >
template<class BlockType , class Alloc >
Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::ParallelOverlappingILU0 ( const Dune::BCRSMatrix< BlockType, Alloc > &  A,
const ParallelInfo &  comm,
const int  n,
const field_type  w,
MILU_VARIANT  milu,
bool  redblack = false,
bool  reorder_sphere = true 
)
inline

Constructor gets all parameters to operate the prec.

Parameters
AThe matrix to operate on.
commcommunication object, e.g. Dune::OwnerOverlapCopyCommunication
nILU fill in level (for testing). This does not work in parallel.
wThe relaxation factor.
miluThe modified ILU variant to use. 0 means traditional ILU.
See also
MILU_VARIANT.
Parameters
redblackWhether to use a red-black ordering.
reorder_sphereIf 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 >
template<class BlockType , class Alloc >
Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::ParallelOverlappingILU0 ( const Dune::BCRSMatrix< BlockType, Alloc > &  A,
const field_type  w,
MILU_VARIANT  milu,
bool  redblack = false,
bool  reorder_sphere = true 
)
inline

Constructor.

Constructor gets all parameters to operate the prec.

Parameters
AThe matrix to operate on.
wThe relaxation factor.
miluThe modified ILU variant to use. 0 means traditional ILU.
See also
MILU_VARIANT.
Parameters
redblackWhether to use a red-black ordering.
reorder_sphereIf 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 >
template<class BlockType , class Alloc >
Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::ParallelOverlappingILU0 ( const Dune::BCRSMatrix< BlockType, Alloc > &  A,
const ParallelInfo &  comm,
const field_type  w,
MILU_VARIANT  milu,
bool  redblack = false,
bool  reorder_sphere = true 
)
inline

Constructor.

Constructor gets all parameters to operate the prec.

Parameters
AThe matrix to operate on.
commcommunication object, e.g. Dune::OwnerOverlapCopyCommunication
wThe relaxation factor.
miluThe modified ILU variant to use. 0 means traditional ILU.
See also
MILU_VARIANT.
Parameters
redblackWhether to use a red-black ordering.
reorder_sphereIf 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 >
template<class BlockType , class Alloc >
Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::ParallelOverlappingILU0 ( const Dune::BCRSMatrix< BlockType, Alloc > &  A,
const ParallelInfo &  comm,
const field_type  w,
MILU_VARIANT  milu,
size_type  interiorSize,
bool  redblack = false,
bool  reorder_sphere = true 
)
inline

Constructor.

Constructor gets all parameters to operate the prec.

Parameters
AThe matrix to operate on.
nILU fill in level (for testing). This does not work in parallel.
wThe relaxation factor.
miluThe modified ILU variant to use. 0 means traditional ILU.
See also
MILU_VARIANT.
Parameters
interiorSizeThe number of interior/owner rows in the matrix.
redblackWhether to use a red-black ordering.
reorder_sphereIf 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.

Member Function Documentation

◆ apply()

template<class Matrix , class Domain , class Range , class ParallelInfoT >
virtual void Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::apply ( Domain &  v,
const Range &  d 
)
inlineoverridevirtual

Apply the preconditoner.

◆ post()

template<class Matrix , class Domain , class Range , class ParallelInfoT >
virtual void Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::post ( Range &  x)
inlineoverridevirtual

Clean up.

◆ pre()

template<class Matrix , class Domain , class Range , class ParallelInfoT >
virtual void Opm::ParallelOverlappingILU0< Matrix, Domain, Range, ParallelInfoT >::pre ( Domain &  x,
Range &  b 
)
inlineoverridevirtual

Prepare the preconditioner.


The documentation for this class was generated from the following file: