dune-fem 2.8.0
Loading...
Searching...
No Matches
solver/parameter.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_SOLVERPARAMETER_HH
2#define DUNE_FEM_SOLVERPARAMETER_HH
3
5
6namespace Dune
7{
8
9 namespace Fem
10 {
12#ifndef DOXYGEN
13 : public LocalParameter< SolverParameter, SolverParameter >
14#endif
15 {
16 protected:
17 // key prefix, default is fem.solver (can be overloaded by user)
18 const std::string keyPrefix_;
19
21
22 public:
23 // identifier for Fem, ISTL and Petsc solvers
24 static const int cg = 1 ; // CG
25 static const int bicgstab = 2 ; // BiCGStab
26 static const int gmres = 3 ; // GMRES
27 static const int minres = 4 ; // MinRes
28 static const int gradient = 5 ; // GradientSolver
29 static const int loop = 6 ; // LoopSolver
30 static const int superlu = 7 ; // SuperLUSolver
31 static const int bicg = 8 ; // BiCG
32 static const int preonly = 9 ; // only preconder
33 static const std::string solverMethodTable(int i)
34 {
35 std::string methods[] =
36 { "cg", "bicgstab", "gmres", "minres", "gradient", "loop", "superlu", "bicg", "preonly" };
37 return methods[i-1]; // starting to count from 1,...
38 }
39
40 static const int none = 1 ; // no preconditioner
41 static const int ssor = 2 ; // SSOR preconditioner
42 static const int sor = 3 ; // SOR preconditioner
43 static const int ilu = 4 ; // ILU preconditioner
44 static const int gauss_seidel = 5 ; // Gauss-Seidel preconditioner
45 static const int jacobi = 6 ; // Jacobi preconditioner
46 static const int amg_ilu = 7 ; // AMG with ILU-0 smoother (deprecated)
47 static const int amg_jacobi = 8 ; // AMG with Jacobi smoother
48 static const int ildl = 9 ; // ILDL from istl
49 static const int oas = 10; // Overlapping Additive Schwarz
50 static const int icc = 11; // Incomplete Cholesky factorization
51 static const std::string preconditionMethodTable(int i)
52 {
53 static const std::string methods[]
54 = { "none", "ssor", "sor", "ilu", "gauss-seidel", "jacobi", "amg-ilu", "amg-jacobi", "ildl", "oas", "icc" };
55 return methods[i-1]; // starting to count from 1,...
56 }
57
59 : keyPrefix_( "fem.solver." ), parameter_( parameter )
60 {}
61
62 explicit SolverParameter ( const std::string keyPrefix, const ParameterReader &parameter = Parameter::container() )
64 {}
65
66 const std::string& keyPrefix() const { return keyPrefix_; }
67
68 const ParameterReader& parameter() const { return parameter_; }
69
70 virtual void reset()
71 {
72 verbose_ = -1;
73 absoluteTol_ = -1;
74 reductionTol_ = -1;
75 maxIterations_ = -1;
76 }
77
78 virtual bool verbose() const
79 {
80 if( verbose_ < 0 )
81 {
82 verbose_ = parameter_.getValue< bool >( keyPrefix_ + "verbose", false ) ? 1 : 0;
83 }
84 return bool(verbose_);
85 }
86
87 virtual void setVerbose( const bool verb )
88 {
89 verbose_ = verb ? 1 : 0;
90 }
91
92 virtual int errorMeasure() const
93 {
94 const std::string errorTypeTable[] =
95 { "absolute", "relative", "residualreduction" };
96 const int errorType = parameter_.getEnum( keyPrefix_ + "errormeasure", errorTypeTable, 0 );
97 return errorType ;
98 }
99
100 virtual double tolerance ( ) const
101 {
102 if( tolerance_ < 0 )
103 {
104 double defaultTol = 1e-8;
105 if(parameter_.exists(keyPrefix_ + "tolerance"))
106 tolerance_ = parameter_.getValue< double >( keyPrefix_ + "tolerance", defaultTol );
107 else
108 {
109 if( parameter_.exists(keyPrefix_ + "absolutetol") ||
110 parameter_.exists(keyPrefix_ + "reductiontol") )
111 // new parameter not used but old parameters exist
112 {
113 int measure = errorMeasure();
114 if (measure == 0)
115 tolerance_ = absoluteTol__();
116 else
117 tolerance_ = reductionTol__();
118 }
119 else tolerance_ = defaultTol; // no parameter set
120 }
121 }
122 return tolerance_;
123 }
124
125 virtual void setTolerance ( const double eps )
126 {
127 assert( eps >= 0.0 );
128 tolerance_ = eps;
129 }
130
131 virtual int maxIterations () const
132 {
133 if( maxIterations_ < 0 )
134 {
135 if(parameter_.exists(keyPrefix_ + "maxlineariterations"))
136 {
137 std::cout << "WARNING: Parameter " + keyPrefix_ +"maxlineariterations is deprecated. Please use " + keyPrefix_ + "maxiterations instead." << std::endl;
138 maxIterations_ = parameter_.getValue< double >(keyPrefix_ + "maxlineariterations");
139 }
140 else
141 maxIterations_ = parameter_.getValue< int >( keyPrefix_ + "maxiterations", std::numeric_limits< int >::max() );
142 }
143 return maxIterations_;
144 }
145
146 virtual void setMaxIterations ( const int maxIter )
147 {
148 assert( maxIter >= 0 );
149 maxIterations_ = maxIter;
150 }
151
152 virtual int solverMethod(
153 const std::vector<int> standardMethods,
154 const std::vector<std::string> &additionalMethods = {},
155 int defaultMethod = 0 // this is the first method passed in
156 ) const
157 {
158 std::vector<std::string> methodTable(standardMethods.size()+additionalMethods.size());
159 for (std::size_t i=0;i<standardMethods.size();++i)
160 methodTable[i] = solverMethodTable(standardMethods[i]);
161 for (std::size_t i=0;i<additionalMethods.size();++i)
162 methodTable[standardMethods.size()+i] = additionalMethods[i];
163 std::size_t method;
164 if( parameter_.exists( keyPrefix_ + "method" ) ||
165 !parameter_.exists( "method" ) )
166 method = parameter_.getEnum( keyPrefix_ + "method", methodTable, defaultMethod );
167 else
168 {
169 method = parameter_.getEnum( "krylovmethod", methodTable, defaultMethod );
170 std::cout << "WARNING: using old parameter name 'krylovmethod' "
171 << "please switch to '" << keyPrefix_ << "method'\n";
172 }
173 if (method < standardMethods.size())
174 return standardMethods[method];
175 else
176 return -(method-standardMethods.size()); // return in [ 0,-additionalMethods.size() )
177 }
178
179 virtual int gmresRestart() const
180 {
181 int defaultRestart = 20;
182 return parameter_.getValue< int >( keyPrefix_ + "gmres.restart", defaultRestart );
183 }
184
186 const std::vector<int> standardMethods,
187 const std::vector<std::string> &additionalMethods = {},
188 int defaultMethod = 0 // this is the first method passed in
189 ) const
190 {
191 std::vector<std::string> methodTable(standardMethods.size()+additionalMethods.size());
192 for (std::size_t i=0;i<standardMethods.size();++i)
193 methodTable[i] = preconditionMethodTable(standardMethods[i]);
194 for (std::size_t i=0;i<additionalMethods.size();++i)
195 methodTable[standardMethods.size()+i] = additionalMethods[i];
196 std::size_t method = parameter_.getEnum( keyPrefix_ + "preconditioning.method", methodTable, defaultMethod );
197 if (method < standardMethods.size())
198 return standardMethods[method];
199 else
200 return -(method-standardMethods.size()); // return in [ 0,-additionalMethods.size() )
201 }
202
203 virtual double relaxation () const
204 {
205 return parameter_.getValue< double >( keyPrefix_ + "preconditioning.relaxation", 1.1 );
206 }
207
208 virtual int preconditionerIteration () const
209 {
210 return parameter_.getValue< int >( keyPrefix_ + "preconditioning.iterations", 1 );
211 }
212
213 virtual int preconditionerLevel () const
214 {
215 return parameter_.getValue< int >( keyPrefix_ + "preconditioning.level", 0 );
216 }
217
218 private:
219 virtual double absoluteTol__ ( ) const
220 {
221 if( absoluteTol_ < 0 )
222 {
223 if(parameter_.exists(keyPrefix_ + "linabstol"))
224 {
225 std::cout << "WARNING: Parameter " + keyPrefix_ + "linabstol is deprecated. Please use " + keyPrefix_ + "absolutetol instead." << std::endl;
226 absoluteTol_ = parameter_.getValue< double >(keyPrefix_ + "linabstol");
227 }
228 else
229 {
230 std::cout << "WARNING: Parameter " + keyPrefix_ + "absolutetol is deprecated. Please use " + keyPrefix_ + "tolerance instead." << std::endl;
231 absoluteTol_ = parameter_.getValue< double >(keyPrefix_ + "absolutetol", 1e-8 );
232 }
233 }
234 return absoluteTol_;
235 }
236
237 virtual double reductionTol__ ( ) const
238 {
239 if( reductionTol_ < 0 )
240 {
241 if(parameter_.exists(keyPrefix_ + "linreduction"))
242 {
243 std::cout << "WARNING: Parameter " + keyPrefix_ +"linreduction is deprecated. Please use " + keyPrefix_ + "reductiontol instead." << std::endl;
244 reductionTol_ = parameter_.getValue< double >(keyPrefix_ + "linreduction");
245 }
246 else
247 {
248 std::cout << "WARNING: Parameter " + keyPrefix_ +"reductiontol is deprecated. Please use " + keyPrefix_ + "tolerance instead." << std::endl;
249 reductionTol_ = parameter_.getValue< double >( keyPrefix_ + "reductiontol", 1e-8 );
250 }
251 }
252 return reductionTol_;
253 }
254
255 mutable int verbose_ = -1;
256 mutable int maxIterations_ = -1;
257 mutable double absoluteTol_ = -1.;
258 mutable double reductionTol_ = -1.;
259 mutable double tolerance_ = -1.;
260 };
261 }
262}
263
264#endif // #ifndef DUNE_FEM_SOLVERPARAMETER_HH
Definition: bindguard.hh:11
static ParameterContainer & container()
Definition: io/parameter.hh:193
Definition: io/parameter.hh:551
int getEnum(const std::string &key, const std::string(&values)[n]) const
Definition: reader.hh:225
T getValue(const std::string &key) const
get mandatory parameter
Definition: reader.hh:159
bool exists(const std::string &key) const
check, whether a parameter is defined
Definition: reader.hh:44
Definition: solver/parameter.hh:15
virtual int errorMeasure() const
Definition: solver/parameter.hh:92
static const int loop
Definition: solver/parameter.hh:29
static const std::string solverMethodTable(int i)
Definition: solver/parameter.hh:33
static const int none
Definition: solver/parameter.hh:40
static const int cg
Definition: solver/parameter.hh:24
static const int icc
Definition: solver/parameter.hh:50
static const std::string preconditionMethodTable(int i)
Definition: solver/parameter.hh:51
virtual double tolerance() const
Definition: solver/parameter.hh:100
virtual void setVerbose(const bool verb)
Definition: solver/parameter.hh:87
static const int gmres
Definition: solver/parameter.hh:26
virtual void setMaxIterations(const int maxIter)
Definition: solver/parameter.hh:146
SolverParameter(const ParameterReader &parameter=Parameter::container())
Definition: solver/parameter.hh:58
static const int bicg
Definition: solver/parameter.hh:31
virtual int preconditionMethod(const std::vector< int > standardMethods, const std::vector< std::string > &additionalMethods={}, int defaultMethod=0) const
Definition: solver/parameter.hh:185
virtual int preconditionerLevel() const
Definition: solver/parameter.hh:213
static const int ssor
Definition: solver/parameter.hh:41
SolverParameter(const std::string keyPrefix, const ParameterReader &parameter=Parameter::container())
Definition: solver/parameter.hh:62
virtual int gmresRestart() const
Definition: solver/parameter.hh:179
virtual void setTolerance(const double eps)
Definition: solver/parameter.hh:125
virtual double relaxation() const
Definition: solver/parameter.hh:203
static const int sor
Definition: solver/parameter.hh:42
static const int ildl
Definition: solver/parameter.hh:48
static const int minres
Definition: solver/parameter.hh:27
static const int gradient
Definition: solver/parameter.hh:28
static const int gauss_seidel
Definition: solver/parameter.hh:44
virtual int preconditionerIteration() const
Definition: solver/parameter.hh:208
static const int ilu
Definition: solver/parameter.hh:43
virtual int solverMethod(const std::vector< int > standardMethods, const std::vector< std::string > &additionalMethods={}, int defaultMethod=0) const
Definition: solver/parameter.hh:152
static const int bicgstab
Definition: solver/parameter.hh:25
static const int superlu
Definition: solver/parameter.hh:30
const ParameterReader & parameter() const
Definition: solver/parameter.hh:68
ParameterReader parameter_
Definition: solver/parameter.hh:20
static const int jacobi
Definition: solver/parameter.hh:45
static const int amg_ilu
Definition: solver/parameter.hh:46
virtual int maxIterations() const
Definition: solver/parameter.hh:131
virtual bool verbose() const
Definition: solver/parameter.hh:78
virtual void reset()
Definition: solver/parameter.hh:70
const std::string keyPrefix_
Definition: solver/parameter.hh:18
static const int amg_jacobi
Definition: solver/parameter.hh:47
const std::string & keyPrefix() const
Definition: solver/parameter.hh:66
static const int preonly
Definition: solver/parameter.hh:32
static const int oas
Definition: solver/parameter.hh:49