// -*- C++ -*- #ifndef __multi_array_operation_h__ #define __multi_array_operation_h__ #include #include #include #include #include #include "base_def.h" #include "mbase_def.h" #include "ranges/rheader.h" #include "pack_num.h" #include "arith.h" namespace MultiArrayTools { namespace { using namespace MultiArrayHelper; } template class OperationBase { public: OperationClass& THIS() { return static_cast(*this); } const OperationClass& THIS() const { return static_cast(*this); } template auto operator+(const Second& in) const -> Operation,OperationClass,Second>; template auto operator-(const Second& in) const -> Operation,OperationClass,Second>; template auto operator*(const Second& in) const -> Operation,OperationClass,Second>; template auto operator/(const Second& in) const -> Operation,OperationClass,Second>; template auto c(std::shared_ptr& ind) const -> Contraction; private: friend OperationClass; OperationTemplate() = default; }; template class OperationTemplate : public OperationBase { /* empty per default; specialize if needed */ }; #include "type_operations.h" template class OperationMaster { public: class AssignmentExpr { private: AssignmentExpr() = default; OperationMaster& mM; const OpClass& mSec; public: static constexpr size_t LAYER = 0; static constexpr size_t SIZE = OpClass::SIZE; typedef decltype(mSec.rootSteps()) ExtType; AssignmentExpr(OperationMaster& m, const OpClass& sec); AssignmentExpr(const AssignmentExpr& in) = default; AssignmentExpr(AssignmentExpr&& in) = default; inline void operator()(size_t start = 0) const; inline void operator()(size_t start, ExtType last) const; auto rootSteps(std::intptr_t iPtrNum = 0) const -> ExtType; }; typedef T value_type; typedef OperationBase OB; typedef ContainerRange CRange; typedef ContainerIndex IndexType; //typedef typename MultiRange::IndexType IndexType; OperationMaster(MutableMultiArrayBase& ma, const OpClass& second, IndexType& index); OperationMaster(T* data, const OpClass& second, IndexType& index); inline void set(size_t pos, T val) { mDataPtr[pos] = val; } inline void add(size_t pos, T val) { mDataPtr[pos] += val; } inline T get(size_t pos) const; private: void performAssignment(std::intptr_t blockIndexNum); OpClass const& mSecond; //MutableMultiArrayBase& mArrayRef; T* mDataPtr; IndexType mIndex; }; template class ConstOperationRoot : public OperationTemplate > { public: typedef T value_type; typedef OperationBase > OT; typedef ContainerRange CRange; typedef ContainerIndex IndexType; static constexpr size_t SIZE = 1; ConstOperationRoot(const MultiArrayBase& ma, const std::shared_ptr&... indices); ConstOperationRoot(const T* data, const IndexType& ind); template inline T get(ET pos) const; MExt rootSteps(std::intptr_t iPtrNum = 0) const; // nullptr for simple usage with decltype template Expr loop(Expr exp) const; private: //MultiArrayBase const& mArrayRef; const T* mDataPtr; IndexType mIndex; }; template class OperationRoot : public OperationTemplate > { public: typedef T value_type; typedef OperationBase > OT; typedef ContainerRange CRange; typedef ContainerIndex IndexType; static constexpr size_t SIZE = 1; OperationRoot(MutableMultiArrayBase& ma, const std::shared_ptr&... indices); OperationRoot(T* data, const IndexType& ind); template OperationMaster operator=(const OpClass& in); template inline T get(ET pos) const; MExt rootSteps(std::intptr_t iPtrNum = 0) const; // nullptr for simple usage with decltype template Expr loop(Expr exp) const; T* data() const; private: //MutableMultiArrayBase& mArrayRef; T* mDataPtr; IndexType mIndex; }; template size_t sumRootNum() { return typename Op::rootNum(); } template size_t sumRootNum() { return typename Op1::rootNum() + sumRootNum(); } template struct RootSumN { template struct rs { static constexpr size_t SIZE = Op1::SIZE + RootSumN::template rs::SIZE; }; }; template <> struct RootSumN<0> { template struct rs { static constexpr size_t SIZE = Op1::SIZE; }; }; template struct RootSum { static constexpr size_t SIZE = RootSumN::template rs::SIZE; }; template class Operation : public OperationTemplate > { public: typedef T value_type; typedef OperationBase > OT; typedef OpFunction F; static constexpr size_t SIZE = RootSum::SIZE; static constexpr bool FISSTATIC = OpFunction::FISSTATIC; private: std::tuple mOps; std::shared_ptr mF; // only if non-static public: typedef decltype(PackNum::template mkSteps(0, mOps)) ETuple; Operation(const Ops&... ops); Operation(std::shared_ptr ff, const Ops&... ops); template inline T get(ET pos) const; auto rootSteps(std::intptr_t iPtrNum = 0) const // nullptr for simple usage with decltype -> decltype(PackNum::mkSteps(iPtrNum, mOps)); template auto loop(Expr exp) const -> decltype(PackNum::mkLoop( mOps, exp)); }; template auto mkOperation(const OpFunction& f, const Ops&... ops) -> Operation { return Operation(ops...); } template class Contraction : public OperationTemplate > { public: typedef T value_type; typedef OperationBase > OT; static constexpr size_t SIZE = Op::SIZE; private: const Op& mOp; std::shared_ptr mInd; public: typedef decltype(mOp.rootSteps(0)) ETuple; Contraction(const Op& op, std::shared_ptr ind); template inline T get(ET pos) const; auto rootSteps(std::intptr_t iPtrNum = 0) const // nullptr for simple usage with decltype -> decltype(mOp.rootSteps(iPtrNum)); template auto loop(Expr exp) const -> decltype(mInd->iforh(exp)); }; template struct operate { static inline OperationRoot apply(const MultiArrayBase& ma, const std::shared_ptr&... indices) { return OperationRoot(ma, indices...); } }; // interchange inheritance // try to implement the 'promote' version template class TypeSpecificOperationSet,OperationClass> { auto operator()() const -> Operation,operate,OperationClass> { return Operation,operate,OperationClass> (/*THIS*/); } }; } /* ========================= * * --- TEMPLATE CODE --- * * ========================= */ namespace MultiArrayTools { namespace { using namespace MultiArrayHelper; } /*************************** * OperationTemplate * ***************************/ template template auto OperationTemplate::operator+(const Second& in) const -> Operation,OperationClass,Second> { return Operation,OperationClass,Second>(THIS(), in); } template template auto OperationTemplate::operator-(const Second& in) const -> Operation,OperationClass,Second> { return Operation,OperationClass,Second>(THIS(), in); } template template auto OperationTemplate::operator*(const Second& in) const -> Operation,OperationClass,Second> { return Operation,OperationClass,Second>(THIS(), in); } template template auto OperationTemplate::operator/(const Second& in) const -> Operation,OperationClass,Second> { return Operation,OperationClass,Second>(THIS(), in); } template template auto OperationTemplate::c(std::shared_ptr& ind) const -> Contraction { return Contraction(THIS(), ind); } /***************************************** * OperationMaster::AssignmentExpr * *****************************************/ template OperationMaster::AssignmentExpr:: AssignmentExpr(OperationMaster& m, const OpClass& sec) : mM(m), mSec(sec) {} template inline void OperationMaster::AssignmentExpr:: operator()(size_t start, ExtType last) const { mM.add(start, mSec.template get(last) ); } template typename OperationMaster::AssignmentExpr::ExtType OperationMaster::AssignmentExpr:: rootSteps(std::intptr_t iPtrNum) const { return mSec.rootSteps(iPtrNum); } /************************* * OperationMaster * *************************/ template OperationMaster:: OperationMaster(MutableMultiArrayBase& ma, const OpClass& second, IndexType& index) : mSecond(second), mDataPtr(ma.data()), mIndex(index) { performAssignment(0); } template OperationMaster:: OperationMaster(T* data, const OpClass& second, IndexType& index) : mSecond(second), mDataPtr(data), mIndex(index) { performAssignment(0); } template void OperationMaster::performAssignment(std::intptr_t blockIndexNum) { AssignmentExpr ae(*this, mSecond); // Expression to be executed within loop const auto loop = mSecond.template loop( mIndex.ifor(ae) ); // hidden Loops outside ! -> auto vectorizable loop(); // execute overall loop(s) and so internal hidden loops and so the inherited expressions } template inline T OperationMaster::get(size_t pos) const { return mDataPtr[pos]; } /**************************** * ConstOperationRoot * ****************************/ template ConstOperationRoot:: ConstOperationRoot(const MultiArrayBase& ma, const std::shared_ptr&... indices) : mDataPtr(ma.data()), mIndex( ma.begin() ) { mIndex(indices...); } template ConstOperationRoot:: ConstOperationRoot(const T* data, const IndexType& ind) : mDataPtr(data), mIndex( ind ) { } template template inline T ConstOperationRoot::get(ET pos) const { return mDataPtr[pos.val()]; } template MExt ConstOperationRoot::rootSteps(std::intptr_t iPtrNum) const { return MExt(getStepSize( mIndex, iPtrNum )); //return MExt(getStepSize( getRootIndices( mIndex->info() ), iPtrNum )); } template template Expr ConstOperationRoot::loop(Expr exp) const { return exp; } /*********************** * OperationRoot * ***********************/ template OperationRoot:: OperationRoot(MutableMultiArrayBase& ma, const std::shared_ptr&... indices) : mDataPtr(ma.data()), mIndex( ma.begin() ) { mIndex(indices...); } template OperationRoot:: OperationRoot(T* data, const IndexType& ind) : mDataPtr(data), mIndex( ind ) { } template template OperationMaster OperationRoot::operator=(const OpClass& in) { return OperationMaster(mDataPtr, in, mIndex); } template template inline T OperationRoot::get(ET pos) const { return mDataPtr[pos.val()]; } template MExt OperationRoot::rootSteps(std::intptr_t iPtrNum) const { return MExt(getStepSize( mIndex, iPtrNum )); //return MExt(getStepSize( mIndex.info(), iPtrNum )); } template template Expr OperationRoot::loop(Expr exp) const { return exp; } template T* OperationRoot::data() const { return mDataPtr + mIndex.pos(); } /******************* * Operation * *******************/ template Operation::Operation(const Ops&... ops) : mOps(ops...) { static_assert( FISSTATIC, "need function instance for non-static function" ); } template Operation::Operation(std::shared_ptr ff, const Ops&... ops) : mOps(ops...), mF(ff) { static_assert( not FISSTATIC, "using instance of function supposed to be static" ); } template template inline T Operation::get(ET pos) const { typedef std::tuple OpTuple; return PackNum:: template mkOpExpr(mF, pos, mOps); } template auto Operation::rootSteps(std::intptr_t iPtrNum) const -> decltype(PackNum::mkSteps(iPtrNum, mOps)) { return PackNum::mkSteps(iPtrNum, mOps); } template template auto Operation::loop(Expr exp) const -> decltype(PackNum::mkLoop( mOps, exp )) { return PackNum::mkLoop( mOps, exp ); } /********************* * Contraction * *********************/ template Contraction::Contraction(const Op& op, std::shared_ptr ind) : mOp(op), mInd(ind) {} // forward loop !!!! template template inline T Contraction::get(ET pos) const { return mOp.template get(pos); } template auto Contraction::rootSteps(std::intptr_t iPtrNum) const -> decltype(mOp.rootSteps(iPtrNum)) { return mOp.rootSteps(iPtrNum); } template template auto Contraction::loop(Expr exp) const -> decltype(mInd->iforh(exp)) { return mInd->iforh(exp); } } #endif