functionalexpression.hh

#ifndef __DUNE_ACFEM_OPERATORS_FUNCTIONALS_FUNCTIONALEXPRESSION_HH__
#define __DUNE_ACFEM_OPERATORS_FUNCTIONALS_FUNCTIONALEXPRESSION_HH__
 
#include "../../common/types.hh"
#include "../../common/literals.hh"
#include "../../expressions/operationtraits.hh"
#include "../../tensors/tensor.hh"
 
#include "linearfunctional.hh"
#include "functionaltraits.hh"
#include "closuretraits.hh"
#include "operations/inner.hh"
 
namespace Dune {
 
  namespace ACFem {
 
    namespace LinearFunctional {
 
      using namespace Literals;
 
      template<class Operation, class T>
      class UnaryLinearFunctionalExpression;
 
      template<class Operation, class T0, class T1>
      class BinaryLinearFunctionalExpression;
 
      template<class T>
      class UnaryLinearFunctionalExpression<MinusOperation, T>
        : public LinearFunctional<typename std::decay_t<T>::DiscreteFunctionSpaceType>
        , public Expressions::Storage<OperationTraits<MinusOperation>, T>
      {
        using BaseType = LinearFunctional<typename std::decay_t<T>::DiscreteFunctionSpaceType>;
        using StorageType = Expressions::Storage<OperationTraits<MinusOperation>, T>;
       public:
        using typename BaseType::DiscreteFunctionSpaceType;
        using typename BaseType::RangeFieldType;
        using StorageType::operand;
        using StorageType::operation;
 
        template<class TArg, std::enable_if_t<std::is_constructible<T, TArg>::value, int> = 0>
        UnaryLinearFunctionalExpression(TArg&& t)
          : StorageType(std::forward<TArg>(t))
        {}
 
        template<class FieldArg, class DiscreteFunction,
                 std::enable_if_t<std::is_constructible<RangeFieldType, FieldArg>::value, int> = 0>
        void coefficients(FieldArg&& s, DiscreteFunction& values) const
        {
          static_assert(std::is_same<DiscreteFunctionSpaceType, typename DiscreteFunction::DiscreteFunctionSpaceType>::value,
                        "DiscreteFunctionSpaceType of functional and discrete storage has to coincide");
 
          operand(0_c).coefficients(-s, values);
        }
 
        template<class DiscreteFunction>
        void coefficients(DiscreteFunction& values) const
        {
          values.clear();
          coefficients(1_f, values);
        }
 
        template<class DiscreteFunction>
        auto operator()(const DiscreteFunction& df) const
        {
          return inner(*this, df);
        }
 
        const DiscreteFunctionSpaceType& space() const
        {
          return arg().space();
        }
 
        std::string name() const
        {
          return operationName(operation(), arg().name());
        }
 
       protected:
        using StorageType::t0_;
        constexpr auto&& arg() && { return t0_; }
        constexpr auto& arg() & { return t0_; }
        constexpr const auto& arg() const& { return t0_; }
      };
 
      template<class Operation, class T0, class T1>
      class BinaryLinearFunctionalExpression
        : public LinearFunctional<typename std::decay_t<T1>::DiscreteFunctionSpaceType>
        , public Expressions::Storage<OperationTraits<Operation>, T0, T1>
      {
        using BaseType = LinearFunctional<typename std::decay_t<T1>::DiscreteFunctionSpaceType>;
        using StorageType = Expressions::Storage<OperationTraits<Operation>, T0, T1>;
       public:
        using typename BaseType::DiscreteFunctionSpaceType;
        using typename BaseType::RangeFieldType;
        using StorageType::operand;
        using StorageType::operation;
 
        template<class T0Arg, class T1Arg,
                 std::enable_if_t<(std::is_constructible<T0, T0Arg>::value
                                   && std::is_constructible<T1, T1Arg>::value
                   ), int> = 0>
        BinaryLinearFunctionalExpression(T0Arg&& t0, T1Arg&& t1)
          : StorageType(std::forward<T0Arg>(t0), std::forward<T1>(t1))
        {}
 
        template<class FieldArg, class DiscreteFunction,
                 std::enable_if_t<std::is_constructible<RangeFieldType, FieldArg>::value, int> = 0>
        void coefficients(FieldArg&& s, DiscreteFunction& values) const
        {
          static_assert(std::is_same<DiscreteFunctionSpaceType, typename DiscreteFunction::DiscreteFunctionSpaceType>::value,
                        "DiscreteFunctionSpaceType of functional and discrete storage has to coincide");
 
          // gurk it together
          if constexpr (std::is_same<Operation, PlusOperation>::value) {
            left().coefficients(s, values);
            right().coefficients(s, values);
          } else if constexpr (std::is_same<Operation, MinusOperation>::value) {
            left().coefficients(s, values);
            right().coefficients(-s, values);
          } else if constexpr (std::is_same<Operation, SMultiplyOperation>::value) {
            // The factor is always the left operand
            right().coefficients(left()*s, values);
          }
        }
 
        template<class DiscreteFunction>
        void coefficients(DiscreteFunction& df) const
        {
          coefficients(1_f, df);
        }
 
        template<class DiscreteFunction>
        auto operator()(const DiscreteFunction& df) const
        {
          return inner(*this, df);
        }
 
        const DiscreteFunctionSpaceType& space() const
        {
          return right().space();
        }
 
        std::string name() const
        {
          return operationName(operation(), left().name(), right().name());
        }
 
       protected:
        using StorageType::t0_;
        using StorageType::t1_;
        constexpr auto&& left() && { return t0_; }
        constexpr auto& left() & { return t0_; }
        constexpr const auto& left() const& { return t0_; }
        constexpr auto&& right() && { return t1_; }
        constexpr auto& right() & { return t1_; }
        constexpr const auto& right() const& { return t1_; }
      };
 
      template<class Op, class T, std::enable_if_t<IsLinearFunctional<T>::value, int> = 0>
      constexpr auto operate(Expressions::DontOptimize, OperationTraits<Op>, T&& t)
      {
        DUNE_ACFEM_RECORD_OPTIMIZATION;
 
        return UnaryLinearFunctionalExpression<Op, T>(std::forward<T>(t));
      }
 
      template<class Op, class T0, class T1,
               std::enable_if_t<((IsLinearFunctional<T0>::value || TensorTraits<T0>::rank == 0)
                                 && IsLinearFunctional<T1>::value
        ), int> = 0>
      constexpr auto operate(Expressions::DontOptimize, OperationTraits<Op>, T0&& t0, T1&& t1)
      {
        DUNE_ACFEM_RECORD_OPTIMIZATION;
 
        return BinaryLinearFunctionalExpression<Op, T0, T1>(std::forward<T0>(t0), std::forward<T1>(t1));
      }
 
    }
 
  } // ACFem::
 
} // Dune::
 
#endif // __DUNE_ACFEM_OPERATORS_FUNCTIONALS_FUNCTIONALEXPRESSION_HH__