fixpointcriterion.hh

Go to the documentation of this file.

// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
  Copyright (C) 2009-2013 by Andreas Lauser

  This file is part of the Open Porous Media project (OPM).

  OPM is free software: you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation, either version 2 of the License, or
  (at your option) any later version.

  OPM is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef EWOMS_ISTL_FIXPOINT_CRITERION_HH
#define EWOMS_ISTL_FIXPOINT_CRITERION_HH

#include "convergencecriterion.hh"

namespace Ewoms {
template <class Vector, class CollectiveCommunication>
class FixPointCriterion : public ConvergenceCriterion<Vector>
{
    typedef typename Vector::field_type Scalar;
    typedef typename Vector::block_type BlockType;

public:
    FixPointCriterion(const CollectiveCommunication &comm) : comm_(comm)
    {}

    FixPointCriterion(const CollectiveCommunication &comm,
                      const Vector &weightVec, Scalar reduction)
        : comm_(comm), weightVec_(weightVec), tolerance_(reduction)
    {}

    void setWeight(const Vector &weightVec)
    { weightVec_ = weightVec; }

    Scalar weight(int outerIdx, int innerIdx) const
    { return (weightVec_.size() == 0) ? 1.0 : weightVec_[outerIdx][innerIdx]; }

    void setTolerance(Scalar tol)
    { tolerance_ = tol; }

    Scalar tolerance() const
    { return tolerance_; }

    void setInitial(const Vector &curSol, const Vector &curResid)
    {
        lastSol_ = curSol;
        delta_ = 1000 * tolerance_;
    }

    void update(const Vector &curSol, const Vector &curResid)
    {
        assert(curSol.size() == lastSol_.size());

        delta_ = 0.0;
        for (size_t i = 0; i < curSol.size(); ++i) {
            for (size_t j = 0; j < BlockType::dimension; ++j) {
                delta_ = std::max<Scalar>(delta_, weight(i, j)
                                                  * std::abs(curSol[i][j]
                                                             - lastSol_[i][j]));
            }
        }

        delta_ = comm_.max(delta_);
        lastSol_ = curSol;
    }

    bool converged() const
    { return delta_ < tolerance(); }

    Scalar accuracy() const
    { return delta_; }

private:
    const CollectiveCommunication &comm_;

    Vector lastSol_;   // solution of the last iteration
    Vector weightVec_; // solution of the last iteration
    Scalar delta_;     // the maximum of the absolute weighted difference of the
                       // last two iterations
    Scalar tolerance_; // the maximum allowed delta for the solution to be
                       // considered converged
};


} // end namespace Ewoms

#endif