All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
fixpointcriterion.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_FIXPOINT_CRITERION_HH
28 #define EWOMS_ISTL_FIXPOINT_CRITERION_HH
29 
30 #include "convergencecriterion.hh"
31 
32 #include <opm/common/Unused.hpp>
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  FixPointCriterion(const CollectiveCommunication& comm) : comm_(comm)
65  {}
66 
67  FixPointCriterion(const CollectiveCommunication& comm,
68  const Vector& weightVec, Scalar reduction)
69  : comm_(comm), weightVec_(weightVec), tolerance_(reduction)
70  {}
71 
88  void setWeight(const Vector& weightVec)
89  { weightVec_ = weightVec; }
90 
107  Scalar weight(int outerIdx, int innerIdx) const
108  { return (weightVec_.size() == 0) ? 1.0 : weightVec_[outerIdx][innerIdx]; }
109 
118  void setTolerance(Scalar tol)
119  { tolerance_ = tol; }
120 
125  Scalar tolerance() const
126  { return tolerance_; }
127 
131  void setInitial(const Vector& curSol, const Vector& curResid OPM_UNUSED)
132  {
133  lastSol_ = curSol;
134  delta_ = 1000 * tolerance_;
135  }
136 
140  void update(const Vector& curSol,
141  const Vector& changeIndicator OPM_UNUSED,
142  const Vector& curResid OPM_UNUSED)
143  {
144  assert(curSol.size() == lastSol_.size());
145 
146  delta_ = 0.0;
147  for (size_t i = 0; i < curSol.size(); ++i) {
148  for (size_t j = 0; j < BlockType::dimension; ++j) {
149  delta_ =
150  std::max(delta_, weight(i, j)*std::abs(curSol[i][j] - lastSol_[i][j]));
151  }
152  }
153 
154  delta_ = comm_.max(delta_);
155  lastSol_ = curSol;
156  }
157 
161  bool converged() const
162  { return accuracy() < tolerance(); }
163 
167  Scalar accuracy() const
168  { return delta_; }
169 
170 private:
171  const CollectiveCommunication& comm_;
172 
173  Vector lastSol_; // solution of the last iteration
174  Vector weightVec_; // solution of the last iteration
175  Scalar delta_; // the maximum of the absolute weighted difference of the
176  // last two iterations
177  Scalar tolerance_; // the maximum allowed delta for the solution to be
178  // considered converged
179 };
180 
182 
183 }} // end namespace Linear, Ewoms
184 
185 #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: fixpointcriterion.hh:161
Scalar tolerance() const
Return the maximum allowed weighted difference between two iterations for the solution considered to ...
Definition: fixpointcriterion.hh:125
void update(const Vector &curSol, const Vector &changeIndicator OPM_UNUSED, const Vector &curResid OPM_UNUSED)
Update the internal members of the convergence criterion with the current solution.
Definition: fixpointcriterion.hh:140
Scalar weight(int outerIdx, int innerIdx) const
Return the relative weight of a primary variable.
Definition: fixpointcriterion.hh:107
void setWeight(const Vector &weightVec)
Sets the relative weight of a primary variable.
Definition: fixpointcriterion.hh:88
Provides a convergence criterion for the linear solvers which looks at the weighted maximum of the di...
Definition: fixpointcriterion.hh:58
Scalar accuracy() const
Returns the accuracy of the solution at the last update.
Definition: fixpointcriterion.hh:167
void setInitial(const Vector &curSol, const Vector &curResid OPM_UNUSED)
Set the initial solution of the linear system of equations.
Definition: fixpointcriterion.hh:131
void setTolerance(Scalar tol)
Set the maximum allowed weighted maximum difference between two iterations.
Definition: fixpointcriterion.hh:118