dune-istl 2.8.0
|
Sequential ILU preconditioner. More...
#include <dune/istl/preconditioners.hh>
Public Types | |
typedef std::remove_const< M >::type | matrix_type |
The matrix type the preconditioner is for. More... | |
typedef matrix_type::block_type | block_type |
block type of matrix More... | |
typedef X | domain_type |
The domain type of the preconditioner. More... | |
typedef Y | range_type |
The range type of the preconditioner. More... | |
typedef X::field_type | field_type |
The field type of the preconditioner. More... | |
typedef Simd::Scalar< field_type > | scalar_field_type |
scalar type underlying the field_type More... | |
typedef FieldTraits< scalar_field_type >::real_type | real_field_type |
real scalar type underlying the field_type More... | |
typedef ILU::CRS< block_type, typename M::allocator_type > | CRS |
type of ILU storage More... | |
Public Member Functions | |
SeqILU (const M &A, real_field_type w, const bool resort=false) | |
Constructor. More... | |
SeqILU (const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration) | |
Constructor. More... | |
SeqILU (const M &A, const ParameterTree &config) | |
Constructor. More... | |
SeqILU (const M &A, int n, real_field_type w, const bool resort=false) | |
Constructor. More... | |
virtual void | pre (X &x, Y &b) |
Prepare the preconditioner. More... | |
virtual void | apply (X &v, const Y &d) |
Apply the preconditioner. More... | |
virtual void | post (X &x) |
Clean up. More... | |
virtual SolverCategory::Category | category () const |
Category of the preconditioner (see SolverCategory::Category) More... | |
Protected Attributes | |
std::unique_ptr< matrix_type > | ILU_ |
The ILU(n) decomposition of the matrix. As storage a BCRSMatrix is used. More... | |
CRS | lower_ |
The ILU(n) decomposition of the matrix. As storage a CRS structure is used. More... | |
CRS | upper_ |
std::vector< block_type, typename matrix_type::allocator_type > | inv_ |
const real_field_type | w_ |
The relaxation factor to use. More... | |
const bool | wNotIdentity_ |
true if w != 1.0 More... | |
Sequential ILU preconditioner.
Wraps the naked ISTL generic ILU preconditioner into the solver framework.
M | The matrix type to operate on |
X | Type of the update |
Y | Type of the defect |
l | Ignored. Just there to have the same number of template arguments as other preconditioners. |
typedef matrix_type::block_type Dune::SeqILU< M, X, Y, l >::block_type |
block type of matrix
typedef ILU::CRS< block_type , typename M::allocator_type> Dune::SeqILU< M, X, Y, l >::CRS |
type of ILU storage
typedef X Dune::SeqILU< M, X, Y, l >::domain_type |
The domain type of the preconditioner.
typedef X::field_type Dune::SeqILU< M, X, Y, l >::field_type |
The field type of the preconditioner.
typedef std::remove_const<M>::type Dune::SeqILU< M, X, Y, l >::matrix_type |
The matrix type the preconditioner is for.
typedef Y Dune::SeqILU< M, X, Y, l >::range_type |
The range type of the preconditioner.
typedef FieldTraits<scalar_field_type>::real_type Dune::SeqILU< M, X, Y, l >::real_field_type |
real scalar type underlying the field_type
typedef Simd::Scalar<field_type> Dune::SeqILU< M, X, Y, l >::scalar_field_type |
scalar type underlying the field_type
|
inline |
Constructor.
Constructor invoking ILU(0) gets all parameters to operate the prec.
A | The matrix to operate on. |
w | The relaxation factor. |
resort | true if a resort of the computed ILU for improved performance should be done. |
|
inline |
Constructor.
A | The assembled linear operator to use. |
configuration | ParameterTree containing preconditioner parameters. |
ParameterTree Key | Meaning |
---|---|
n | The order of the ILU decomposition. default=0 |
relaxation | The relaxation factor. default=1.0 |
resort | True if a resort of the computed ILU for improved performance should be done. default=false |
See ISTL_Factory for the ParameterTree layout and examples.
|
inline |
Constructor.
A | The matrix to operate on. |
configuration | ParameterTree containing preconditioner parameters. |
ParameterTree Key | Meaning |
---|---|
n | The order of the ILU decomposition. default=0 |
relaxation | The relaxation factor. default=1.0 |
resort | True if a resort of the computed ILU for improved performance should be done. default=false |
See ISTL_Factory for the ParameterTree layout and examples.
|
inline |
|
inlinevirtual |
Apply the preconditioner.
Apply one step of the preconditioner to the system A(v)=d.
On entry v=0 and d=b-A(x) (although this might not be computed in that way. On exit v contains the update, i.e one step computes where
is the approximate inverse of the operator
characterizing the preconditioner.
[out] | v | The update to be computed |
d | The current defect. |
Implements Dune::Preconditioner< X, Y >.
|
inlinevirtual |
Category of the preconditioner (see SolverCategory::Category)
Implements Dune::Preconditioner< X, Y >.
|
inlinevirtual |
Clean up.
Clean up.
This method is called after the last apply call for the linear system to be solved. Memory may be deallocated safely here. x is the solution of the linear equation.
x | The right hand side of the equation. |
Implements Dune::Preconditioner< X, Y >.
|
inlinevirtual |
Prepare the preconditioner.
Prepare the preconditioner.
A solver solves a linear operator equation A(x)=b by applying one or several steps of the preconditioner. The method pre() is called before the first apply operation. b and x are right hand side and solution vector of the linear system respectively. It may. e.g., scale the system, allocate memory or compute a (I)LU decomposition. Note: The ILU decomposition could also be computed in the constructor or with a separate method of the derived method if several linear systems with the same matrix are to be solved.
x | The left hand side of the equation. |
b | The right hand side of the equation. |
Implements Dune::Preconditioner< X, Y >.
|
protected |
The ILU(n) decomposition of the matrix. As storage a BCRSMatrix is used.
|
protected |
|
protected |
The ILU(n) decomposition of the matrix. As storage a CRS structure is used.
|
protected |
|
protected |
The relaxation factor to use.
|
protected |
true if w != 1.0