Marker.hpp

#pragma once

#include <limits>
#include <string>
#include <utility>

#include <dune/common/std/optional.hh>
#include <dune/grid/common/grid.hh>

#include <amdis/common/ConceptsBase.hpp>
#include <amdis/common/TypeTraits.hpp>

#include <amdis/gridfunctions/GridFunction.hpp>

#include <amdis/AdaptInfo.hpp>
#include <amdis/Flag.hpp>
#include <amdis/Initfile.hpp>

namespace AMDiS
{
  /**
   * \ingroup Adaption
   *
   * \brief
   * Base class for all markers.
   */
  template <class Grid>
  class Marker
  {
  protected:
    using GridView = typename Grid::LeafGridView;
    using Element  = typename GridView::template Codim<0>::Entity;

  public:
    /// Constructor.
    Marker() = default;

    /// Constructor.
    Marker(std::string const& name, std::shared_ptr<Grid> const& grid)
      : name_(name)
      , grid_(grid)
    {
      Parameters::get(name + "->info", info_);
      Parameters::get(name + "->max refinement level", maxRefineLevel_);
      Parameters::get(name + "->min refinement level", minRefineLevel_);
    }

    /// Destructor.
    virtual ~Marker() = default;

    /// Marks element with newMark. If \ref maximumMarking_ is set, the element
    /// is marked only if the new mark is bigger than the old one. The return
    /// value specifies whether the element has been marked, or not.
    void mark(Element const& elem, int newMark);

    /// Called before traversal.
    virtual void initMarking(AdaptInfo& adaptInfo);

    /// Called after traversal.
    virtual void finishMarking(AdaptInfo& adaptInfo);

    /// Marks one element.
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem) = 0;

    /// Marking of the whole grid.
    virtual Flag markGrid(AdaptInfo& adaptInfo);

    /// Returns \ref elMarkRefine_ of the Marker
    int elMarkRefine() const
    {
      return elMarkRefine_;
    }

    /// Returns \ref elMarkCoarsen_ of the Marker
    int elMarkCoarsen() const
    {
      return elMarkCoarsen_;
    }

    /// Returns \ref name_ of the Marker
    std::string const& name() const
    {
      return name_;
    }

    /// Sets \ref maximumMarking_.
    void setMaximumMarking(bool maxMark)
    {
      maximumMarking_ = maxMark;
    }

    /// Sets \ref refineAllowed_.
    void allowRefinement(bool allow)
    {
      refineAllowed_ = allow;
    }

    /// Sets \ref coarsenAllowed_.
    void allowCoarsening(bool allow)
    {
      coarsenAllowed_ = allow;
    }

  protected:
    /// Name of the marker.
    std::string name_;

    /// Pointer to the grid
    std::shared_ptr<Grid> grid_;

    /// If true, elements are marked only if the new value is bigger than
    /// the current marking.
    bool maximumMarking_ = false;

    /// Info level.
    int info_ = 0;

    /// Counter for elements marked for refinement
    int elMarkRefine_ = 0;

    /// Counter for elements marked for coarsening
    int elMarkCoarsen_ = 0;

    /// Maximal level of all elements.
    int maxRefineLevel_ = std::numeric_limits<int>::max();

    /// Minimal level of all elements.
    int minRefineLevel_ = 0;

    /// Allow elements to be marked for refinement
    bool refineAllowed_ = true;

    /// Allow elements to be marked for coarsening
    bool coarsenAllowed_ = true;
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Base class for all estimator-based markers.
   */
  template <class Grid>
  class EstimatorMarker
      : public Marker<Grid>
  {
  protected:
    using Super = Marker<Grid>;
    using Element   = typename Super::Element;
    using Estimates = std::vector<double>;

  public:
    /// Constructor.
    EstimatorMarker() = default;

    /// Constructor.
    EstimatorMarker(std::string name, std::string component, Estimates const& est,
                    std::shared_ptr<Grid> const& grid)
      : Super{std::move(name), grid}
      , component_(std::move(component))
      , est_(est)
    {
      Parameters::get(this->name_ + "->p", p_);
    }

    /// Can be used by sub classes. Called before traversal.
    void initMarking(AdaptInfo& adaptInfo) override;

    /// Marks one element.
    void markElement(AdaptInfo& adaptInfo, Element const& elem) override;

    /// Creates a scalar marker depending on the strategy set in parameters.
    static std::unique_ptr<EstimatorMarker<Grid> > createMarker(std::string const& name,
      std::string const& component, Estimates const& est, std::shared_ptr<Grid> const& grid);

  protected:
    /// Problem component to work on
    std::string component_;

    /// Pointer to the local error estimates
    Estimates est_;

    /// Estimation sum
    double estSum_ = 0.0;

    /// Estmation maximum
    double estMax_ = 0.0;

    /// Lower limit for error estimation, from which an element is marked for
    /// refinement
    double markRLimit_;

    /// Upper limit for error estimation, from which an element is marked for
    /// coarsening
    double markCLimit_;

    /// Power in estimator norm
    double p_ = 2.0;
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Global refinement.
   */
  template <class Grid>
  class GRMarker
      : public EstimatorMarker<Grid>
  {
    using Super = EstimatorMarker<Grid>;
    using Element   = typename Super::Element;
    using Estimates = typename Super::Estimates;

  public:
    /// Constructor.
    GRMarker(std::string const& name, std::string const& component, Estimates const& est,
             std::shared_ptr<Grid> const& grid)
      : Super{name, component, est, grid}
    {}

    void markElement(AdaptInfo& adaptInfo, Element const& elem) override
    {
      if (this->refineAllowed_)
        this->mark(elem, 1);
    }
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Maximum strategy.
   */

  template <class Grid>
  class MSMarker
      : public EstimatorMarker<Grid>
  {
    using Super = EstimatorMarker<Grid>;
    using Estimates = typename Super::Estimates;

  public:
    /// Constructor.
    MSMarker(std::string const& name, std::string component, Estimates const& est,
             std::shared_ptr<Grid> const& grid)
      : Super{name, std::move(component), est, grid}
    {
      Parameters::get(name + "->MSGamma", msGamma_);
      Parameters::get(name + "->MSGammaC", msGammaC_);
    }

    void initMarking(AdaptInfo& adaptInfo) override;

  protected:
    /// Marking parameter.
    double msGamma_ = 0.5;

    /// Marking parameter.
    double msGammaC_ = 0.1;
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Equidistribution strategy
   */

  template <class Grid>
  class ESMarker
      : public EstimatorMarker<Grid>
  {
    using Super = EstimatorMarker<Grid>;
    using Estimates = typename Super::Estimates;

  public:
    /// Constructor.
    ESMarker(std::string const& name, std::string component, Estimates const& est,
             std::shared_ptr<Grid> const& grid)
      : Super{name, std::move(component), est, grid}
    {
      Parameters::get(name + "->ESTheta", esTheta_);
      Parameters::get(name + "->ESThetaC", esThetaC_);
    }

    void initMarking(AdaptInfo& adaptInfo) override;

  protected:
    /// Marking parameter.
    double esTheta_ = 0.9;

    /// Marking parameter.
    double esThetaC_ = 0.2;
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Guaranteed error reduction strategy
   */

  template <class Grid>
  class GERSMarker
      : public EstimatorMarker<Grid>
  {
    using Super = EstimatorMarker<Grid>;
    using Element   = typename Super::Element;
    using Estimates = typename Super::Estimates;

  public:
    /// Constructor.
    GERSMarker(std::string const& name, std::string component, Estimates const& est,
               std::shared_ptr<Grid> const& grid)
      : Super{name, std::move(component), est, grid}
    {
      Parameters::get(name + "->GERSThetaStar", gersThetaStar_);
      Parameters::get(name + "->GERSNu", gersNu_);
      Parameters::get(name + "->GERSThetaC", gersThetaC_);
    }

    Flag markGrid(AdaptInfo& adaptInfo) override;

  protected:
    /// Refinement marking function.
    void markElementForRefinement(AdaptInfo& adaptInfo, Element const& elem);

    /// Coarsening marking function.
    void markElementForCoarsening(AdaptInfo& adaptInfo, Element const& elem);

  protected:
    /// Marking parameter.
    double gersSum_ = 0.0;

    /// Marking parameter.
    double oldErrSum_ = 0.0;

    /// Marking parameter.
    double gersThetaStar_ = 0.6;

    /// Marking parameter.
    double gersNu_ = 0.1;

    /// Marking parameter.
    double gersThetaC_ = 0.1;
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Marker based on an indicator given as grid-function.
   *
   * On each element the grid-function is evaluated in the barycenter. The returned
   * values must be an integer (or convertible to an integer) indicating the
   * desired refinement level of this element. The element is marked for refinement
   * if the current level is < the desired level or coarsened, if >.
   *
   * \tparam Grid     An Implementation of the \ref Dune::Grid interface
   * \tparam GridFct  A grid-function with `Range` convertible to `int`
   */
  template <class Grid, class GridFct>
  class GridFunctionMarker
      : public Marker<Grid>
  {
    using Super = Marker<Grid>;
    using Element = typename Super::Element;

    template <class GF>
    using IsGridFunction = decltype(localFunction(std::declval<GF>()));

    static_assert(Dune::Std::is_detected<IsGridFunction,GridFct>::value, "");

  public:
    /// Constructor.
    template <class GF>
    GridFunctionMarker(std::string const& name, std::shared_ptr<Grid> const& grid, GF&& gf)
      : Super{name, grid}
      , gridFct_{makeGridFunction(FWD(gf), grid->leafGridView())}
    {}

    /// \brief Implementation of \ref Marker::markElement. Does nothing since marking is
    /// done in \ref markGrid().
    void markElement(AdaptInfo& adaptInfo, Element const& elem) final {}

    /// Mark element for refinement, if grid-function \ref gridFct_ evaluates to
    /// a value larger than the current level and is marked for coarsening of
    /// the result is less than the current value.
    Flag markGrid(AdaptInfo& adaptInfo) override;

  private:
    /// Indicator function for adaptation
    GridFct gridFct_;
  };


#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
  // Deduction guide for GridFunctionMarker class
  template <class Grid, class GF>
  GridFunctionMarker(std::string const& name, std::shared_ptr<Grid> const& grid, GF&& gf)
    -> GridFunctionMarker<Grid,
          TYPEOF( makeGridFunction(FWD(gf), grid->leafGridView()) )>;
#endif

  // Generator function for GridFunctionMarker class
  template <class Grid, class PreGridFct>
  auto makeGridFunctionMarker(std::string const& name, std::shared_ptr<Grid> const& grid,
                              PreGridFct&& preGridFct)
  {
    auto gridFct = makeGridFunction(FWD(preGridFct), grid->leafGridView());
    return GridFunctionMarker<Grid,TYPEOF(gridFct)>{name, grid, std::move(gridFct)};
  }

} // end namespace AMDiS

#include "Marker.inc.hpp"