restriction.hh

#ifndef __DUNE_ACFEM_TENSORS_OPTIMIZATION_RESTRICTION_HH__
#define __DUNE_ACFEM_TENSORS_OPTIMIZATION_RESTRICTION_HH__
 
#include "../../expressions/optimizegeneral.hh"
#include "../expressionoperations.hh"
#include "../modules.hh"
#include "../operations/restriction.hh"
 
namespace Dune {
 
  namespace ACFem {
 
    namespace Expressions {
 
      template<class DefectPositions, class PivotSequence, class T>
      constexpr inline bool IsUnaryZeroV<
        OperationTraits<RestrictionOperation<DefectPositions, PivotSequence> >,
        T> = (DefectPositions::size() == PivotSequence::size() // only compile-time constant restrictions
              && TensorTraits<T>::isZero(PivotSequence{}, DefectPositions{}));
 
    }
 
    namespace Tensor {
 
      namespace Optimization::Restriction {
 
        using RestrictionTag = OptimizeNext<Policy::DefaultOptimizationTag>;
 
        template<class DefectPositions, class PivotSequence, class T>
        constexpr auto zero(OperationTraits<RestrictionOperation<DefectPositions, PivotSequence> >, T&&)
        {
          using Signature = typename TensorTraits<T>::Signature;
          using DefectSignature = SequenceSliceComplement<Signature, DefectPositions>;
 
          DUNE_ACFEM_EXPRESSION_RESULT(
            zeros(DefectSignature{}, Disclosure{})
            , "zero restriction"
            );
        }
 
        template<class F, class T>
        constexpr inline bool IsConstRestrictionOfConstRestrictionV = false;
 
        template<class DefectPositions, class PivotSequence, class T>
        constexpr inline bool IsConstRestrictionOfConstRestrictionV<
          OperationTraits<RestrictionOperation<DefectPositions, PivotSequence> >,
          T> = (DefectPositions::size() == PivotSequence::size()
                && IsConstRestriction<T>::value);
 
        template<class F, class T, std::enable_if_t<IsConstRestrictionOfConstRestrictionV<std::decay_t<F>, T>, int> = 0>
        constexpr decltype(auto) operate(RestrictionTag, F&&, T&& t)
        {
          DUNE_ACFEM_RECORD_OPTIMIZATION;
 
          using OuterTraits = RestrictionTraits<std::decay_t<F> >;
          using InnerTraits = RestrictionTraits<Operation<T> >;
          using OuterDefects = typename OuterTraits::DefectPositions;
          using OuterPivots = typename OuterTraits::PivotSequence;
          using InnerDefects = typename InnerTraits::DefectPositions;
          using InnerPivots = typename InnerTraits::PivotSequence;
 
          using MergedDefects = JoinedDefects<InnerDefects, OuterDefects>;
          using Injects = JoinedInjections<InnerDefects, OuterDefects>;
          using MergedPivots = InsertAt<InnerPivots, OuterPivots, Injects>;
 
          DUNE_ACFEM_EXPRESSION_RESULT(
            (operate<RestrictionOperation<MergedDefects, MergedPivots> >(
              std::forward<T>(t).operand(0_c)
              ))
             , "merge restrictions"
            );
        }
 
        template<class F, class T>
        constexpr inline bool IsNonConstRestrictionOfNonConstRestrictionV = false;
 
        template<std::size_t D0, std::size_t... DRest, class T>
        constexpr inline bool IsNonConstRestrictionOfNonConstRestrictionV<
          OperationTraits<RestrictionOperation<Seq<D0, DRest...>, Seq<> > >,
          T> = (IsDynamicRestriction<T>::value
                && !std::is_reference<Operand<0, T> >::value);
 
        template<class F, class T, std::enable_if_t<IsNonConstRestrictionOfNonConstRestrictionV<F, T>, int> = 0>
        constexpr decltype(auto) operate(RestrictionTag, F&& f, T&& t)
        {
          DUNE_ACFEM_RECORD_OPTIMIZATION;
 
          using OuterTraits = RestrictionTraits<std::decay_t<F> >;
          using InnerTraits = RestrictionTraits<Operation<T> >;
          using InnerDefects = typename InnerTraits::DefectPositions;
          using OuterDefects = typename OuterTraits::DefectPositions;
 
          using MergedDefects = JoinedDefects<InnerDefects, OuterDefects>;
          using Injects = JoinedInjections<InnerDefects, OuterDefects>;
 
          auto at = insertAt(std::forward<T>(t).lookAt(), f.pivotIndices_, Injects{});
 
          DUNE_ACFEM_EXPRESSION_RESULT(
            operate(OperationTraits<RestrictionOperation<MergedDefects, Seq<> > >(at), std::move(std::forward<T>(t)).operand(0_c))
            , "merge non-const restrictions"
            );
        }
 
        template<class F, class T>
        constexpr inline bool IsConstRestrictionOfEyeV = false;
 
        template<class DefectPositions, class PivotSequence, class T>
        constexpr inline bool IsConstRestrictionOfEyeV<
          OperationTraits<RestrictionOperation<DefectPositions, PivotSequence> >,
          T> = (DefectPositions::size() == PivotSequence::size()
                && IsEye<T>::value
                && isConstant(PivotSequence{}));
 
        template<class F, class T, std::enable_if_t<IsConstRestrictionOfEyeV<F, T>, int> = 0>
        constexpr auto operate(OptimizeTerminal1, F&&, T&& t)
        {
          DUNE_ACFEM_RECORD_OPTIMIZATION;
 
          using Traits = RestrictionTraits<std::decay_t<F> >;
          using DefectPositions = typename Traits::DefectPositions;
          using PivotSequence = typename Traits::PivotSequence;
          using Signature = typename TensorTraits<T>::Signature;
          using DefectSignature = SequenceSliceComplement<Signature, DefectPositions>;
          using ValueType = typename std::decay_t<T>::ValueType;
          constexpr std::size_t defectRank = DefectSignature::size();
 
          DUNE_ACFEM_EXPRESSION_RESULT(
            (KroneckerDelta<DefectSignature, MakeConstantSequence<defectRank, Head<PivotSequence>::value>, ValueType>{})
            , "non-zero restriction of eye"
            );
        }
 
        template<class F, class T>
        constexpr inline bool IsConstRestrictionOfBlockEyeV = false;
 
        template<class DefectPositions, class PivotSequence, class T>
        constexpr inline bool IsConstRestrictionOfBlockEyeV<
          OperationTraits<RestrictionOperation<DefectPositions, PivotSequence> >,
          T> = (DefectPositions::size() == PivotSequence::size()
                && IsBlockEye<T>::value
                && IsBlockEyePositionsV<DefectPositions, BlockEyeTraits<T>::blockRank_, typename BlockEyeTraits<T>::BlockSignature>);
 
        template<class F, class T, std::enable_if_t<IsConstRestrictionOfBlockEyeV<F, T>, int> = 0>
        constexpr auto operate(OptimizeTerminal1, F&&, T&& t)
        {
          DUNE_ACFEM_RECORD_OPTIMIZATION;
 
          using OperandType = std::decay_t<T>;
          using Signature = typename OperandType::Signature;
          using Traits = RestrictionTraits<std::decay_t<F> >;
          using DefectPositions = typename Traits::DefectPositions;
          using PivotSequence = typename Traits::PivotSequence;
          using DefectSignature = SequenceSliceComplement<Signature, DefectPositions>;
 
          using ValueType = typename OperandType::ValueType;
          constexpr std::size_t defectRank = DefectSignature::size();
 
          using BlockSignature = typename BlockEyeTraits<OperandType>::BlockSignature;
          constexpr std::size_t blockSize = BlockSignature::size();
 
          static_assert(defectRank % blockSize == 0,
                        "Internal template specialization error: rank of resulting KroneckerDelta must be a multiple of the BlockSignature size.");
 
          constexpr std::size_t numBlocks = defectRank / blockSize;
 
          // One building block for the Kronecker-delta
          using PivotBlock = HeadPart<blockSize, PivotSequence>;
 
          // Generate the Pivot sequence
          using KroneckerPivots = SequenceProd<numBlocks, PivotBlock>;
 
          DUNE_ACFEM_EXPRESSION_RESULT(
            (KroneckerDelta<DefectSignature, KroneckerPivots, ValueType>{})
            , "non-zero restriction of block-eye"
            );
        }
 
      } // Optimization::Restriction::
 
      using Optimization::Restriction::zero;
 
    } // Tensor::
 
    namespace Expressions {
 
      using Tensor::Optimization::Restriction::operate;
 
    } // Expressions
 
  } // ACFem
 
} // Dune
 
#endif // __DUNE_ACFEM_TENSORS_OPTIMIZATION_RESTRICTION_HH__