high level operation works. TODO: shift code to lib/include files

This commit is contained in:
Christian Zimmermann 2020-09-02 23:55:08 +02:00
parent 705f6351a5
commit 4beb352c0e
2 changed files with 122 additions and 95 deletions

View file

@ -280,10 +280,7 @@ namespace
auto di4 = getIndex(dr4); auto di4 = getIndex(dr4);
auto di4a = getIndex(dr4a); auto di4a = getIndex(dr4a);
//(*di1)({imap["i2_1"],imap["i2_2"],imap["i3_1"],imap["i4_1"]}); (*di1)({"ia_1","ia_2","ib_1","ic_1"});
//(*di2)({imap["i3_1"],imap["i3_1"],imap["i4_2"]});
//(*di4)({imap["i2_1"],imap["i3_1"],imap["i4_1"],imap["i4_2"]});
(*di1)({"ia_1","ia_2","ib_1","ic_1"});
(*di2)({"ib_1","ib_1","ic_2"}); (*di2)({"ib_1","ib_1","ic_2"});
(*di4)({"ia_1","ib_1","ic_1","ic_2"}); (*di4)({"ia_1","ib_1","ic_1","ic_2"});
(*di4a)(svec({"ia_1","ib_1"})); (*di4a)(svec({"ia_1","ib_1"}));
@ -299,7 +296,6 @@ namespace
resx2(i1,di4) = mkDynOp(ma1(i1,di1) * exp(ma2(i1,di2))); resx2(i1,di4) = mkDynOp(ma1(i1,di1) * exp(ma2(i1,di2)));
resx3(i1,di4) = mkDynOp(mkDynOp(ma1(i1,di1)) * mkDynOp(exp(mkDynOp(ma2(i1,di2))))); resx3(i1,di4) = mkDynOp(mkDynOp(ma1(i1,di1)) * mkDynOp(exp(mkDynOp(ma2(i1,di2)))));
//auto xx = std::make_shared<decltype(resx4)>(resx4);
auto xx = mkArrayPtr<double>(nullr()); auto xx = mkArrayPtr<double>(nullr());
auto mi = mkMIndex(i1,di4a); auto mi = mkMIndex(i1,di4a);
@ -307,19 +303,8 @@ namespace
auto hop3 = mkHLO(ma1(i1,di1)); auto hop3 = mkHLO(ma1(i1,di1));
auto hop2 = exp(hop1); auto hop2 = exp(hop1);
auto hop4 = hop3 * hop2; auto hop4 = hop3 * hop2;
auto hopr = mkHLO(resx4(i1,di4));
auto opr = resx4(i1,di4); hopr.assign( hop4, mi, ic_1, ic_2 );
auto loop = mkPILoop
( [&opr,&hop4,&xx,&ic_1,&ic_2,this](){
auto hop4x = hop4;
auto dop2 = hop4x.create(ic_1,ic_2);
auto gexp = mkDynOp1<size_t>(mkMOp<size_t>(dop2.outer,dop2.op));
auto xloop = mkILoop(std::make_tuple(*dop2.op.data()->mOp), std::make_tuple(ic_1, ic_2),
std::make_tuple(xx),
std::make_tuple(opr.assign( *dop2.op.data()->mOp, mkMIndex(ic_1, ic_2) )),
std::array<size_t,1>({1}), std::array<size_t,1>({0}));
return mkGetExpr(gexp, xloop); });
mi->pifor(1,loop)();
auto i2_1 = imap.at("i2_1"); auto i2_1 = imap.at("i2_1");
auto i2_2 = imap.at("i2_2"); auto i2_2 = imap.at("i2_2");

View file

@ -13,6 +13,7 @@ DynamicO<T> mkDynOp1(const Op& op)
return DynamicO<T>(op); return DynamicO<T>(op);
} }
template <class ROP>
class HighLevelOpBase class HighLevelOpBase
{ {
public: public:
@ -66,65 +67,37 @@ public:
virtual RetT<CI,CI> create(const std::shared_ptr<CI> ind1, virtual RetT<CI,CI> create(const std::shared_ptr<CI> ind1,
const std::shared_ptr<CI> ind2) = 0; const std::shared_ptr<CI> ind2) = 0;
virtual const OperationRoot<double,DR>* get1() const = 0; //virtual const OperationRoot<double,CR,DR>* get2() const = 0;
virtual const OperationRoot<double,CR,DR>* get2() const = 0; virtual const ROP* get() const = 0;
}; };
template <bool SAME> template <class ROP>
struct Fwd class HighLevelOpRoot : public HighLevelOpBase<ROP>
{
template <class O1, class O2>
static inline const O1* fwd(const O2* in)
{
assert(0);
return nullptr;
}
};
template <>
struct Fwd<true>
{
template <class O1, class O2>
static inline const O1* fwd(const O2* in)
{
return in;
}
};
template <class OR>
class HighLevelOpRoot : public HighLevelOpBase
{ {
private: private:
typedef OperationRoot<double,DR> OType1; typedef HighLevelOpBase<ROP> B;
typedef OperationRoot<double,CR,DR> OType2;
typedef HighLevelOpBase B;
OR mOp; ROP mOp;
public: public:
HighLevelOpRoot(const OR& op) : mOp(op) {} HighLevelOpRoot(const ROP& op) : mOp(op) {}
virtual bool root() const override final virtual bool root() const override final
{ {
return true; return true;
} }
virtual B::RetT<CI,CI> create(const std::shared_ptr<CI> ind1, virtual typename B::RetT<CI,CI> create(const std::shared_ptr<CI> ind1,
const std::shared_ptr<CI> ind2) override final const std::shared_ptr<CI> ind2) override final
{ {
assert(0); assert(0);
return B::RetT<CI,CI>(); return typename B::RetT<CI,CI>();
} }
virtual const OType1* get1() const override final virtual const ROP* get() const override final
{ {
return Fwd<std::is_same<OType1,OR>::value>::template fwd<OType1>(&mOp); return &mOp;
}
virtual const OType2* get2() const override final
{
return Fwd<std::is_same<OType2,OR>::value>::template fwd<OType2>(&mOp);
} }
@ -143,13 +116,13 @@ struct Create
template <class... Indices> template <class... Indices>
struct cx struct cx
{ {
template <class OpF, class... OPs> template <class ROP, class OpF, class... OPs>
struct ccx struct ccx
{ {
template <size_t M, class... DOPs> template <size_t M, class... DOPs>
static inline void static inline void
cccx(HighLevelOpBase::RetT<Indices...>& res, cccx(typename HighLevelOpBase<ROP>::RetT<Indices...>& res,
const std::array<std::shared_ptr<HighLevelOpBase>,M>& in, const std::array<std::shared_ptr<HighLevelOpBase<ROP>>,M>& in,
const std::shared_ptr<Indices>&... inds, const std::shared_ptr<Indices>&... inds,
const OPs&... ops, const OPs&... ops,
const DOPs&... dops) const DOPs&... dops)
@ -161,13 +134,13 @@ struct Create
auto op = *dop.op.data()->mOp; auto op = *dop.op.data()->mOp;
typedef decltype(op) OP; typedef decltype(op) OP;
res.appendOuter(dop); res.appendOuter(dop);
Create<N-1>::template cx<Indices...>::template ccx<OpF,OP,OPs...>::template cccx<M> Create<N-1>::template cx<Indices...>::template ccx<ROP,OpF,OP,OPs...>::template cccx<M>
(res, in, inds..., op, ops..., dop, dops...); (res, in, inds..., op, ops..., dop, dops...);
} }
else { else {
auto& op = *inn->get2(); auto& op = *inn->get();
typedef typename std::remove_reference<decltype(op)>::type OP; typedef typename std::remove_reference<decltype(op)>::type OP;
Create<N-1>::template cx<Indices...>::template ccx<OpF,OP,OPs...>::template cccx<M> Create<N-1>::template cx<Indices...>::template ccx<ROP,OpF,OP,OPs...>::template cccx<M>
(res, in, inds..., op, ops..., dops...); (res, in, inds..., op, ops..., dops...);
} }
} }
@ -181,13 +154,13 @@ struct Create<0>
template <class... Indices> template <class... Indices>
struct cx struct cx
{ {
template <class OpF, class... OPs> template <class ROP, class OpF, class... OPs>
struct ccx struct ccx
{ {
template <size_t M, class... DOPs> template <size_t M, class... DOPs>
static inline void static inline void
cccx(HighLevelOpBase::RetT<Indices...>& res, cccx(typename HighLevelOpBase<ROP>::RetT<Indices...>& res,
const std::array<std::shared_ptr<HighLevelOpBase>,M>& in, const std::array<std::shared_ptr<HighLevelOpBase<ROP>>,M>& in,
const std::shared_ptr<Indices>&... inds, const std::shared_ptr<Indices>&... inds,
const OPs&... ops, const OPs&... ops,
const DOPs&... dops) const DOPs&... dops)
@ -201,7 +174,7 @@ struct Create<0>
res.appendOuterM(dop.op,dops.op...); res.appendOuterM(dop.op,dops.op...);
} }
else { else {
auto& op = *inn->get2(); auto& op = *inn->get();
res.op = mkDynOutOp(mkFOp<OpF>(op,ops...), inds...); res.op = mkDynOutOp(mkFOp<OpF>(op,ops...), inds...);
res.appendOuterM(dops.op...); res.appendOuterM(dops.op...);
} }
@ -210,42 +183,41 @@ struct Create<0>
}; };
}; };
template <class OpF, size_t N> template <class ROP, class OpF, size_t N>
class HighLevelOp : public HighLevelOpBase class HighLevelOp : public HighLevelOpBase<ROP>
{ {
private: private:
std::array<std::shared_ptr<HighLevelOpBase>,N> mIn; std::array<std::shared_ptr<HighLevelOpBase<ROP>>,N> mIn;
public: public:
typedef HighLevelOpBase B; typedef HighLevelOpBase<ROP> B;
HighLevelOp(std::array<std::shared_ptr<HighLevelOpBase>,N> in) : mIn(in) {} HighLevelOp(std::array<std::shared_ptr<HighLevelOpBase<ROP>>,N> in) : mIn(in) {}
virtual bool root() const override final virtual bool root() const override final
{ {
return false; return false;
} }
virtual const OperationRoot<double,DR>* get1() const override final virtual const ROP* get() const override final
{ assert(0); return nullptr; }
virtual const OperationRoot<double,CR,DR>* get2() const override final
{ assert(0); return nullptr; } { assert(0); return nullptr; }
virtual B::RetT<CI,CI> create(const std::shared_ptr<CI> ind1, virtual typename B::RetT<CI,CI> create(const std::shared_ptr<CI> ind1,
const std::shared_ptr<CI> ind2) override final const std::shared_ptr<CI> ind2) override final
{ {
B::RetT<CI,CI> res; typename B::RetT<CI,CI> res;
Create<N-1>::template cx<CI,CI>::template ccx<OpF>::template cccx<N>(res,mIn,ind1,ind2); Create<N-1>::template cx<CI,CI>::template ccx<ROP,OpF>::template cccx<N>(res,mIn,ind1,ind2);
return res; return res;
} }
}; };
template <class ROP>
class HighLevelOpHolder class HighLevelOpHolder
{ {
private: private:
std::shared_ptr<HighLevelOpBase> mOp; std::shared_ptr<HighLevelOpBase<ROP>> mOp;
public: public:
HighLevelOpHolder() = default; HighLevelOpHolder() = default;
@ -254,34 +226,104 @@ public:
HighLevelOpHolder& operator=(const HighLevelOpHolder& in) = default; HighLevelOpHolder& operator=(const HighLevelOpHolder& in) = default;
HighLevelOpHolder& operator=(HighLevelOpHolder&& in) = default; HighLevelOpHolder& operator=(HighLevelOpHolder&& in) = default;
HighLevelOpHolder(const std::shared_ptr<HighLevelOpBase>& op) : mOp(op) {} HighLevelOpHolder(const std::shared_ptr<HighLevelOpBase<ROP>>& op) : mOp(op) {}
bool root() const { return mOp->root(); } bool root() const { return mOp->root(); }
auto create(const std::shared_ptr<CI> ind1, auto create(const std::shared_ptr<CI> ind1,
const std::shared_ptr<CI> ind2) const { return mOp->create(ind1,ind2); } const std::shared_ptr<CI> ind2) const { return mOp->create(ind1,ind2); }
auto get1() const { return mOp->get1(); } auto get() const { return mOp->get(); }
auto get2() const { return mOp->get2(); }
std::shared_ptr<HighLevelOpBase> op() const { return mOp; } std::shared_ptr<HighLevelOpBase<ROP>> op() const { return mOp; }
HighLevelOpHolder operator*(const HighLevelOpHolder in) const HighLevelOpHolder operator*(const HighLevelOpHolder& in) const
{ {
return HighLevelOpHolder return HighLevelOpHolder<ROP>
( std::make_shared<HighLevelOp<multipliesx<double,double>,2>> ( std::make_shared<HighLevelOp<ROP,multipliesx<double,double>,2>>
( std::array<std::shared_ptr<HighLevelOpBase>,2>({mOp, in.mOp}) ) ); ( std::array<std::shared_ptr<HighLevelOpBase<ROP>>,2>({mOp, in.mOp}) ) );
}
HighLevelOpHolder operator+(const HighLevelOpHolder& in) const
{
return HighLevelOpHolder<ROP>
( std::make_shared<HighLevelOp<ROP,plusx<double,double>,2>>
( std::array<std::shared_ptr<HighLevelOpBase<ROP>>,2>({mOp, in.mOp}) ) );
}
HighLevelOpHolder operator-(const HighLevelOpHolder& in) const
{
return HighLevelOpHolder<ROP>
( std::make_shared<HighLevelOp<ROP,minusx<double,double>,2>>
( std::array<std::shared_ptr<HighLevelOpBase<ROP>>,2>({mOp, in.mOp}) ) );
}
HighLevelOpHolder operator/(const HighLevelOpHolder& in) const
{
return HighLevelOpHolder<ROP>
( std::make_shared<HighLevelOp<ROP,dividesx<double,double>,2>>
( std::array<std::shared_ptr<HighLevelOpBase<ROP>>,2>({mOp, in.mOp}) ) );
}
template <class MIndex, class... Indices>
HighLevelOpHolder& assign(const HighLevelOpHolder& in,
const std::shared_ptr<MIndex>& mi,
const std::shared_ptr<Indices>&... inds)
{
auto xx = mkArrayPtr<double>(nullr());
auto& opr = *mOp->get();
auto loop = mkPILoop
( [&opr,&in,&xx,&inds...,this](){
auto inx = in;
auto dop = inx.create(inds...);
auto gexp = mkDynOp1<size_t>(mkMOp<size_t>(dop.outer,dop.op));
auto xloop = mkILoop(std::make_tuple(*dop.op.data()->mOp),
std::make_tuple(inds...),
std::make_tuple(xx),
std::make_tuple(opr.assign( *dop.op.data()->mOp,
mkMIndex(inds...) )),
std::array<size_t,1>({1}), std::array<size_t,1>({0}));
return mkGetExpr(gexp, xloop); });
mi->pifor(1,loop)();
return *this;
}
template <class MIndex, class... Indices>
HighLevelOpHolder& plus(const HighLevelOpHolder& in,
const std::shared_ptr<MIndex>& mi,
const std::shared_ptr<Indices>&... inds)
{
auto xx = mkArrayPtr<double>(nullr());
auto& opr = *mOp->get();
auto loop = mkPILoop
( [&opr,&in,&xx,&inds...,this](){
auto inx = in;
auto dop = inx.create(inds...);
auto gexp = mkDynOp1<size_t>(mkMOp<size_t>(dop.outer,dop.op));
auto xloop = mkILoop(std::make_tuple(*dop.op.data()->mOp),
std::make_tuple(inds...),
std::make_tuple(xx),
std::make_tuple(opr.plus( *dop.op.data()->mOp,
mkMIndex(inds...) )),
std::array<size_t,1>({1}), std::array<size_t,1>({0}));
return mkGetExpr(gexp, xloop); });
mi->pifor(1,loop)();
return *this;
} }
}; };
template <class OR> typedef OperationRoot<double,CR,DR> OpType1;
HighLevelOpHolder mkHLO(const OR& op)
template <class ROP>
HighLevelOpHolder<ROP> mkHLO(const ROP& op)
{ {
return HighLevelOpHolder(std::make_shared<HighLevelOpRoot<OR>>( op ) ); return HighLevelOpHolder<ROP>(std::make_shared<HighLevelOpRoot<ROP>>( op ) );
} }
HighLevelOpHolder exp(const HighLevelOpHolder& in) template <class ROP>
HighLevelOpHolder<ROP> exp(const HighLevelOpHolder<ROP>& in)
{ {
return HighLevelOpHolder( std::make_shared<HighLevelOp<x_exp<double>,1>> return HighLevelOpHolder<ROP>
( std::array<std::shared_ptr<HighLevelOpBase>,1>( {in.op()} ) ) ); ( std::make_shared<HighLevelOp<ROP,x_exp<double>,1>>
( std::array<std::shared_ptr<HighLevelOpBase<ROP>>,1>( {in.op()} ) ) );
} }