DUNE PDELab (git)

preconditioners.hh
Go to the documentation of this file.
1// SPDX-FileCopyrightText: Copyright © DUNE Project contributors, see file LICENSE.md in module root
2// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4// vi: set et ts=4 sw=2 sts=2:
5#ifndef DUNE_ISTL_PRECONDITIONERS_HH
6#define DUNE_ISTL_PRECONDITIONERS_HH
7
8#include <cmath>
9#include <complex>
10#include <iostream>
11#include <iomanip>
12#include <memory>
13#include <string>
14
17
18#include <dune/istl/solverregistry.hh>
19#include "preconditioner.hh"
20#include "solver.hh"
21#include "solvercategory.hh"
22#include "istlexception.hh"
23#include "matrixutils.hh"
24#include "gsetc.hh"
25#include "dilu.hh"
26#include "ildl.hh"
27#include "ilu.hh"
28
29
30namespace Dune {
73 template<class O, int c = -1>
75 public Preconditioner<typename O::domain_type, typename O::range_type>
76 {
77 public:
79 typedef typename O::domain_type domain_type;
81 typedef typename O::range_type range_type;
83 typedef typename range_type::field_type field_type;
87 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
89 typedef O InverseOperator;
90
96 : inverse_operator_(inverse_operator)
97 {
98 if(c != -1 && SolverCategory::category(inverse_operator_) != c)
99 DUNE_THROW(InvalidStateException, "User-supplied solver category does not match that of the given inverse operator");
100 }
101
102 void pre([[maybe_unused]] domain_type&,[[maybe_unused]] range_type&) override
103 {}
104
105 void apply(domain_type& v, const range_type& d) override
106 {
108 range_type copy(d);
109 inverse_operator_.apply(v, copy, res);
110 }
111
112 void post([[maybe_unused]] domain_type&) override
113 {}
114
117 {
118 return SolverCategory::category(inverse_operator_);
119 }
120
121 private:
122 InverseOperator& inverse_operator_;
123 };
124
125 //=====================================================================
126 // Implementation of this interface for sequential ISTL-preconditioners
127 //=====================================================================
128
129
141 template<class M, class X, class Y, int l=1>
142 class SeqSSOR : public Preconditioner<X,Y> {
143 public:
145 typedef M matrix_type;
147 typedef X domain_type;
149 typedef Y range_type;
151 typedef typename X::field_type field_type;
155 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
156
164 SeqSSOR (const M& A, int n, real_field_type w)
165 : _A_(A), _n(n), _w(w)
166 {
168 }
169
183 SeqSSOR (const std::shared_ptr<const AssembledLinearOperator<M,X,Y>>& A, const ParameterTree& configuration)
184 : SeqSSOR(A->getmat(), configuration)
185 {}
186
200 SeqSSOR (const M& A, const ParameterTree& configuration)
201 : SeqSSOR(A, configuration.get<int>("iterations",1), configuration.get<real_field_type>("relaxation",1.0))
202 {}
203
209 void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b) override
210 {}
211
217 void apply (X& v, const Y& d) override
218 {
219 for (int i=0; i<_n; i++) {
220 bsorf(_A_,v,d,_w,BL<l>());
221 bsorb(_A_,v,d,_w,BL<l>());
222 }
223 }
224
230 void post ([[maybe_unused]] X& x) override
231 {}
232
235 {
237 }
238
239 private:
241 const M& _A_;
243 int _n;
246 };
247 DUNE_REGISTER_PRECONDITIONER("ssor", defaultPreconditionerBlockLevelCreator<Dune::SeqSSOR>());
248
249
261 template<class M, class X, class Y, int l=1>
262 class SeqSOR : public Preconditioner<X,Y> {
263 public:
265 typedef M matrix_type;
267 typedef X domain_type;
269 typedef Y range_type;
271 typedef typename X::field_type field_type;
275 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
276
284 SeqSOR (const M& A, int n, real_field_type w)
285 : _A_(A), _n(n), _w(w)
286 {
288 }
289
303 SeqSOR (const std::shared_ptr<const AssembledLinearOperator<M,X,Y>>& A, const ParameterTree& configuration)
304 : SeqSOR(A->getmat(), configuration)
305 {}
306
320 SeqSOR (const M& A, const ParameterTree& configuration)
321 : SeqSOR(A, configuration.get<int>("iterations",1), configuration.get<real_field_type>("relaxation",1.0))
322 {}
323
329 void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b) override
330 {}
331
337 void apply (X& v, const Y& d) override
338 {
339 this->template apply<true>(v,d);
340 }
341
350 template<bool forward>
351 void apply(X& v, const Y& d)
352 {
353 if(forward)
354 for (int i=0; i<_n; i++) {
355 bsorf(_A_,v,d,_w,BL<l>());
356 }
357 else
358 for (int i=0; i<_n; i++) {
359 bsorb(_A_,v,d,_w,BL<l>());
360 }
361 }
362
368 void post ([[maybe_unused]] X& x) override
369 {}
370
373 {
375 }
376
377 private:
379 const M& _A_;
381 int _n;
384 };
385 DUNE_REGISTER_PRECONDITIONER("sor", defaultPreconditionerBlockLevelCreator<Dune::SeqSOR>());
386
387
398 template<class M, class X, class Y, int l=1>
400 DUNE_REGISTER_PRECONDITIONER("gs", defaultPreconditionerBlockLevelCreator<Dune::SeqGS>());
401
412 template<class M, class X, class Y, int l=1>
414 public:
416 typedef M matrix_type;
418 typedef X domain_type;
420 typedef Y range_type;
422 typedef typename X::field_type field_type;
426 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
427
435 SeqJac (const M& A, int n, real_field_type w)
436 : _A_(A), _n(n), _w(w)
437 {
439 }
440
454 SeqJac (const std::shared_ptr<const AssembledLinearOperator<M,X,Y>>& A, const ParameterTree& configuration)
455 : SeqJac(A->getmat(), configuration)
456 {}
457
471 SeqJac (const M& A, const ParameterTree& configuration)
472 : SeqJac(A, configuration.get<int>("iterations",1), configuration.get<real_field_type>("relaxation",1.0))
473 {}
474
480 void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b) override
481 {}
482
488 void apply (X& v, const Y& d) override
489 {
490 for (int i=0; i<_n; i++) {
491 dbjac(_A_,v,d,_w,BL<l>());
492 }
493 }
494
500 void post ([[maybe_unused]] X& x) override
501 {}
502
505 {
507 }
508
509 private:
511 const M& _A_;
513 int _n;
515 real_field_type _w;
516 };
517 DUNE_REGISTER_PRECONDITIONER("jac", defaultPreconditionerBlockLevelCreator<Dune::SeqJac>());
518
562 template <class M, class X, class Y, int l = 1>
564 {
565 public:
567 using matrix_type = M;
569 using block_type = typename matrix_type::block_type;
571 using domain_type = X;
573 using range_type = Y;
574
576 using field_type = typename X::field_type;
577
581 using real_field_type = typename FieldTraits<scalar_field_type>::real_type;
582
590 : _A_(A),
591 _w(w),
592 wNotIdentity_([w]
593 {using std::abs; return abs(w - real_field_type(1)) > 1e-15; }())
594 {
595 Dinv_.resize(_A_.N());
597 DILU::blockDILUDecomposition(_A_, Dinv_);
598 }
599
612 SeqDILU(const std::shared_ptr<const AssembledLinearOperator<M, X, Y>> &A, const ParameterTree &configuration)
613 : SeqDILU(A->getmat(), configuration)
614 {
615 }
616
629 SeqDILU(const M &A, const ParameterTree &config)
630 : SeqDILU(A, config.get<real_field_type>("relaxation", 1.0))
631 {
632 }
633
639 void pre([[maybe_unused]] X &x, [[maybe_unused]] Y &b) override
640 {
641 }
642
648 void apply(X &v, const Y &d) override
649 {
650
651 DILU::blockDILUBacksolve(_A_, Dinv_, v, d);
652
653 if (wNotIdentity_)
654 {
655 v *= _w;
656 }
657 }
658
664 void post([[maybe_unused]] X &x) override
665 {
666 }
667
670 {
672 }
673
674 protected:
675 std::vector<block_type> Dinv_;
677 const M &_A_;
681 const bool wNotIdentity_;
682 };
683 DUNE_REGISTER_PRECONDITIONER("dilu", defaultPreconditionerBlockLevelCreator<Dune::SeqDILU>());
684
696 template<class M, class X, class Y, int l=1>
697 class SeqILU : public Preconditioner<X,Y> {
698 public:
700 typedef typename std::remove_const<M>::type matrix_type;
702 typedef typename matrix_type :: block_type block_type;
704 typedef X domain_type;
706 typedef Y range_type;
707
709 typedef typename X::field_type field_type;
710
714 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
715
718
726 SeqILU (const M& A, real_field_type w, const bool resort = false )
727 : SeqILU( A, 0, w, resort ) // construct ILU(0)
728 {
729 }
730
745 SeqILU (const std::shared_ptr<const AssembledLinearOperator<M,X,Y>>& A, const ParameterTree& configuration)
746 : SeqILU(A->getmat(), configuration)
747 {}
748
763 SeqILU(const M& A, const ParameterTree& config)
764 : SeqILU(A, config.get("n", 0),
765 config.get<real_field_type>("relaxation", 1.0),
766 config.get("resort", false))
767 {}
768
777 SeqILU (const M& A, int n, real_field_type w, const bool resort = false )
778 : ILU_(),
779 lower_(),
780 upper_(),
781 inv_(),
782 w_(w),
783 wNotIdentity_([w]{using std::abs; return abs(w - real_field_type(1)) > 1e-15;}() )
784 {
785 if( n == 0 )
786 {
787 // copy A
788 ILU_.reset( new matrix_type( A ) );
789 // create ILU(0) decomposition
790 ILU::blockILU0Decomposition( *ILU_ );
791 }
792 else
793 {
794 // create matrix in build mode
795 ILU_.reset( new matrix_type( A.N(), A.M(), matrix_type::row_wise) );
796 // create ILU(n) decomposition
797 ILU::blockILUDecomposition( A, n, *ILU_ );
798 }
799
800 if( resort )
801 {
802 // store ILU in simple CRS format
803 ILU::convertToCRS( *ILU_, lower_, upper_, inv_ );
804 ILU_.reset();
805 }
806 }
807
813 void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b) override
814 {}
815
821 void apply (X& v, const Y& d) override
822 {
823 if( ILU_ )
824 {
825 ILU::blockILUBacksolve( *ILU_, v, d);
826 }
827 else
828 {
829 ILU::blockILUBacksolve(lower_, upper_, inv_, v, d);
830 }
831
832 if( wNotIdentity_ )
833 {
834 v *= w_;
835 }
836 }
837
843 void post ([[maybe_unused]] X& x) override
844 {}
845
848 {
850 }
851
852 protected:
854 std::unique_ptr< matrix_type > ILU_;
855
858 CRS upper_;
859 std::vector< block_type, typename matrix_type::allocator_type > inv_;
860
864 const bool wNotIdentity_;
865 };
866 DUNE_REGISTER_PRECONDITIONER("ilu", defaultPreconditionerBlockLevelCreator<Dune::SeqILU>());
867
868
877 template<class X, class Y>
878 class Richardson : public Preconditioner<X,Y> {
879 public:
881 typedef X domain_type;
883 typedef Y range_type;
885 typedef typename X::field_type field_type;
889 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
890
897 _w(w)
898 {}
899
911 Richardson (const ParameterTree& configuration)
912 : Richardson(configuration.get<real_field_type>("relaxation", 1.0))
913 {}
914
920 void pre ([[maybe_unused]] X& x, [[maybe_unused]] Y& b) override
921 {}
922
928 void apply (X& v, const Y& d) override
929 {
930 v = d;
931 v *= _w;
932 }
933
939 void post ([[maybe_unused]] X& x) override
940 {}
941
944 {
946 }
947
948 private:
951 };
952 DUNE_REGISTER_PRECONDITIONER("richardson", [](auto opTraits, const auto& op, const ParameterTree& config){
953 using OpTraits = std::decay_t<decltype(opTraits)>;
954 using D = typename OpTraits::domain_type;
955 using R = typename OpTraits::range_type;
956 return std::make_shared<Richardson<D,R>>(config);
957 });
958
959
970 template< class M, class X, class Y >
972 : public Preconditioner< X, Y >
973 {
974 typedef SeqILDL< M, X, Y > This;
976
977 public:
979 typedef std::remove_const_t< M > matrix_type;
981 typedef X domain_type;
983 typedef Y range_type;
985 typedef typename X::field_type field_type;
989 typedef typename FieldTraits<scalar_field_type>::real_type real_field_type;
990
1003 SeqILDL (const std::shared_ptr<const AssembledLinearOperator<M,X,Y>>& A, const ParameterTree& configuration)
1004 : SeqILDL(A->getmat(), configuration)
1005 {}
1006
1019 SeqILDL(const matrix_type& A, const ParameterTree& config)
1020 : SeqILDL(A, config.get<real_field_type>("relaxation", 1.0))
1021 {}
1022
1031 explicit SeqILDL ( const matrix_type &A, real_field_type relax = real_field_type( 1 ) )
1032 : decomposition_( A.N(), A.M(), matrix_type::random ),
1033 relax_( relax )
1034 {
1035 // setup row sizes for lower triangular matrix
1036 for( auto i = A.begin(), iend = A.end(); i != iend; ++i )
1037 {
1038 const auto &A_i = *i;
1039 const auto ij = A_i.find( i.index() );
1040 if( ij != A_i.end() )
1041 decomposition_.setrowsize( i.index(), ij.offset()+1 );
1042 else
1043 DUNE_THROW( ISTLError, "diagonal entry missing" );
1044 }
1045 decomposition_.endrowsizes();
1046
1047 // setup row indices for lower triangular matrix
1048 for( auto i = A.begin(), iend = A.end(); i != iend; ++i )
1049 {
1050 const auto &A_i = *i;
1051 for( auto ij = A_i.begin(); ij.index() < i.index() ; ++ij )
1052 decomposition_.addindex( i.index(), ij.index() );
1053 decomposition_.addindex( i.index(), i.index() );
1054 }
1055 decomposition_.endindices();
1056
1057 // copy values of lower triangular matrix
1058 auto i = A.begin();
1059 for( auto row = decomposition_.begin(), rowend = decomposition_.end(); row != rowend; ++row, ++i )
1060 {
1061 auto ij = i->begin();
1062 for( auto col = row->begin(), colend = row->end(); col != colend; ++col, ++ij )
1063 *col = *ij;
1064 }
1065
1066 // perform ILDL decomposition
1067 bildl_decompose( decomposition_ );
1068 }
1069
1071 void pre ([[maybe_unused]] X &x, [[maybe_unused]] Y &b) override
1072 {}
1073
1075 void apply ( X &v, const Y &d ) override
1076 {
1077 bildl_backsolve( decomposition_, v, d, true );
1078 v *= relax_;
1079 }
1080
1082 void post ([[maybe_unused]] X &x) override
1083 {}
1084
1087
1088 private:
1089 matrix_type decomposition_;
1090 real_field_type relax_;
1091 };
1092 DUNE_REGISTER_PRECONDITIONER("ildl", defaultPreconditionerCreator<Dune::SeqILDL>());
1093
1096} // end namespace
1097
1098
1099#endif
A linear operator exporting itself in matrix form.
Definition: operators.hh:110
derive error class from the base class in common
Definition: istlexception.hh:19
Default exception if a function was called while the object is not in a valid state for that function...
Definition: exceptions.hh:281
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
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:116
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
void post(domain_type &) override
Clean up.
Definition: preconditioners.hh:112
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
void pre(domain_type &, range_type &) override
Prepare the preconditioner.
Definition: preconditioners.hh:102
void apply(domain_type &v, const range_type &d) override
Apply one step of the preconditioner to the system A(v)=d.
Definition: preconditioners.hh:105
Abstract base class for all solvers.
Definition: solver.hh:101
Hierarchical structure of string parameters.
Definition: parametertree.hh:37
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:33
Y range_type
The range type of the preconditioner.
Definition: preconditioner.hh:38
X domain_type
The domain type of the preconditioner.
Definition: preconditioner.hh:36
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioner.hh:40
Richardson preconditioner.
Definition: preconditioners.hh:878
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:885
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:883
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:920
Richardson(real_field_type w=1.0)
Constructor.
Definition: preconditioners.hh:896
void apply(X &v, const Y &d) override
Apply the precondioner.
Definition: preconditioners.hh:928
void post(X &x) override
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
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:943
Sequential DILU preconditioner.
Definition: preconditioners.hh:564
const bool wNotIdentity_
true if w != 1.0
Definition: preconditioners.hh:681
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
typename FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition: preconditioners.hh:581
void apply(X &v, const Y &d) override
Apply the preconditioner.
Definition: preconditioners.hh:648
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:567
typename matrix_type::block_type block_type
block type of matrix
Definition: preconditioners.hh:569
SeqDILU(const M &A, real_field_type w)
Constructor.
Definition: preconditioners.hh:589
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:639
SeqDILU(const M &A, const ParameterTree &config)
Constructor.
Definition: preconditioners.hh:629
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:669
void post(X &x) override
Clean up.
Definition: preconditioners.hh:664
sequential ILDL preconditioner
Definition: preconditioners.hh:973
SeqILDL(const matrix_type &A, const ParameterTree &config)
Constructor.
Definition: preconditioners.hh:1019
SeqILDL(const matrix_type &A, real_field_type relax=real_field_type(1))
constructor
Definition: preconditioners.hh:1031
X domain_type
domain type of the preconditioner
Definition: preconditioners.hh:981
void post(X &x) override
Clean up.
Definition: preconditioners.hh:1082
Y range_type
range type of the preconditioner
Definition: preconditioners.hh:983
std::remove_const_t< M > matrix_type
type of matrix the preconditioner is for
Definition: preconditioners.hh:979
void apply(X &v, const Y &d) override
Apply one step of the preconditioner to the system A(v)=d.
Definition: preconditioners.hh:1075
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition: preconditioners.hh:989
SeqILDL(const std::shared_ptr< const AssembledLinearOperator< M, X, Y > > &A, const ParameterTree &configuration)
Constructor.
Definition: preconditioners.hh:1003
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:1071
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:987
X::field_type field_type
field type of the preconditioner
Definition: preconditioners.hh:985
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:1086
Sequential ILU preconditioner.
Definition: preconditioners.hh:697
SeqILU(const M &A, int n, real_field_type w, const bool resort=false)
Constructor.
Definition: preconditioners.hh:777
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
void post(X &x) override
Clean up.
Definition: preconditioners.hh:843
const bool wNotIdentity_
true if w != 1.0
Definition: preconditioners.hh:864
SeqILU(const M &A, const ParameterTree &config)
Constructor.
Definition: preconditioners.hh:763
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:813
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
SeqILU(const M &A, real_field_type w, const bool resort=false)
Constructor.
Definition: preconditioners.hh:726
void apply(X &v, const Y &d) override
Apply the preconditioner.
Definition: preconditioners.hh:821
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:847
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
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
The sequential jacobian preconditioner.
Definition: preconditioners.hh:413
SeqJac(const M &A, const ParameterTree &configuration)
Constructor.
Definition: preconditioners.hh:471
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
void apply(X &v, const Y &d) override
Apply the preconditioner.
Definition: preconditioners.hh:488
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:422
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:480
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:418
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:504
void post(X &x) override
Clean up.
Definition: preconditioners.hh:500
FieldTraits< scalar_field_type >::real_type real_field_type
real scalar type underlying the field_type
Definition: preconditioners.hh:426
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 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
void post(X &x) override
Clean up.
Definition: preconditioners.hh:368
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
Simd::Scalar< field_type > scalar_field_type
scalar type underlying the field_type
Definition: preconditioners.hh:273
void apply(X &v, const Y &d) override
Apply the preconditioner.
Definition: preconditioners.hh:337
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:329
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
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:372
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
Sequential SSOR preconditioner.
Definition: preconditioners.hh:142
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition: preconditioners.hh:209
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
SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition: preconditioners.hh:234
void post(X &x) override
Clean up.
Definition: preconditioners.hh:230
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
void apply(X &v, const Y &d) override
Apply the preconditioner.
Definition: preconditioners.hh:217
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
The diagonal incomplete LU factorization kernels.
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
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
typename Overloads::ScalarType< std::decay_t< V > >::type Scalar
Element type of some SIMD type.
Definition: interface.hh:235
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
Incomplete LDL decomposition.
The incomplete LU factorization kernels.
Some handy generic functions for ISTL matrices.
Dune namespace.
Definition: alignedallocator.hh:13
void bildl_decompose(Matrix &A)
compute ILDL decomposition of a symmetric matrix A
Definition: ildl.hh:88
constexpr auto get(std::integer_sequence< T, II... >, std::integral_constant< std::size_t, pos >={})
Return the entry at position pos of the given sequence.
Definition: integersequence.hh:22
A hierarchical structure of string parameters.
Include file for users of the SIMD abstraction layer.
compile-time parameter for block recursion depth
Definition: gsetc.hh:45
static void check(const Matrix &mat)
Check whether the a matrix has diagonal values on blocklevel recursion levels.
Definition: matrixutils.hh:53
Statistics about the application of an inverse operator.
Definition: solver.hh:50
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
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)