residreductioncriterion.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_RESID_REDUCTION_CRITERION_HH
28 #define EWOMS_ISTL_RESID_REDUCTION_CRITERION_HH
29 
30 #include "convergencecriterion.hh"
31 
32 #include <opm/common/Unused.hpp>
33 
34 #include <dune/istl/scalarproducts.hh>
35 
36 namespace Ewoms {
37 namespace Linear {
38 
52 template <class Vector>
54 {
55  typedef typename Vector::field_type Scalar;
56 
57 public:
58  ResidReductionCriterion(Dune::ScalarProduct<Vector>& scalarProduct,
59  Scalar tolerance = 1e-6)
60  : scalarProduct_(scalarProduct), tolerance_(tolerance)
61  {}
62 
67  void setTolerance(Scalar tol)
68  { tolerance_ = tol; }
69 
73  Scalar tolerance() const
74  { return tolerance_; }
75 
79  void setInitial(const Vector& curSol OPM_UNUSED, const Vector& curResid)
80  {
81  static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1e10;
82 
83  // make sure that we don't allow an initial error of 0 to avoid
84  // divisions by zero
85  curDefect_ = scalarProduct_.norm(curResid);
86  lastDefect_ = curDefect_;
87  initialDefect_ = std::max(curDefect_, eps);
88  }
89 
93  void update(const Vector& curSol OPM_UNUSED,
94  const Vector& changeIndicator OPM_UNUSED,
95  const Vector& curResid)
96  {
97  lastDefect_ = curDefect_;
98  curDefect_ = scalarProduct_.norm(curResid);
99  }
100 
104  bool converged() const
105  { return accuracy() <= tolerance(); }
106 
110  Scalar accuracy() const
111  { return curDefect_/initialDefect_; }
112 
116  void printInitial(std::ostream& os=std::cout) const
117  {
118  os << std::setw(20) << "iteration ";
119  os << std::setw(20) << "residual ";
120  os << std::setw(20) << "accuracy ";
121  os << std::setw(20) << "rate ";
122  os << std::endl;
123  }
124 
128  void print(Scalar iter, std::ostream& os=std::cout) const
129  {
130  static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1e10;
131 
132  os << std::setw(20) << iter << " ";
133  os << std::setw(20) << curDefect_ << " ";
134  os << std::setw(20) << accuracy() << " ";
135  os << std::setw(20) << (lastDefect_/std::max(eps, curDefect_)) << " ";
136  os << std::endl;
137  }
138 
139 private:
140  Dune::ScalarProduct<Vector>& scalarProduct_;
141 
142  Scalar tolerance_;
143  Scalar initialDefect_;
144  Scalar curDefect_;
145  Scalar lastDefect_;
146 };
147 
149 
150 }} // end namespace Linear,Ewoms
151 
152 #endif
Base class for all convergence criteria which only defines an virtual API.
Definition: convergencecriterion.hh:57
Definition: baseauxiliarymodule.hh:37
Scalar accuracy() const
Returns the accuracy of the solution at the last update.
Definition: residreductioncriterion.hh:110
Provides a convergence criterion which looks at the reduction of the two-norm of the residual for the...
Definition: residreductioncriterion.hh:53
void printInitial(std::ostream &os=std::cout) const
Prints the initial information about the convergence behaviour.
Definition: residreductioncriterion.hh:116
void update(const Vector &curSol OPM_UNUSED, const Vector &changeIndicator OPM_UNUSED, const Vector &curResid)
Definition: residreductioncriterion.hh:93
Scalar tolerance() const
Return the maximum allowed weighted maximum of the reduction of the linear residual.
Definition: residreductioncriterion.hh:73
void setTolerance(Scalar tol)
Set the maximum allowed weighted maximum of the reduction of the linear residual. ...
Definition: residreductioncriterion.hh:67
void setInitial(const Vector &curSol OPM_UNUSED, const Vector &curResid)
Set the initial solution of the linear system of equations.
Definition: residreductioncriterion.hh:79
bool converged() const
Returns true if and only if the convergence criterion is met.
Definition: residreductioncriterion.hh:104
void print(Scalar iter, std::ostream &os=std::cout) const
Prints the information about the convergence behaviour for the current iteration. ...
Definition: residreductioncriterion.hh:128