combinedcriterion.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 
54 template <class Vector, class CollectiveCommunication>
56 {
57  typedef typename Vector::field_type Scalar;
58  typedef typename Vector::block_type BlockType;
59 
60 public:
61  CombinedCriterion(const CollectiveCommunication& comm)
62  : comm_(comm)
63  {}
64 
65  CombinedCriterion(const CollectiveCommunication& comm,
67  Scalar absResidualTolerance = 0.0,
68  Scalar maxResidual = 0.0)
69  : comm_(comm),
70  residualReductionTolerance_(residualReductionTolerance),
71  absResidualTolerance_(absResidualTolerance),
72  maxResidual_(maxResidual)
73  { }
74 
79  { residualReductionTolerance_ = tol; }
80 
85  { return residualReductionTolerance_; }
86 
91  Scalar residualReduction() const
92  { return residualError_/std::max<Scalar>(1e-20, initialResidualError_); }
93 
97  void setAbsResidualTolerance(Scalar tol)
98  { absResidualTolerance_ = tol; }
99 
104  Scalar absResidualTolerance() const
105  { return absResidualTolerance_; }
106 
110  Scalar absResidual() const
111  { return residualError_; }
112 
116  void setInitial(const Vector& curSol, const Vector& curResid) override
117  {
118  updateErrors_(curSol, curSol, curResid);
119  stagnates_ = false;
120 
121  // to avoid divisions by zero, make sure that we don't use an initial error of 0
122  residualError_ = std::max<Scalar>(residualError_,
123  std::numeric_limits<Scalar>::min()*1e10);
124  initialResidualError_ = residualError_;
125  lastResidualError_ = residualError_;
126  }
127 
131  void update(const Vector& curSol, const Vector& changeIndicator, const Vector& curResid) override
132  { updateErrors_(curSol, changeIndicator, curResid); }
133 
137  bool converged() const override
138  {
139  // we're converged if the solution is better than the tolerance
140  // fix-point and residual tolerance.
141  return
144  }
145 
149  bool failed() const override
150  { return !converged() && (stagnates_ || residualError_ > maxResidual_); }
151 
157  Scalar accuracy() const override
158  { return residualError_/initialResidualError_; }
159 
163  void printInitial(std::ostream& os = std::cout) const override
164  {
165  os << std::setw(20) << "iteration ";
166  os << std::setw(20) << "residual ";
167  os << std::setw(20) << "reduction ";
168  os << std::setw(20) << "rate ";
169  os << std::endl;
170  }
171 
175  void print(Scalar iter, std::ostream& os = std::cout) const override
176  {
177  const Scalar eps = std::numeric_limits<Scalar>::min()*1e10;
178 
179  os << std::setw(20) << iter << " ";
180  os << std::setw(20) << absResidual() << " ";
181  os << std::setw(20) << accuracy() << " ";
182  os << std::setw(20) << lastResidualError_/std::max<Scalar>(residualError_, eps) << " ";
183  os << std::endl << std::flush;
184  }
185 
186 private:
187  // update the weighted absolute residual
188  void updateErrors_(const Vector& curSol OPM_UNUSED, const Vector& changeIndicator, const Vector& curResid)
189  {
190  lastResidualError_ = residualError_;
191  residualError_ = 0.0;
192  stagnates_ = true;
193  for (size_t i = 0; i < curResid.size(); ++i) {
194  for (unsigned j = 0; j < BlockType::dimension; ++j) {
195  residualError_ =
196  std::max<Scalar>(residualError_,
197  std::abs(curResid[i][j]));
198 
199  if (stagnates_ && changeIndicator[i][j] != 0.0)
200  // only stagnation means that we've failed!
201  stagnates_ = false;
202  }
203  }
204 
205  residualError_ = comm_.max(residualError_);
206 
207  // the linear solver only stagnates if all processes stagnate
208  stagnates_ = comm_.min(stagnates_);
209  }
210 
211  const CollectiveCommunication& comm_;
212 
213  // the infinity norm of the residual of the last iteration
214  Scalar lastResidualError_;
215 
216  // the infinity norm of the residual of the current iteration
217  Scalar residualError_;
218 
219  // the infinity norm of the residual of the initial solution
220  Scalar initialResidualError_;
221 
222  // the minimum reduction of the residual norm where the solution is to be considered
223  // converged
224  Scalar residualReductionTolerance_;
225 
226  // the maximum residual norm for the residual for the solution to be considered to be
227  // converged
228  Scalar absResidualTolerance_;
229 
230  // The maximum error which is tolerated before we fail.
231  Scalar maxResidual_;
232 
233  // does the linear solver seem to stagnate, i.e. were the last two solutions
234  // identical?
235  bool stagnates_;
236 };
237 
239 
240 }} // end namespace Linear,Ewoms
241 
242 #endif
Base class for all convergence criteria which only defines an virtual API.
Definition: convergencecriterion.hh:57
Scalar absResidual() const
Returns the infinity norm of the absolute residual.
Definition: combinedcriterion.hh:110
bool failed() const override
Returns true if the convergence criterion cannot be met anymore because the solver has broken down...
Definition: combinedcriterion.hh:149
Definition: baseauxiliarymodule.hh:37
Scalar residualReduction() const
Returns the reduction of the maximum of the residual compared to the initial solution.
Definition: combinedcriterion.hh:91
void update(const Vector &curSol, const Vector &changeIndicator, const Vector &curResid) override
Update the internal members of the convergence criterion with the current solution.
Definition: combinedcriterion.hh:131
void setInitial(const Vector &curSol, const Vector &curResid) override
Set the initial solution of the linear system of equations.
Definition: combinedcriterion.hh:116
Scalar absResidualTolerance() const
Returns the tolerated maximum of the the infinity norm of the absolute residual.
Definition: combinedcriterion.hh:104
Convergence criterion which looks at the absolute value of the residual and fails if the linear solve...
Definition: combinedcriterion.hh:55
Scalar accuracy() const override
Returns the accuracy of the solution at the last update.
Definition: combinedcriterion.hh:157
void setAbsResidualTolerance(Scalar tol)
Sets the maximum absolute tolerated residual.
Definition: combinedcriterion.hh:97
void setResidualReductionTolerance(Scalar tol)
Sets the residual reduction tolerance.
Definition: combinedcriterion.hh:78
void print(Scalar iter, std::ostream &os=std::cout) const override
Prints the information about the convergence behaviour for the current iteration. ...
Definition: combinedcriterion.hh:175
bool converged() const override
Returns true if and only if the convergence criterion is met.
Definition: combinedcriterion.hh:137
Scalar residualReductionTolerance() const
Returns the tolerance of the residual reduction of the solution.
Definition: combinedcriterion.hh:84
void printInitial(std::ostream &os=std::cout) const override
Prints the initial information about the convergence behaviour.
Definition: combinedcriterion.hh:163