3 #ifndef DUNE_ISTL_PRECONDITIONERS_HH 4 #define DUNE_ISTL_PRECONDITIONERS_HH 13 #include <dune/common/unused.hh> 68 template<
class O,
int c = -1>
70 public Preconditioner<typename O::domain_type, typename O::range_type>
89 : inverse_operator_(inverse_operator)
92 DUNE_THROW(InvalidStateException,
"User supplied solver category does not match that of the supplied iverser operator");
95 virtual void pre(domain_type&,range_type&)
98 virtual void apply(domain_type& v,
const range_type& d)
102 inverse_operator_.apply(v, copy, res);
105 virtual void post(domain_type&)
115 InverseOperator& inverse_operator_;
134 template<
class M,
class X,
class Y,
int l=1>
155 SeqSSOR (
const M& A,
int n, scalar_field_type w)
156 : _A_(A), _n(n), _w(w)
166 virtual void pre (X& x, Y& b)
168 DUNE_UNUSED_PARAMETER(x);
169 DUNE_UNUSED_PARAMETER(b);
178 virtual void apply (X& v,
const Y& d)
180 for (
int i=0; i<_n; i++) {
193 DUNE_UNUSED_PARAMETER(x);
208 scalar_field_type _w;
224 template<
class M,
class X,
class Y,
int l=1>
245 SeqSOR (
const M& A,
int n, scalar_field_type w)
246 : _A_(A), _n(n), _w(w)
256 virtual void pre (X& x, Y& b)
258 DUNE_UNUSED_PARAMETER(x);
259 DUNE_UNUSED_PARAMETER(b);
267 virtual void apply (X& v,
const Y& d)
269 this->
template apply<true>(v,d);
280 template<
bool forward>
284 for (
int i=0; i<_n; i++) {
288 for (
int i=0; i<_n; i++) {
300 DUNE_UNUSED_PARAMETER(x);
315 scalar_field_type _w;
329 template<
class M,
class X,
class Y,
int l=1>
350 SeqGS (
const M& A,
int n, scalar_field_type w)
351 : _A_(A), _n(n), _w(w)
361 virtual void pre (X& x, Y& b)
363 DUNE_UNUSED_PARAMETER(x);
364 DUNE_UNUSED_PARAMETER(b);
372 virtual void apply (X& v,
const Y& d)
374 for (
int i=0; i<_n; i++) {
386 DUNE_UNUSED_PARAMETER(x);
401 scalar_field_type _w;
415 template<
class M,
class X,
class Y,
int l=1>
436 SeqJac (
const M& A,
int n, scalar_field_type w)
437 : _A_(A), _n(n), _w(w)
447 virtual void pre (X& x, Y& b)
449 DUNE_UNUSED_PARAMETER(x);
450 DUNE_UNUSED_PARAMETER(b);
458 virtual void apply (X& v,
const Y& d)
460 for (
int i=0; i<_n; i++) {
472 DUNE_UNUSED_PARAMETER(x);
487 scalar_field_type _w;
503 template<
class M,
class X,
class Y,
int l=1>
531 SeqILU (
const M& A, scalar_field_type w,
const bool resort =
false )
532 :
SeqILU( A, 0, w, resort )
544 SeqILU (
const M& A,
int n, scalar_field_type w,
const bool resort =
false )
550 wNotIdentity_(
std::abs( w_ - scalar_field_type(1) ) > 1e-15 )
555 ILU_.reset(
new matrix_type( A ) );
562 ILU_.reset(
new matrix_type( A.N(), A.M(), matrix_type::row_wise) );
580 virtual void pre (X& x, Y& b)
582 DUNE_UNUSED_PARAMETER(x);
583 DUNE_UNUSED_PARAMETER(b);
591 virtual void apply (X& v,
const Y& d)
615 DUNE_UNUSED_PARAMETER(x);
626 std::unique_ptr< matrix_type >
ILU_;
631 std::vector< block_type >
inv_;
634 const scalar_field_type
w_;
651 template<
class M,
class X,
class Y,
int l=1>
684 virtual void pre (X& x, Y& b)
686 DUNE_UNUSED_PARAMETER(x);
687 DUNE_UNUSED_PARAMETER(b);
695 virtual void apply (X& v,
const Y& d)
708 DUNE_UNUSED_PARAMETER(x);
719 scalar_field_type _w;
738 template<
class M,
class X,
class Y,
int l=1>
759 SeqILUn (
const M& A,
int n, scalar_field_type w)
760 : ILU(A.N(),A.M(),M::row_wise),
772 virtual void pre (X& x, Y& b)
774 DUNE_UNUSED_PARAMETER(x);
775 DUNE_UNUSED_PARAMETER(b);
783 virtual void apply (X& v,
const Y& d)
796 DUNE_UNUSED_PARAMETER(x);
811 scalar_field_type _w;
824 template<
class X,
class Y>
850 virtual void pre (X& x, Y& b)
852 DUNE_UNUSED_PARAMETER(x);
853 DUNE_UNUSED_PARAMETER(b);
861 virtual void apply (X& v,
const Y& d)
874 DUNE_UNUSED_PARAMETER(x);
885 scalar_field_type _w;
899 template<
class M,
class X,
class Y >
927 : decomposition_( A.N(), A.M(), matrix_type::random ),
931 for(
auto i = A.begin(), iend = A.end(); i != iend; ++i )
933 const auto &A_i = *i;
934 const auto ij = A_i.find( i.index() );
935 if( ij != A_i.end() )
936 decomposition_.setrowsize( i.index(), ij.offset()+1 );
938 DUNE_THROW(
ISTLError,
"diagonal entry missing" );
940 decomposition_.endrowsizes();
943 for(
auto i = A.begin(), iend = A.end(); i != iend; ++i )
945 const auto &A_i = *i;
946 for(
auto ij = A_i.begin(); ij.index() < i.index() ; ++ij )
947 decomposition_.addindex( i.index(), ij.index() );
948 decomposition_.addindex( i.index(), i.index() );
950 decomposition_.endindices();
954 for(
auto row = decomposition_.begin(), rowend = decomposition_.end(); row != rowend; ++row, ++i )
956 auto ij = i->begin();
957 for(
auto col = row->begin(), colend = row->end();
col != colend; ++
col, ++ij )
966 void pre ( X &x, Y &b )
override 968 DUNE_UNUSED_PARAMETER( x );
969 DUNE_UNUSED_PARAMETER( b );
973 void apply ( X &v,
const Y &d )
override 982 DUNE_UNUSED_PARAMETER( x );
989 matrix_type decomposition_;
990 scalar_field_type relax_;
void apply(X &v, const Y &d)
Apply the preconditioner in a special direction.
Definition: preconditioners.hh:281
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition: preconditioners.hh:783
Sequential SSOR preconditioner.
Definition: preconditioners.hh:135
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:712
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:511
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:361
matrix_type ::block_type block_type
block type of matrix
Definition: preconditioners.hh:509
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:706
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:423
SeqGS(const M &A, int n, scalar_field_type w)
Constructor.
Definition: preconditioners.hh:350
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:298
const bool wNotIdentity_
true if w != 1.0
Definition: preconditioners.hh:636
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:663
void bsorb(const M &A, X &x, const Y &b, const K &w)
SSOR step.
Definition: gsetc.hh:591
SeqJac(const M &A, int n, scalar_field_type w)
Constructor.
Definition: preconditioners.hh:436
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:613
Sequential Gauss Seidel preconditioner.
Definition: preconditioners.hh:330
Definition: allocator.hh:7
void bilu0_decomposition(M &A)
compute ILU decomposition of A. A is overwritten by its decomposition
Definition: ilu.hh:36
SeqILUn(const M &A, int n, scalar_field_type w)
Constructor.
Definition: preconditioners.hh:759
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:470
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:872
void post(X &x) override
Clean up.
Definition: preconditioners.hh:980
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:142
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:427
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:234
virtual void pre(domain_type &, range_type &)
Prepare the preconditioner.
Definition: preconditioners.hh:95
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:191
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:32
Define general, extensible interface for inverse operators.
Richardson preconditioner.
Definition: preconditioners.hh:825
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:830
std::remove_const_t< M > matrix_type
type of matrix the preconditioner is for
Definition: preconditioners.hh:908
std::unique_ptr< matrix_type > ILU_
The ILU(n) decomposition of the matrix. As storage a BCRSMatrix is used.
Definition: preconditioners.hh:626
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
virtual void post(domain_type &)
Clean up.
Definition: preconditioners.hh:105
Sequential ILU preconditioner.
Definition: preconditioners.hh:504
void bilu_decomposition(const M &A, int n, M &ILU)
Definition: ilu.hh:157
void bildl_decompose(Matrix &A)
compute ILDL decomposition of a symmetric matrix A
Definition: ildl.hh:76
void dbgs(const M &A, X &x, const Y &b, const K &w)
GS step.
Definition: gsetc.hh:567
SeqSOR(const M &A, int n, scalar_field_type w)
Constructor.
Definition: preconditioners.hh:245
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:421
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition: preconditioners.hh:861
SeqILU0(const M &A, scalar_field_type w)
Constructor.
Definition: preconditioners.hh:671
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:419
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:878
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:655
compile-time parameter for block recursion depth
Definition: gsetc.hh:40
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:138
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:390
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:197
range_type::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:78
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:828
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:230
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:657
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:966
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:746
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:140
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:304
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:476
ILU::CRS< block_type > CRS
type of ILU storage
Definition: preconditioners.hh:522
O::domain_type domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:74
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:337
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:580
O InverseOperator
type of the wrapped inverse operator
Definition: preconditioners.hh:82
X domain_type
domain type of the preconditioner
Definition: preconditioners.hh:910
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:333
void apply(X &v, const Y &d) override
Apply one step of the preconditioner to the system A(v)=d.
Definition: preconditioners.hh:973
Incomplete LDL decomposition.
virtual void apply(X &v, const Y &d)
Apply the preconditoner.
Definition: preconditioners.hh:591
SeqSSOR(const M &A, int n, scalar_field_type w)
Constructor.
Definition: preconditioners.hh:155
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:109
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:619
const scalar_field_type w_
The relaxation factor to use.
Definition: preconditioners.hh:634
O::range_type range_type
The range type of the preconditioner.
Definition: preconditioners.hh:76
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:519
Richardson(scalar_field_type w=1.0)
Constructor.
Definition: preconditioners.hh:841
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:339
Sequential ILU(n) preconditioner.
Definition: preconditioners.hh:739
SeqILU(const M &A, scalar_field_type w, const bool resort=false)
Constructor.
Definition: preconditioners.hh:531
std::vector< block_type > inv_
Definition: preconditioners.hh:631
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:507
The sequential jacobian preconditioner.
Definition: preconditioners.hh:416
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:794
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:447
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:166
void bildl_backsolve(const Matrix &A, X &v, const Y &d, bool isLowerTriangular=false)
Definition: ildl.hh:137
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:146
CRS upper_
Definition: preconditioners.hh:630
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:832
derive error class from the base class in common
Definition: istlexception.hh:16
X::field_type field_type
field type of the preconditioner
Definition: preconditioners.hh:914
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:916
Turns an InverseOperator into a Preconditioner.
Definition: preconditioners.hh:69
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:98
sequential ILDL preconditioner
Definition: preconditioners.hh:900
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:748
Sequential ILU0 preconditioner.
Definition: preconditioners.hh:652
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:516
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:267
Y range_type
range type of the preconditioner
Definition: preconditioners.hh:912
CRS lower_
The ILU(n) decomposition of the matrix. As storage a CRS structure is used.
Definition: preconditioners.hh:629
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:659
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:335
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:850
void dbjac(const M &A, X &x, const Y &b, const K &w)
Jacobi step.
Definition: gsetc.hh:603
Category
Definition: solvercategory.hh:21
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:228
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:144
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:800
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
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:750
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:772
virtual void apply(X &v, const Y &d)
Apply the preconditoner.
Definition: preconditioners.hh:695
Col col
Definition: matrixmatrix.hh:349
void bilu_backsolve(const CRS &lower, const CRS &upper, const std::vector< mblock > &inv, X &v, const Y &d)
LU backsolve with stored inverse in CRS format for lower and upper triangular.
Definition: ilu.hh:374
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:341
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:458
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:425
InverseOperator2Preconditioner(InverseOperator &inverse_operator)
Construct the preconditioner from the solver.
Definition: preconditioners.hh:88
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:236
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:372
void bilu_backsolve(const M &A, X &v, const Y &d)
LU backsolve with stored inverse.
Definition: ilu.hh:97
void bsorf(const M &A, X &x, const Y &b, const K &w)
SOR step.
Definition: gsetc.hh:579
static void check(const Matrix &mat)
Check whether the a matrix has diagonal values on blocklevel recursion levels.
Definition: matrixutils.hh:88
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:742
SeqILU(const M &A, int n, scalar_field_type w, const bool resort=false)
Constructor.
Definition: preconditioners.hh:544
SeqILDL(const matrix_type &A, scalar_field_type relax=scalar_field_type(1))
constructor
Definition: preconditioners.hh:926
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:684
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:178
Sequential SOR preconditioner.
Definition: preconditioners.hh:225
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:661
Statistics about the application of an inverse operator.
Definition: solver.hh:40
Category for sequential solvers.
Definition: solvercategory.hh:23
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:384
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:232
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:986
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:80
Some handy generic functions for ISTL matrices.
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:744
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:256
SimdScalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:834
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:513