multistagetimestepper.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:
/*****************************************************************************
 *   See the file COPYING for full copying permissions.                      *
 *                                                                           *
 *   This program 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 3 of the License, or       *
 *   (at your option) any later version.                                     *
 *                                                                           *
 *   This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.   *
 *****************************************************************************/
#ifndef DUMUX_TIMESTEPPING_MULTISTAGE_TIMESTEPPER_HH
#define DUMUX_TIMESTEPPING_MULTISTAGE_TIMESTEPPER_HH
 
#include <memory>
#include <vector>
#include <cmath>
 
namespace Dumux::Experimental {
 
template<class Scalar>
class MultiStageMethod;
 
template<class Scalar>
class MultiStageParams
{
    struct Params {
        Scalar alpha, betaDt, timeAtStage, dtFraction;
        bool skipTemporal, skipSpatial;
    };
public:
    MultiStageParams(const MultiStageMethod<Scalar>& m, std::size_t i, const Scalar t, const Scalar dt)
    : size_(i+1)
    {
        params_.resize(size_);
        for (std::size_t k = 0; k < size_; ++k)
        {
            auto& p = params_[k];
            p.alpha = m.temporalWeight(i, k);
            p.betaDt = m.spatialWeight(i, k)*dt;
            p.timeAtStage = t + m.timeStepWeight(k)*dt;
            p.dtFraction = m.timeStepWeight(k);
 
            using std::abs;
            p.skipTemporal = (abs(p.alpha) < 1e-6);
            p.skipSpatial = (abs(p.betaDt) < 1e-6);
        }
    }
 
    std::size_t size () const
    { return size_; }
 
    Scalar temporalWeight (std::size_t k) const
    { return params_[k].alpha; }
 
    Scalar spatialWeight (std::size_t k) const
    { return params_[k].betaDt; }
 
    Scalar timeAtStage (std::size_t k) const
    { return params_[k].timeAtStage; }
 
    Scalar timeStepFraction (std::size_t k) const
    { return params_[k].dtFraction; }
 
    Scalar skipTemporal (std::size_t k) const
    { return params_[k].skipTemporal; }
 
    Scalar skipSpatial (std::size_t k) const
    { return params_[k].skipSpatial; }
 
private:
    std::size_t size_;
    std::vector<Params> params_;
};
 
template<class PDESolver>
class MultiStageTimeStepper
{
    using Variables = typename PDESolver::Variables;
    using Scalar = typename Variables::Scalar;
    using StageParams = MultiStageParams<Scalar>;
 
public:
 
    MultiStageTimeStepper(std::shared_ptr<PDESolver> pdeSolver,
                          std::shared_ptr<const MultiStageMethod<Scalar>> msMethod)
    : pdeSolver_(pdeSolver)
    , msMethod_(msMethod)
    {}
 
    void step(Variables& vars, const Scalar t, const Scalar dt)
    {
        // make sure there are no traces of previous stages
        pdeSolver_->assembler().clearStages();
 
        for (auto stageIdx = 1UL; stageIdx <= msMethod_->numStages(); ++stageIdx)
        {
            // extract parameters for this stage from the time stepping method
            auto stageParams = std::make_shared<StageParams>(*msMethod_, stageIdx, t, dt);
 
            // prepare the assembler for this stage
            pdeSolver_->assembler().prepareStage(vars, stageParams);
 
            // assemble & solve
            pdeSolver_->solve(vars);
        }
 
        // clear traces of previously registered stages
        pdeSolver_->assembler().clearStages();
    }
 
    void setMethod(std::shared_ptr<const MultiStageMethod<Scalar>> msMethod)
    { msMethod_ = msMethod; }
 
private:
    std::shared_ptr<PDESolver> pdeSolver_;
    std::shared_ptr<const MultiStageMethod<Scalar>> msMethod_;
};
 
} // end namespace Dumux::Experimental
 
#endif