All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
weightedresidreductioncriterion.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef EWOMS_ISTL_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH
28 #define EWOMS_ISTL_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH
29 
30 #include "convergencecriterion.hh"
31 
32 #include <iostream>
33 
34 namespace Ewoms {
35 namespace Linear {
36 
57 template <class Vector, class CollectiveCommunication>
59 {
60  typedef typename Vector::field_type Scalar;
61  typedef typename Vector::block_type BlockType;
62 
63 public:
64  WeightedResidualReductionCriterion(const CollectiveCommunication& comm)
65  : comm_(comm)
66  {}
67 
68  WeightedResidualReductionCriterion(const CollectiveCommunication& comm,
69  const Vector& residWeights,
71  Scalar fixPointTolerance,
72  Scalar absResidualTolerance = 0.0,
73  Scalar maxError = 0.0)
74  : comm_(comm),
75  residWeightVec_(residWeights),
76  fixPointTolerance_(fixPointTolerance),
77  residualReductionTolerance_(residualReductionTolerance),
78  absResidualTolerance_(absResidualTolerance),
79  maxError_(maxError)
80  { }
81 
98  void setResidualWeight(const Vector& residWeightVec)
99  { residWeightVec_ = residWeightVec; }
100 
107  Scalar residualWeight(size_t outerIdx, unsigned innerIdx) const
108  {
109  return (residWeightVec_.size() == 0)
110  ? 1.0
111  : residWeightVec_[outerIdx][innerIdx];
112  }
113 
118  { residualReductionTolerance_ = tol; }
119 
124  { return residualReductionTolerance_; }
125 
129  void setResidualTolerance(Scalar tol)
130  { absResidualTolerance_ = tol; }
131 
135  Scalar absResidualTolerance() const
136  { return absResidualTolerance_; }
137 
142  Scalar residualAccuracy() const
143  { return residualError_/std::max<Scalar>(1e-20, initialResidualError_); }
144 
148  void setFixPointTolerance(Scalar tol)
149  { fixPointTolerance_ = tol; }
150 
156  Scalar fixPointTolerance() const
157  { return fixPointTolerance_; }
158 
163  Scalar fixPointAccuracy() const
164  { return fixPointError_; }
165 
169  void setInitial(const Vector& curSol, const Vector& curResid)
170  {
171  lastResidualError_ = 1e100;
172 
173  lastSolVec_ = curSol;
174  updateErrors_(curSol, curResid);
175  // the fix-point error is not applicable for the initial solution!
176  fixPointError_ = 1e100;
177 
178  // make sure that we don't allow an initial error of 0 to avoid
179  // divisions by zero
180  residualError_ = std::max<Scalar>(residualError_, 1e-20);
181  initialResidualError_ = residualError_;
182  }
183 
187  void update(const Vector& curSol, const Vector& curResid)
188  {
189  lastResidualError_ = residualError_;
190  updateErrors_(curSol, curResid);
191  }
192 
196  bool converged() const
197  {
198  // we're converged if the solution is better than the tolerance
199  // fix-point and residual tolerance.
200  return
202  residualError_ <= absResidualTolerance_;
203  }
204 
208  bool failed() const
209  {
210  return
211  (!converged() && fixPointError_ <= fixPointTolerance_)
212  || residualError_ > maxError_;
213  }
214 
220  Scalar accuracy() const
221  { return residualError_; }
222 
226  void printInitial(std::ostream& os = std::cout) const
227  {
228  os << std::setw(20) << " Iter ";
229  os << std::setw(20) << " Delta ";
230  os << std::setw(20) << " Residual ";
231  os << std::setw(20) << " ResidRed ";
232  os << std::setw(20) << " Rate ";
233  os << std::endl;
234  }
235 
239  void print(Scalar iter, std::ostream& os = std::cout) const
240  {
241  static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1e10;
242 
243  os << std::setw(20) << iter << " ";
244  os << std::setw(20) << fixPointAccuracy() << " ";
245  os << std::setw(20) << residualError_ << " ";
246  os << std::setw(20) << 1.0/residualAccuracy() << " ";
247  os << std::setw(20) << lastResidualError_ / std::max<Scalar>(residualError_, eps) << " ";
248  os << std::endl << std::flush;
249  }
250 
251 private:
252  // update the weighted absolute residual
253  void updateErrors_(const Vector& curSol, const Vector& curResid)
254  {
255  residualError_ = 0.0;
256  fixPointError_ = 0.0;
257  for (size_t i = 0; i < curResid.size(); ++i) {
258  for (unsigned j = 0; j < BlockType::dimension; ++j) {
259  residualError_ =
260  std::max<Scalar>(residualError_,
261  residualWeight(i, j)*std::abs(curResid[i][j]));
262  fixPointError_ =
263  std::max<Scalar>(fixPointError_,
264  std::abs(curSol[i][j] - lastSolVec_[i][j])
265  /std::max<Scalar>(1.0, curSol[i][j]));
266  }
267  }
268  lastSolVec_ = curSol;
269 
270  residualError_ = comm_.max(residualError_);
271  fixPointError_ = comm_.max(fixPointError_);
272  }
273 
274  const CollectiveCommunication& comm_;
275 
276  // the weights of the components of the residual
277  Vector residWeightVec_;
278 
279  // solution vector of the last iteration
280  Vector lastSolVec_;
281 
282  // the maximum of the weighted difference between the last two
283  // iterations
284  Scalar fixPointError_;
285 
286  // the maximum allowed relative tolerance for difference of the
287  // solution of two iterations
288  Scalar fixPointTolerance_;
289 
290  // the maximum of the absolute weighted residual of the last
291  // iteration
292  Scalar residualError_;
293 
294  // the maximum of the absolute weighted difference of the last
295  // iteration
296  Scalar lastResidualError_;
297 
298  // the maximum of the absolute weighted residual of the initial
299  // solution
300  Scalar initialResidualError_;
301 
302  // the maximum allowed relative tolerance of the residual for the
303  // solution to be considered converged
304  Scalar residualReductionTolerance_;
305 
306  // the maximum allowed absolute tolerance of the residual for the
307  // solution to be considered converged
308  Scalar absResidualTolerance_;
309 
310  // The maximum error which is tolerated before we fail.
311  Scalar maxError_;
312 };
313 
315 
316 }} // end namespace Linear,Ewoms
317 
318 #endif
Base class for all convergence criteria which only defines an virtual API.
Definition: convergencecriterion.hh:57
bool converged() const
Returns true if and only if the convergence criterion is met.
Definition: weightedresidreductioncriterion.hh:196
void setFixPointTolerance(Scalar tol)
Sets the fix-point tolerance.
Definition: weightedresidreductioncriterion.hh:148
void setResidualReductionTolerance(Scalar tol)
Sets the residual reduction tolerance.
Definition: weightedresidreductioncriterion.hh:117
Scalar absResidualTolerance() const
Returns the maximum absolute tolerated residual.
Definition: weightedresidreductioncriterion.hh:135
Scalar residualAccuracy() const
Returns the reduction of the weighted maximum of the residual compared to the initial solution...
Definition: weightedresidreductioncriterion.hh:142
Scalar fixPointAccuracy() const
Returns the weighted maximum of the difference between the last two iterative solutions.
Definition: weightedresidreductioncriterion.hh:163
bool failed() const
Returns true if the convergence criterion cannot be met anymore because the solver has broken down...
Definition: weightedresidreductioncriterion.hh:208
Scalar fixPointTolerance() const
Returns the maximum tolerated difference between two iterations to be met before a solution is consid...
Definition: weightedresidreductioncriterion.hh:156
Scalar accuracy() const
Returns the accuracy of the solution at the last update.
Definition: weightedresidreductioncriterion.hh:220
void setInitial(const Vector &curSol, const Vector &curResid)
Set the initial solution of the linear system of equations.
Definition: weightedresidreductioncriterion.hh:169
void setResidualWeight(const Vector &residWeightVec)
Sets the relative weight of each row of the residual.
Definition: weightedresidreductioncriterion.hh:98
void printInitial(std::ostream &os=std::cout) const
Prints the initial information about the convergence behaviour.
Definition: weightedresidreductioncriterion.hh:226
void print(Scalar iter, std::ostream &os=std::cout) const
Prints the information about the convergence behaviour for the current iteration. ...
Definition: weightedresidreductioncriterion.hh:239
Scalar residualReductionTolerance() const
Returns the tolerance of the residual reduction of the solution.
Definition: weightedresidreductioncriterion.hh:123
Convergence criterion which looks at the weighted absolute value of the residual. ...
Definition: weightedresidreductioncriterion.hh:58
void setResidualTolerance(Scalar tol)
Sets the maximum absolute tolerated residual.
Definition: weightedresidreductioncriterion.hh:129
void update(const Vector &curSol, const Vector &curResid)
Definition: weightedresidreductioncriterion.hh:187
Scalar residualWeight(size_t outerIdx, unsigned innerIdx) const
Return the relative weight of a row of the residual.
Definition: weightedresidreductioncriterion.hh:107