DiscreteFunction.hpp

#pragma once

#include <vector>

#include <dune/common/std/optional.hh>
#include <dune/functions/common/defaultderivativetraits.hh>
#include <dune/functions/functionspacebases/defaultnodetorangemap.hh>
#include <dune/functions/functionspacebases/flatvectorview.hh>
#include <dune/functions/gridfunctions/gridviewentityset.hh>
#include <dune/typetree/childextraction.hh>

#include <amdis/LinearAlgebra.hpp>
#include <amdis/typetree/FiniteElementType.hpp>
#include <amdis/typetree/RangeType.hpp>
#include <amdis/typetree/TreePath.hpp>

namespace AMDiS
{
  template <class GB, class VT>
  class DOFVector;
  
  /// \class DiscreteFunction
  /// \brief A view on a subspace of a \ref DOFVector
  /**
  * \ingroup GridFunctions
  *
  * \tparam GB  Type of the global basis
  * \tparam VT  Coefficient type of the DOFVector
  * \tparam TP  A realization of \ref Dune::TypeTree::HybridTreePath
  * \tparam is_const  Specifies whether a const or mutable view is implemented.
  *
  * **Requirements:**
  * - GB models \ref Dune::Functions::Concept::GlobalBasis
  **/
  template <class GB, class VT, class TP, bool is_const>
  class DiscreteFunction;

#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
  // Deduction guide for DiscreteFunction class
  template <class GB, class VT>
  DiscreteFunction(DOFVector<GB, VT> const& dofVector)
    -> DiscreteFunction<GB, VT, Dune::TypeTree::HybridTreePath<>, true>;

  template <class GB, class VT>
  DiscreteFunction(DOFVector<GB, VT>& dofVector)
    -> DiscreteFunction<GB, VT, Dune::TypeTree::HybridTreePath<>, false>;
#endif


  /// A Generator for a const \ref DiscreteFunction
  template <class GlobalBasis, class ValueType,
            class PreTreePath = Dune::TypeTree::HybridTreePath<>>
  auto makeDiscreteFunction(DOFVector<GlobalBasis, ValueType> const& dofVector,
                            PreTreePath const& preTreePath = {})
  {
    auto treePath = makeTreePath(preTreePath);
    return DiscreteFunction<GlobalBasis, ValueType, decltype(treePath), true>{dofVector, treePath};
  }

  /// A Generator for a mutable \ref DiscreteFunction
  template <class GlobalBasis, class ValueType,
            class PreTreePath = Dune::TypeTree::HybridTreePath<>>
  auto makeDiscreteFunction(DOFVector<GlobalBasis, ValueType>& dofVector,
                            PreTreePath const& preTreePath = {})
  {
    auto treePath = makeTreePath(preTreePath);
    return DiscreteFunction<GlobalBasis, ValueType, decltype(treePath), false>{dofVector, treePath};
  }



  /// A Const DiscreteFunction
  template <class GB, class VT, class TP>
  class DiscreteFunction<GB,VT,TP,true>
  {
  private:
    using GlobalBasis = GB;
    using TreePath = TP;

    using Tree = typename GlobalBasis::LocalView::Tree;
    using SubTree = typename Dune::TypeTree::ChildForTreePath<Tree, TreePath>;
    using NodeToRangeEntry = Dune::Functions::DefaultNodeToRangeMap<SubTree>;

    using GridView = typename GlobalBasis::GridView;

  public:
    /// Set of entities the DiscreteFunction is defined on
    using EntitySet = Dune::Functions::GridViewEntitySet<GridView, 0>;

    /// Global coordinates of the EntitySet
    using Domain = typename EntitySet::GlobalCoordinate;

    /// Range type of this DiscreteFunction
    using Range = RangeType_t<SubTree,VT>;

    /// \brief This GridFunction has no derivative function, it can be created
    /// by \ref DiscreteGridFunction.
    enum { hasDerivative = false };

  public:
    /// A LocalFunction representing the derivative of the DOFVector on a bound element
    template <class Type>
    class DerivativeLocalFunctionBase;

    class GradientLocalFunction;
    class PartialLocalFunction;
    class DivergenceLocalFunction;

    /// A LocalFunction representing the value the DOFVector on a bound element
    class LocalFunction;

  public:
    /// Constructor. Stores a pointer to the dofVector and a copy of the treePath.
    DiscreteFunction(DOFVector<GB,VT> const& dofVector, TP const& treePath = {})
      : dofVector_(&dofVector)
      , treePath_(treePath)
      , entitySet_(dofVector.basis()->gridView())
      , nodeToRangeEntry_(Dune::Functions::makeDefaultNodeToRangeMap(*dofVector.basis(), treePath))
    {}

    /// \brief Evaluate DiscreteFunction in global coordinates. NOTE: expensive
    Range operator()(Domain const& x) const;

    /// \brief Create a local function for this view on the DOFVector. \relates LocalFunction
    LocalFunction makeLocalFunction() const
    {
      return LocalFunction{*this};
    }

    /// \brief Return a \ref Dune::Functions::GridViewEntitySet
    EntitySet const& entitySet() const
    {
      return entitySet_;
    }

    /// \brief Return global basis bound to the DOFVector
    std::shared_ptr<GlobalBasis const> basis() const
    {
      return dofVector_->basis();
    }

    /// \brief Return treePath associated with this discrete function
    TreePath const& treePath() const
    {
      return treePath_;
    }

    /// \brief Return const coefficient vector
    DOFVector<GB,VT> const& coefficients() const
    {
      return *dofVector_;
    }

  protected:
    DOFVector<GB,VT> const* dofVector_;
    TreePath treePath_;
    EntitySet entitySet_;
    NodeToRangeEntry nodeToRangeEntry_;
  };



  /// A mutable view on the subspace of a DOFVector, \relates DiscreteFunction
  template <class GB, class VT, class TP>
  class DiscreteFunction<GB,VT,TP,false>
      : public DiscreteFunction<GB, VT, TP, true>
  {
    using Self = DiscreteFunction<GB, VT, TP, false>;
    using Super = DiscreteFunction<GB, VT, TP, true>;

    using GlobalBasis = GB;
    using TreePath = TP;

  public:
    /// Constructor. Stores a pointer to the mutable `dofvector`.
    DiscreteFunction(DOFVector<GB,VT>& dofVector, TP const& treePath = {})
      : Super(dofVector, treePath)
      , mutableDofVector_(&dofVector)
    {}

  public:
    /// \brief Interpolation of GridFunction to DOFVector, assuming that there is no
    /// reference to this DOFVector in the expression.
    /**
     * **Example:**
     * ```
     * auto v = makeDiscreteFunction(prob.solutionVector(),0);
     * v.interpolate_noalias([](auto const& x) { return x[0]; });
     * ```
     **/
    template <class Expr, class Tag = tag::average>
    void interpolate_noalias(Expr&& expr, Tag strategy = {});

    /// \brief Interpolation of GridFunction to DOFVector
    /**
     * **Example:**
     * ```
     * auto v = makeDiscreteFunction(prob.solutionVector(),0);
     * v.interpolate(v + [](auto const& x) { return x[0]; });
     * ```
     * Allows to have a reference to the DOFVector in the expression, e.g. as
     * \ref DiscreteFunction or \ref gradientAtQP() of a DiscreteFunction.
     **/
    template <class Expr, class Tag = tag::average>
    void interpolate(Expr&& expr, Tag strategy = {});

    /// \brief Interpolation of GridFunction to DOFVector, alias to \ref interpolate()
    template <class Expr>
    Self& operator<<(Expr&& expr)
    {
      interpolate(FWD(expr));
      return *this;
    }

    /// \brief interpolate `(*this) + expr` to DOFVector
    template <class Expr>
    Self& operator+=(Expr&& expr)
    {
      interpolate((*this) + expr);
      return *this;
    }

    /// \brief interpolate `(*this) - expr` to DOFVector
    template <class Expr>
    Self& operator-=(Expr&& expr)
    {
      interpolate((*this) - expr);
      return *this;
    }

    /// Return the mutable DOFVector
    DOFVector<GB,VT>& coefficients() { return *mutableDofVector_; }

    /// Return the const DOFVector
    using Super::coefficients;

  protected:
    DOFVector<GB,VT>* mutableDofVector_;
  };

} // end namespace AMDiS

#include "DiscreteLocalFunction.inc.hpp"
#include "DiscreteFunction.inc.hpp"