5#ifndef DUNE_ISTL_PRECONDITIONERS_HH
6#define DUNE_ISTL_PRECONDITIONERS_HH
15#include <dune/common/simd/simd.hh>
16#include <dune/common/parametertree.hh>
73 template<
class O,
int c = -1>
75 public Preconditioner<typename O::domain_type, typename O::range_type>
96 : inverse_operator_(inverse_operator)
99 DUNE_THROW(InvalidStateException,
"User-supplied solver category does not match that of the given inverse operator");
109 inverse_operator_.apply(v, copy, res);
141 template<
class M,
class X,
class Y,
int l=1>
165 : _A_(A), _n(n), _w(w)
184 :
SeqSSOR(A->getmat(), configuration)
200 SeqSSOR (
const M& A,
const ParameterTree& configuration)
209 virtual void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b)
217 virtual void apply (X& v,
const Y& d)
219 for (
int i=0; i<_n; i++) {
230 virtual void post ([[maybe_unused]] X& x)
261 template<
class M,
class X,
class Y,
int l=1>
285 : _A_(A), _n(n), _w(w)
304 :
SeqSOR(A->getmat(), configuration)
320 SeqSOR (
const M& A,
const ParameterTree& configuration)
329 virtual void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b)
337 virtual void apply (X& v,
const Y& d)
339 this->
template apply<true>(v,d);
350 template<
bool forward>
354 for (
int i=0; i<_n; i++) {
358 for (
int i=0; i<_n; i++) {
368 virtual void post ([[maybe_unused]] X& x)
398 template<
class M,
class X,
class Y,
int l=1>
412 template<
class M,
class X,
class Y,
int l=1>
436 : _A_(A), _n(n), _w(w)
455 :
SeqJac(A->getmat(), configuration)
471 SeqJac (
const M& A,
const ParameterTree& configuration)
480 virtual void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b)
488 virtual void apply (X& v,
const Y& d)
490 for (
int i=0; i<_n; i++) {
500 virtual void post ([[maybe_unused]] X& x)
562 template <
class M,
class X,
class Y,
int l = 1>
595 Dinv_.resize(_A_.N());
613 :
SeqDILU(A->getmat(), configuration)
629 SeqDILU(
const M &A,
const ParameterTree &config)
639 virtual void pre([[maybe_unused]] X &x, [[maybe_unused]] Y &b)
648 virtual void apply(X &v,
const Y &d)
664 virtual void post([[maybe_unused]] X &x)
696 template<
class M,
class X,
class Y,
int l=1>
727 :
SeqILU( A, 0, w, resort )
746 :
SeqILU(A->getmat(), configuration)
763 SeqILU(
const M& A,
const ParameterTree& config)
766 config.
get(
"resort", false))
813 virtual void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b)
821 virtual void apply (X& v,
const Y& d)
843 virtual void post ([[maybe_unused]] X& x)
854 std::unique_ptr< matrix_type >
ILU_;
859 std::vector< block_type, typename matrix_type::allocator_type >
inv_;
877 template<
class X,
class Y>
920 virtual void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b)
928 virtual void apply (X& v,
const Y& d)
939 virtual void post ([[maybe_unused]] X& x)
953 using D =
typename Dune::TypeListElement<1,
decltype(tl)>::type;
954 using R =
typename Dune::TypeListElement<2,
decltype(tl)>::type;
955 return std::make_shared<Richardson<D,R>>(config);
969 template<
class M,
class X,
class Y >
1003 :
SeqILDL(A->getmat(), configuration)
1031 : decomposition_( A.N(), A.M(),
matrix_type::random ),
1035 for(
auto i = A.begin(), iend = A.end(); i != iend; ++i )
1037 const auto &A_i = *i;
1038 const auto ij = A_i.find( i.index() );
1039 if( ij != A_i.end() )
1040 decomposition_.setrowsize( i.index(), ij.offset()+1 );
1042 DUNE_THROW(
ISTLError,
"diagonal entry missing" );
1044 decomposition_.endrowsizes();
1047 for(
auto i = A.begin(), iend = A.end(); i != iend; ++i )
1049 const auto &A_i = *i;
1050 for(
auto ij = A_i.begin(); ij.index() < i.index() ; ++ij )
1051 decomposition_.addindex( i.index(), ij.index() );
1052 decomposition_.addindex( i.index(), i.index() );
1054 decomposition_.endindices();
1058 for(
auto row = decomposition_.begin(), rowend = decomposition_.end(); row != rowend; ++row, ++i )
1060 auto ij = i->begin();
1061 for(
auto col = row->begin(), colend = row->end();
col != colend; ++
col, ++ij )
1070 void pre ([[maybe_unused]] X &x, [[maybe_unused]] Y &b)
override
1074 void apply ( X &v,
const Y &d )
override
1081 void post ([[maybe_unused]] X &x)
override
The diagonal incomplete LU factorization kernels.
Some handy generic functions for ISTL matrices.
#define DUNE_REGISTER_PRECONDITIONER(name,...)
Definition solverregistry.hh:16
Define general, extensible interface for inverse operators.
Incomplete LDL decomposition.
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
The incomplete LU factorization kernels.
Col col
Definition matrixmatrix.hh:351
void bsorb(const M &A, X &x, const Y &b, const K &w)
SSOR step.
Definition gsetc.hh:646
void dbjac(const M &A, X &x, const Y &b, const K &w)
Jacobi step.
Definition gsetc.hh:658
void bsorf(const M &A, X &x, const Y &b, const K &w)
SOR step.
Definition gsetc.hh:634
Definition allocator.hh:11
void bildl_decompose(Matrix &A)
compute ILDL decomposition of a symmetric matrix A
Definition ildl.hh:88
PropertyMapTypeSelector< Amg::VertexVisitedTag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > >::Type get(const Amg::VertexVisitedTag &tag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > &graph)
Definition dependency.hh:293
void bildl_backsolve(const Matrix &A, X &v, const Y &d, bool isLowerTriangular=false)
Definition ildl.hh:149
void blockDILUDecomposition(M &A, std::vector< typename M::block_type > &Dinv_)
Definition dilu.hh:51
void blockDILUBacksolve(const M &A, const std::vector< typename M::block_type > Dinv_, X &v, const Y &d)
Definition dilu.hh:107
void convertToCRS(const M &A, CRS &lower, CRS &upper, InvVector &inv)
convert ILU decomposition into CRS format for lower and upper triangular and inverse.
Definition ilu.hh:307
void blockILUBacksolve(const M &A, X &v, const Y &d)
LU backsolve with stored inverse.
Definition ilu.hh:94
void blockILU0Decomposition(M &A)
compute ILU decomposition of A. A is overwritten by its decomposition
Definition ilu.hh:33
void blockILUDecomposition(const M &A, int n, M &ILU)
Definition ilu.hh:167
compile-time parameter for block recursion depth
Definition gsetc.hh:45
a simple compressed row storage matrix class
Definition ilu.hh:259
derive error class from the base class in common
Definition istlexception.hh:19
static void check(const Matrix &mat)
Check whether the a matrix has diagonal values on blocklevel recursion levels.
Definition matrixutils.hh:53
A linear operator exporting itself in matrix form.
Definition operators.hh:111
Base class for matrix free definition of preconditioners.
Definition preconditioner.hh:33
Turns an InverseOperator into a Preconditioner.
Definition preconditioners.hh:76
O::range_type range_type
The range type of the preconditioner.
Definition preconditioners.hh:81
O::domain_type domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:79
virtual void post(domain_type &)
Clean up.
Definition preconditioners.hh:112
range_type::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:83
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:87
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:85
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:116
virtual void pre(domain_type &, range_type &)
Prepare the preconditioner.
Definition preconditioners.hh:102
InverseOperator2Preconditioner(InverseOperator &inverse_operator)
Construct the preconditioner from the solver.
Definition preconditioners.hh:95
O InverseOperator
type of the wrapped inverse operator
Definition preconditioners.hh:89
virtual void apply(domain_type &v, const range_type &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition preconditioners.hh:105
Sequential SSOR preconditioner.
Definition preconditioners.hh:142
virtual void post(X &x)
Clean up.
Definition preconditioners.hh:230
SeqSSOR(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:183
SeqSSOR(const M &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:200
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:234
X::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:151
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:153
X domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:147
M matrix_type
The matrix type the preconditioner is for.
Definition preconditioners.hh:145
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition preconditioners.hh:217
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition preconditioners.hh:209
Y range_type
The range type of the preconditioner.
Definition preconditioners.hh:149
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:155
SeqSSOR(const M &A, int n, real_field_type w)
Constructor.
Definition preconditioners.hh:164
Sequential SOR preconditioner.
Definition preconditioners.hh:262
SeqSOR(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:303
M matrix_type
The matrix type the preconditioner is for.
Definition preconditioners.hh:265
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:275
void apply(X &v, const Y &d)
Apply the preconditioner in a special direction.
Definition preconditioners.hh:351
X domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:267
virtual void post(X &x)
Clean up.
Definition preconditioners.hh:368
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition preconditioners.hh:329
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:273
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:372
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition preconditioners.hh:337
Y range_type
The range type of the preconditioner.
Definition preconditioners.hh:269
SeqSOR(const M &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:320
X::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:271
SeqSOR(const M &A, int n, real_field_type w)
Constructor.
Definition preconditioners.hh:284
The sequential jacobian preconditioner.
Definition preconditioners.hh:413
virtual void post(X &x)
Clean up.
Definition preconditioners.hh:500
SeqJac(const M &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:471
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition preconditioners.hh:488
M matrix_type
The matrix type the preconditioner is for.
Definition preconditioners.hh:416
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:424
SeqJac(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:454
X::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:422
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition preconditioners.hh:480
X domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:418
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:426
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:504
SeqJac(const M &A, int n, real_field_type w)
Constructor.
Definition preconditioners.hh:435
Y range_type
The range type of the preconditioner.
Definition preconditioners.hh:420
Sequential DILU preconditioner.
Definition preconditioners.hh:564
const bool wNotIdentity_
true if w != 1.0
Definition preconditioners.hh:681
typename X::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:576
const real_field_type _w
The relaxation factor to use.
Definition preconditioners.hh:679
const M & _A_
The matrix we operate on.
Definition preconditioners.hh:677
SeqDILU(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:612
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:579
virtual void post(X &x)
Clean up.
Definition preconditioners.hh:664
typename FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:581
M matrix_type
The matrix type the preconditioner is for.
Definition preconditioners.hh:567
X domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:571
typename matrix_type::block_type block_type
block type of matrix
Definition preconditioners.hh:569
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition preconditioners.hh:648
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition preconditioners.hh:639
SeqDILU(const M &A, real_field_type w)
Constructor.
Definition preconditioners.hh:589
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:669
Y range_type
The range type of the preconditioner.
Definition preconditioners.hh:573
SeqDILU(const M &A, const ParameterTree &config)
Constructor.
Definition preconditioners.hh:629
std::vector< block_type > Dinv_
Definition preconditioners.hh:675
Sequential ILU preconditioner.
Definition preconditioners.hh:697
virtual void post(X &x)
Clean up.
Definition preconditioners.hh:843
SeqILU(const M &A, int n, real_field_type w, const bool resort=false)
Constructor.
Definition preconditioners.hh:777
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition preconditioners.hh:813
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition preconditioners.hh:821
ILU::CRS< block_type, typename M::allocator_type > CRS
type of ILU storage
Definition preconditioners.hh:717
Y range_type
The range type of the preconditioner.
Definition preconditioners.hh:706
CRS lower_
The ILU(n) decomposition of the matrix. As storage a CRS structure is used.
Definition preconditioners.hh:857
const bool wNotIdentity_
true if w != 1.0
Definition preconditioners.hh:864
SeqILU(const M &A, const ParameterTree &config)
Constructor.
Definition preconditioners.hh:763
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition preconditioners.hh:700
matrix_type::block_type block_type
block type of matrix
Definition preconditioners.hh:702
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:714
X::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:709
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:847
SeqILU(const M &A, real_field_type w, const bool resort=false)
Constructor.
Definition preconditioners.hh:726
const real_field_type w_
The relaxation factor to use.
Definition preconditioners.hh:862
SeqILU(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:745
X domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:704
std::vector< block_type, typename matrix_type::allocator_type > inv_
Definition preconditioners.hh:859
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:712
std::unique_ptr< matrix_type > ILU_
The ILU(n) decomposition of the matrix. As storage a BCRSMatrix is used.
Definition preconditioners.hh:854
CRS upper_
Definition preconditioners.hh:858
Richardson preconditioner.
Definition preconditioners.hh:878
X::field_type field_type
The field type of the preconditioner.
Definition preconditioners.hh:885
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:943
Y range_type
The range type of the preconditioner.
Definition preconditioners.hh:883
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition preconditioners.hh:920
Richardson(real_field_type w=1.0)
Constructor.
Definition preconditioners.hh:896
virtual void post(X &x)
Clean up.
Definition preconditioners.hh:939
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:889
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:887
Richardson(const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:911
X domain_type
The domain type of the preconditioner.
Definition preconditioners.hh:881
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition preconditioners.hh:928
sequential ILDL preconditioner
Definition preconditioners.hh:972
SeqILDL(const matrix_type &A, const ParameterTree &config)
Constructor.
Definition preconditioners.hh:1018
SeqILDL(const matrix_type &A, real_field_type relax=real_field_type(1))
constructor
Definition preconditioners.hh:1030
X domain_type
domain type of the preconditioner
Definition preconditioners.hh:980
void post(X &x) override
Clean up.
Definition preconditioners.hh:1081
Y range_type
range type of the preconditioner
Definition preconditioners.hh:982
std::remove_const_t< M > matrix_type
type of matrix the preconditioner is for
Definition preconditioners.hh:978
void apply(X &v, const Y &d) override
Apply one step of the preconditioner to the system A(v)=d.
Definition preconditioners.hh:1074
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition preconditioners.hh:988
SeqILDL(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition preconditioners.hh:1002
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition preconditioners.hh:1070
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition preconditioners.hh:986
X::field_type field_type
field type of the preconditioner
Definition preconditioners.hh:984
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition preconditioners.hh:1085
Statistics about the application of an inverse operator.
Definition solver.hh:50
Abstract base class for all solvers.
Definition solver.hh:101
Category
Definition solvercategory.hh:23
@ sequential
Category for sequential solvers.
Definition solvercategory.hh:25
static Category category(const OP &op, decltype(op.category()) *=nullptr)
Helperfunction to extract the solver category either from an enum, or from the newly introduced virtu...
Definition solvercategory.hh:34