fieldtensor.hh

#ifndef __DUNE_ACFEM_TENSORS_BINDINGS_DUNE_FIELDTENSOR_HH__
#define __DUNE_ACFEM_TENSORS_BINDINGS_DUNE_FIELDTENSOR_HH__
 
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <dune/common/hybridutilities.hh>
 
#include "../../../expressions/terminal.hh"
#include "../../tensorbase.hh"
#include "../../operations/assign.hh"
#include "../../optimization/policy.hh"
#include "fieldvectortraits.hh"
#include "densevectorview.hh"
 
namespace Dune {
 
  namespace ACFem {
 
    namespace Tensor {
 
      template<class FV>
      class FieldVectorTensor
        : public TensorBase<typename FieldTraits<std::decay_t<FV> >::field_type,
                            FieldVectorSignature<FV>,
                            FieldVectorTensor<FV> >
        , public Expressions::SelfExpression<FieldVectorTensor<FV> >
        , public MPL::UniqueTags<ConditionalType<
                                   !std::is_reference<FV>::value,
                                   IndependentExpression,
                                   void> >
      {
        using ThisType = FieldVectorTensor;
        using DataType = FV;
        using ArrayType = std::decay_t<FV>;
        using BaseType = TensorBase<typename FieldTraits<ArrayType>::field_type,
                                    FieldVectorSignature<ArrayType>,
                                    FieldVectorTensor<FV> >;
       public:
        using typename BaseType::Signature;
        using typename BaseType::FieldType;
        using BaseType::rank;
 
        FieldVectorTensor()
        {}
 
        template<class T, std::enable_if_t<SameDecay<T, ArrayType>::value, int> = 0>
        FieldVectorTensor(T&& t)
          : data_(std::forward<T>(t))
        {}
 
        template<class T,
                 std::enable_if_t<(std::is_constructible<DataType, T>::value
                                   && !SameDecay<T, ArrayType>::value
                   ), int> = 0>
        explicit FieldVectorTensor(T&& t)
          : data_(std::forward<T>(t))
        {}
 
        template<class Tensor,
                 std::enable_if_t<(IsTensor<Tensor>::value
                                   && !Expressions::IsClosure<Tensor>::value
                                   && !std::is_const<DataType>::value
                                   && !std::is_same<Tensor, FieldVectorTensor>::value
                                   && std::is_same<typename Tensor::Signature, Signature>::value
                   ), int> = 0>
        explicit FieldVectorTensor(const Tensor& other)
          : data_(DenseVectorView<const Tensor&>(other))
        {}
 
        template<class T,
                 std::enable_if_t<(IsTensor<T>::value
                                   && !RefersConst<DataType>::value
                                   && std::is_same<Signature, typename TensorTraits<T>::Signature>::value),
                                  int> = 0>
        FieldVectorTensor& operator=(T&& other)
        {
          //std::clog << "Assign from " << typeString(other) << std::endl;
          return assign(*this, std::forward<T>(other), [](auto& dst, auto&& src) { dst = src; });
        }
 
        FieldVectorTensor& operator=(const ArrayType& other)
        {
          data_ = other;
          return *this;
        }
 
        template<class... Dims,
                 std::enable_if_t<(sizeof...(Dims) == rank
                                   &&
                                   IsIntegralPack<Dims...>::value), int> = 0>
        decltype(auto) operator()(Dims... indices)
        {
          return value(data_, indices...);
        }
 
        template<class... Dims,
                 std::enable_if_t<(sizeof...(Dims) == rank
                                   &&
                                   IsIntegralPack<Dims...>::value), int> = 0>
        const FieldType& operator()(Dims... indices) const
        {
          return value(data_, indices...);
        }
 
        template<std::size_t... Indices, std::enable_if_t<sizeof...(Indices) == rank, int> = 0>
        decltype(auto) constexpr operator()(Seq<Indices...>)
        {
          return (*this)(Indices...);
        }
 
        template<std::size_t... Indices, std::enable_if_t<sizeof...(Indices) == rank, int> = 0>
        decltype(auto) constexpr operator()(Seq<Indices...>) const
        {
          return (*this)(Indices...);
        }
 
        operator DataType () &&
        {
          return data_;
        }
 
        operator ConditionalType<RefersConst<DataType>::value, DataType, ArrayType>& () &
        {
          return data_;
        }
 
        operator const ArrayType& () const&
        {
          return data_;
        }
 
        std::string name() const
        {
          return "FT(" + refPrefix() + name<Signature>() + ")";
        }
 
      private:
        std::string refPrefix() const
        {
          constexpr bool isRef = std::is_lvalue_reference<DataType>::value;
          constexpr bool isConst = std::is_const<std::remove_reference_t<DataType> >::value;
 
          return isRef ? (isConst ? "cref" : "ref") : "";
        }
 
        template<class Sig, std::enable_if_t<Sig::size() == 0, int> = 0>
        std::string name() const
        {
          return toString(data_);
        }
 
        template<class Sig, std::enable_if_t<Sig::size() == 1, int> = 0>
        std::string name() const
        {
          return "[" + toString(data_) + "]";
        }
 
        template<class Sig, std::enable_if_t<(Sig::size() > 1), int> = 0>
        std::string name() const
        {
          return "<" + toString(Sig{}) + ">";
        }
 
        template<class Data>
        FieldType& value(Data& data)
        {
          return data;
        }
        template<class Data>
        const FieldType& value(const Data& data)
        {
          return data;
        }
        template<class Data>
        const FieldType& value(const Data& data) const
        {
          return data;
        }
 
        template<int Dim>
        FieldType& value(FieldVector<FieldType, Dim>& data, const std::size_t& index)
        {
          return data[index];
        }
        template<int Dim>
        const FieldType& value(const FieldVector<FieldType, Dim>& data, const std::size_t& index)
        {
          return data[index];
        }
        template<int Dim>
        const FieldType& value(const FieldVector<FieldType, Dim>& data, const std::size_t& index) const
        {
          return data[index];
        }
 
        template<int rows, int cols>
        FieldType& value(FieldMatrix<FieldType, rows, cols>& data, const std::size_t& i, const std::size_t& j)
        {
          return data[i][j];
        }
        template<int rows, int cols>
        const FieldType& value(const FieldMatrix<FieldType, rows, cols>& data, const std::size_t& i, const std::size_t& j)
        {
          return data[i][j];
        }
        template<int rows, int cols>
        const FieldType& value(const FieldMatrix<FieldType, rows, cols>& data, const std::size_t& i, const std::size_t& j) const
        {
          return data[i][j];
        }
 
        template<class Data, class... Dims>
        FieldType& value(Data& data, const std::size_t& index0, Dims... rest)
        {
          return value(data[index0], rest...);
        }
        template<class Data, class... Dims>
        const FieldType& value(const Data& data, const std::size_t& index0, Dims... rest)
        {
          return value(data[index0], rest...);
        }
        template<class Data, class... Dims>
        const FieldType& value(const Data& data, const std::size_t& index0, Dims... rest) const
        {
          return value(data[index0], rest...);
        }
 
        DataType data_;
      };
 
      template<class Field, std::size_t... Dimensions>
      using FieldTensor = FieldVectorTensor<FieldVectorType<Field, Dimensions...> >;
 
      namespace {
 
        template<class Field, class Seq>
        struct MatchingFieldTensorHelper;
 
        template<class Field, std::size_t... Dimensions>
        struct MatchingFieldTensorHelper<Field, Seq<Dimensions...> >
        {
          using Type = FieldTensor<Field, Dimensions...>;
        };
      }
 
      template<class T>
      using MatchingFieldTensor = typename MatchingFieldTensorHelper<FloatingPointClosure<typename TensorTraits<T>::FieldType>, typename TensorTraits<T>::Signature>::Type;
 
      namespace {
 
        template<class FV>
        struct ScalarDecayFieldTensorHelper
        {
          using Type = FieldVectorTensor<FV>;
        };
 
        template<class Data>
        struct ScalarDecayFieldTensorHelper<FieldVector<Data, 1> >
        {
          using Type = FieldVectorTensor<Data>;
        };
 
        template<class Data>
        struct ScalarDecayFieldTensorHelper<Fem::ExplicitFieldVector<Data, 1> >
        {
          using Type = FieldVectorTensor<Data>;
        };
 
        template<int NCols, class Data>
        struct ScalarDecayFieldTensorHelper<FieldMatrix<Data, 1, NCols> >
        {
          using Type = FieldVectorTensor<FieldVector<Data, NCols> >;
        };
 
      }
 
      template<class FV>
      using ScalarDecayFieldTensor = typename ScalarDecayFieldTensorHelper<FV>::Type;
 
      template<class T, class SFINAE = void>
      struct IsFieldVectorTensor
        : FalseType
      {};
 
      template<class T>
      struct IsFieldVectorTensor<T, std::enable_if_t<!IsDecay<T>::value> >
        : IsFieldVectorTensor<std::decay_t<T> >
      {};
 
      template<class FV>
      struct IsFieldVectorTensor<FieldVectorTensor<FV> >
        : TrueType
      {};
 
      template<class T>
      using IsFieldTensor = IsFieldVectorTensor<T>;
 
      template<class FV>
      struct TensorTraits<FV, std::enable_if_t<(IsFieldObject<FV>::value && IsDecay<FV>::value)> >
        : public TensorTraits<FieldVectorTensor<FV> >
      {
        template<class T>
        static constexpr auto toTensor(T&& t)
        {
          return FieldVectorTensor<T>(std::forward<T>(t));
        }
      };
 
      template<class FV>
      struct IsTensorOperand<FV, std::enable_if_t<(IsFieldObject<FV>::value && IsDecay<FV>::value)> >
        : TrueType
      {};
 
      template<class T>
      struct IsUnaryTensorOperand;
 
      template<class T1, class T2, class SFINAE>
      struct AreBinaryTensorOperands;
 
      template<class T1, class T2>
      struct AreBinaryTensorOperands<
        T1, T2,
        std::enable_if_t<(IsFieldObject<T1>::value
                          && IsFieldObject<T2>::value
                          && IsUnaryTensorOperand<T1>::value
                          && IsUnaryTensorOperand<T2>::value
        )> >
        : TrueType
      {};
 
    } // NS Tensor
 
  } // NS ACFem
 
  template<class FV>
  struct FieldTraits<ACFem::Tensor::FieldVectorTensor<FV> >
    : FieldTraits<std::decay_t<FV> >
  {};
 
} // NS Dune
 
#endif // __DUNE_ACFEM_TENSORS_BINDINGS_DUNE_FIELDTENSOR_HH__