add contraction class -> tests work

This commit is contained in:
Christian Zimmermann 2017-03-21 15:18:34 +01:00
parent 7f6762f159
commit a3183f5ced
4 changed files with 234 additions and 42 deletions

View file

@ -113,6 +113,10 @@ namespace MultiArrayTools
template <typename T, class Operation, class... MAOps>
class MultiArrayOperation;
// multi_array_operation.h
template <typename T, class Operation, class Range, class MAOp>
class MultiArrayContraction;
// slice.h
template <typename T, class Range, class MARange>
class Slice;

View file

@ -34,7 +34,7 @@ namespace MultiArrayTools
template <typename T, class Range>
void MultiArrayOperationRoot<T,Range>::performAssignment(const MultiArrayOperationBase<T>& in)
{
#error "WRITE MAOR INTRINSIC CONTRACT FUNCTION"
//#error "WRITE MAOR INTRINSIC CONTRACT FUNCTION"
//CHECK;
in.linkIndicesTo(MAOB::mIibPtr);
//CHECK;
@ -42,6 +42,7 @@ namespace MultiArrayTools
//CHECK;
const size_t endPos = mArrayRef.endIndex().pos();
std::cout << "assignment: " << endPos << " elements" << std::endl;
// assignment loop
for(iref = mArrayRef.beginIndex().pos(); iref != mArrayRef.endIndex(); ++iref){
std::cout << iref.pos() << '\r' << std::flush;
get() = in.get();
@ -50,6 +51,7 @@ namespace MultiArrayTools
MAOB::mIibPtr->freeLinked();
}
/*
template <typename T, class Range>
template <class RangeX>
MultiArrayOperationRoot<T,Range>&
@ -70,7 +72,7 @@ namespace MultiArrayTools
sl.setConst(in.mArrayRef, name(), dynamic_cast<const typename RangeX::IndexType&>( in.index() ), in.name());
return *this;
}
*/
// CONST SLICE !!!!!
@ -102,10 +104,9 @@ namespace MultiArrayTools
MultiArrayOperationRoot<T,Range>::operator=(MultiArrayOperationRoot<T,Range>& in)
{
//CHECK;
maketurnSlice(in);
if(mArrayRef.isSlice() and not mArrayRef.isInit()){
return makeSlice(in);
}
// if(mArrayRef.isSlice() and not mArrayRef.isInit()){
// return makeSlice(in);
//}
performAssignment(in);
freeIndex();
return *this;
@ -117,9 +118,9 @@ namespace MultiArrayTools
MultiArrayOperationRoot<T,Range>::operator=(MultiArrayOperationRoot<T,Range2>& in)
{
//CHECK;
if(mArrayRef.isSlice() and not mArrayRef.isInit()){
return makeSlice(in);
}
//if(mArrayRef.isSlice() and not mArrayRef.isInit()){
// return makeSlice(in);
//}
performAssignment(in);
freeIndex();
return *this;
@ -132,10 +133,10 @@ namespace MultiArrayTools
MultiArrayOperationRoot<T,Range>::operator=(const MultiArrayOperationRoot<T,Range2>& in)
{
//CHECK;
if(mArrayRef.isSlice() and not mArrayRef.isInit()){
//if(mArrayRef.isSlice() and not mArrayRef.isInit()){
//CHECK;
return makeConstSlice(in);
}
// return makeConstSlice(in);
//}
performAssignment(in);
freeIndex();
return *this;
@ -147,11 +148,11 @@ namespace MultiArrayTools
MultiArrayOperationRoot<T,Range>::operator=(const MultiArrayOperation<T,Operation,MAOps...>& in)
{
//CHECK;
if(mArrayRef.isSlice() and not mArrayRef.isInit()){
//if(mArrayRef.isSlice() and not mArrayRef.isInit()){
// NO SLICE CREATION !!! (total array not initialized!!)
// throw !
assert(0);
}
// assert(0);
//}
performAssignment(in);
freeIndex();
return *this;
@ -166,6 +167,31 @@ namespace MultiArrayTools
return MultiArrayOperation<T,Operation,MultiArrayOperationRoot<T,Range>, MAOps...>(op, *this, secs...);
}
template <typename T, class Range>
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,MultiArrayOperationRoot<T,Range> >
MultiArrayOperationRoot<T,Range>::contract(const ContractOperation& cop,
const std::string& indexName) const
{
typename Range2::IndexType* ind = dynamic_cast<typename Range2::IndexType*>( mIndex.getLinked(indexName) );
return MultiArrayContraction<T,ContractOperation,Range2,
MultiArrayOperationRoot<T,Range> >(cop, *this, *ind);
}
template <typename T, class Range>
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,MultiArrayOperationRoot<T,Range> >
MultiArrayOperationRoot<T,Range>::contract(const ContractOperation& cop,
const std::string& indexName,
const typename Range2::IndexType& begin,
const typename Range2::IndexType& end) const
{
typename Range2::IndexType* ind = dynamic_cast<typename Range2::IndexType*>( mIndex.getLinked(indexName) );
return MultiArrayContraction<T,ContractOperation,Range2,
MultiArrayOperationRoot<T,Range> >(cop, *this, *ind, begin, end);
}
template <typename T, class Range>
template <class MAOp>
auto MultiArrayOperationRoot<T,Range>::operator+(const MAOp& sec)
@ -411,6 +437,31 @@ namespace MultiArrayTools
return MultiArrayOperation<T,Operation,ConstMultiArrayOperationRoot<T,Range>, MAOps...>(op, *this, secs...);
}
template <typename T, class Range>
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,ConstMultiArrayOperationRoot<T,Range> >
ConstMultiArrayOperationRoot<T,Range>::contract(const ContractOperation& cop,
const std::string& indexName) const
{
typename Range2::IndexType* ind = dynamic_cast<typename Range2::IndexType*>( mIndex.getLinked(indexName) );
return MultiArrayContraction<T,ContractOperation,Range2,
ConstMultiArrayOperationRoot<T,Range> >(cop, *this, *ind);
}
template <typename T, class Range>
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,ConstMultiArrayOperationRoot<T,Range> >
ConstMultiArrayOperationRoot<T,Range>::contract(const ContractOperation& cop,
const std::string& indexName,
const typename Range2::IndexType& begin,
const typename Range2::IndexType& end) const
{
typename Range2::IndexType* ind = dynamic_cast<typename Range2::IndexType*>( mIndex.getLinked(indexName) );
return MultiArrayContraction<T,ContractOperation,Range2,
ConstMultiArrayOperationRoot<T,Range> >(cop, *this, *ind, begin, end);
}
template <typename T, class Range>
template <class MAOp>
auto ConstMultiArrayOperationRoot<T,Range>::operator+(const MAOp& sec) const
@ -639,20 +690,55 @@ namespace MultiArrayTools
TupleIndicesLinker<sizeof...(MAOps)-1>::linkTupleIndicesTo(mArgs, target);
}
/*
template <typename T, class Operation, class... MAOps>
T& MultiArrayOperation<T,Operation,MAOps...>::get()
{
mVal = OperationCall<sizeof...(MAOps)-1>::
template callOperation(mOp, mArgs);
return mVal;
}*/
template <typename T, class Operation, class... MAOps>
template <typename T, class Operation, class... MAOps>
const T& MultiArrayOperation<T,Operation,MAOps...>::get() const
{
mVal = OperationCall<sizeof...(MAOps)-1>::
template callOperation(mOp, mArgs);
return mVal;
}
/*******************************
* MultiArrayContraction *
*******************************/
template <typename T, class ContractOperation, class Range, class MAOp>
MultiArrayContraction<T,ContractOperation,Range,MAOp>::
MultiArrayContraction(const ContractOperation& cop, const MAOp& mao,
const typename Range::IndexType& runIndex) :
MultiArrayOperation<T,ContractOperation,MAOp>(cop, mao),
mBeginIndex(runIndex), mEndIndex(runIndex),
mRunIndex(runIndex)
{
mBeginIndex.setPos(0);
mEndIndex.setPos(mRunIndex.max());
MAO::linkIndicesTo(&mRunIndex);
}
template <typename T, class ContractOperation, class Range, class MAOp>
MultiArrayContraction<T,ContractOperation,Range,MAOp>::
MultiArrayContraction(const ContractOperation& cop, const MAOp& mao,
const typename Range::IndexType& runIndex,
const typename Range::IndexType& beginIndex,
const typename Range::IndexType& endIndex) :
MultiArrayOperation<T,ContractOperation,MAOp>(cop, mao),
mBeginIndex(beginIndex), mEndIndex(endIndex),
mRunIndex(runIndex)
{
MAO::linkIndicesTo(&mRunIndex);
}
// for the moment simplest case only:
template <typename T, class ContractOperation, class Range, class MAOp>
const T& MultiArrayContraction<T,ContractOperation,Range,MAOp>::get() const
{
MAO::mOp.reset();
for(mRunIndex.copyPos( mBeginIndex ); mRunIndex.pos() != mEndIndex.pos(); ++mRunIndex){
MAO::mOp(std::get<0>(MAO::mArgs).get() );
}
MAO::mOp.endOp(MAO::mVal);
return MAO::mOp();
}
}

View file

@ -73,6 +73,17 @@ namespace MultiArrayTools
template <class Operation, class... MAOps>
MultiArrayOperation<T,Operation,MultiArrayOperationRoot<T,Range>, MAOps...>
operator()(const Operation& op, const MAOps&... secs) const;
template < class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,MultiArrayOperationRoot<T,Range> >
contract(const ContractOperation& cop, const std::string& indexName) const;
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,MultiArrayOperationRoot<T,Range> >
contract(const ContractOperation& cop, const std::string& indexName,
const typename Range2::IndexType& begin,
const typename Range2::IndexType& end) const;
template <class MAOp>
auto operator+(const MAOp& sec) -> decltype(operator()(std::plus<T>(), sec));
@ -131,13 +142,13 @@ namespace MultiArrayTools
void performAssignment(const MultiArrayOperationBase<T>& in);
/*
template <class RangeX>
MultiArrayOperationRoot& makeSlice(MultiArrayOperationRoot<T,RangeX>& in);
template <class RangeX>
const MultiArrayOperationRoot& makeConstSlice(const MultiArrayOperationRoot<T,RangeX>& in);
*/
MutableMultiArrayBase<T,Range>& mArrayRef;
mutable IndexType mIndex;
Name mNm;
@ -155,23 +166,20 @@ namespace MultiArrayTools
ConstMultiArrayOperationRoot(const MultiArrayBase<T,Range>& ma, const Name& nm);
ConstMultiArrayOperationRoot(const MultiArrayOperationRoot<T,Range>& in);
/*
const ConstMultiArrayOperationRoot& operator=(const ConstMultiArrayOperationRoot& in);
template <class Range2>
const ConstMultiArrayOperationRoot& operator=(const ConstMultiArrayOperationRoot<T,Range2>& in);
template <class Range2>
const ConstMultiArrayOperationRoot& operator=(const MultiArrayOperationRoot<T,Range2>& in);
*/
//template <class Operation, class... MAOps>
//MultiArrayOperation<T,Operation,MultiArrayOperationRoot<T,Range>, MAOps...>
//operator()(Operation& op, const MAOps&... secs) const;
template <class Operation, class... MAOps>
MultiArrayOperation<T,Operation,ConstMultiArrayOperationRoot<T,Range>, MAOps...>
operator()(const Operation& op, const MAOps&... secs) const;
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,ConstMultiArrayOperationRoot<T,Range> >
contract(const ContractOperation& cop, const std::string& indexName) const;
template <class Range2, class ContractOperation>
MultiArrayContraction<T,ContractOperation,Range2,ConstMultiArrayOperationRoot<T,Range> >
contract(const ContractOperation& cop, const std::string& indexName,
const typename Range2::IndexType& begin,
const typename Range2::IndexType& end) const;
template <class MAOp>
auto operator+(const MAOp& sec) const -> decltype(operator()(std::plus<T>(), sec));
@ -262,7 +270,6 @@ namespace MultiArrayTools
virtual void linkIndicesTo(IndefinitIndexBase* target) const override;
//virtual T& get() override;
virtual const T& get() const override;
protected:
@ -272,8 +279,30 @@ namespace MultiArrayTools
OBT mArgs; // include first arg also here !!!
};
template <typename T, class ContractOperation, class Range, class MAOp>
class MultiArrayContraction : public MultiArrayOperation<T,ContractOperation,MAOp>
{
public:
typedef MultiArrayOperationBase<T> MAOB;
typedef MultiArrayOperation<T,ContractOperation,MAOp> MAO;
MultiArrayContraction(const ContractOperation& cop, const MAOp& mao,
const typename Range::IndexType& runIndex);
MultiArrayContraction(const ContractOperation& cop, const MAOp& mao,
const typename Range::IndexType& runIndex,
const typename Range::IndexType& beginIndex,
const typename Range::IndexType& endIndex);
virtual const T& get() const override;
protected:
typename Range::IndexType mBeginIndex;
typename Range::IndexType mEndIndex;
mutable typename Range::IndexType mRunIndex;
};
}
#include "multi_array_operation.cc"

View file

@ -9,6 +9,33 @@
namespace MAT = MultiArrayTools;
namespace {
template <typename T>
struct sum
{
public:
sum() = default;
T& operator()() const
{
return res;
}
T& operator()(const T& a) const
{
return res += a;
}
void endOp(T& res) const {}
void reset() const
{
res = static_cast<T>(0);
}
private:
mutable T res = static_cast<T>(0);
};
class OneDimTest : public ::testing::Test
{
@ -135,6 +162,32 @@ namespace {
MultiArray3dAny ma;
//Slice2d3dAny sl;
};
class ContractionTest : public ::testing::Test
{
protected:
typedef MAT::SingleRange<char,MAT::RangeType::ANY> Range1dAny;
typedef MAT::MultiRange<Range1dAny,Range1dAny> Range2dAny;
typedef MAT::MultiRange<Range1dAny,Range1dAny,Range1dAny> Range3dAny;
typedef MAT::MultiArray<int,Range3dAny> MultiArray3dAny;
typedef MAT::Slice<int,Range2dAny,Range3dAny> Slice2d3dAny;
typedef MAT::MultiArray<int,Range2dAny> MultiArray2dAny;
ContractionTest() : r1({'a','b','c'}), r2({'a','b','c','d'}), r3({'a','b'}),
ra(r1,r3),
rb(r1,r2),
r3d(r1,r2,r3),
ma(r3d, {-5,6,2,1,9,54,27,-7,-13,32,90,-67,
-10,16,-2,101,39,-64,81,-22,14,34,95,-62}) {}
Range1dAny r1;
Range1dAny r2;
Range1dAny r3;
Range2dAny ra;
Range2dAny rb;
Range3dAny r3d;
MultiArray3dAny ma;
};
TEST_F(OneDimTest, CorrectExtensions)
{
@ -454,6 +507,26 @@ namespace {
EXPECT_EQ(sl[j(j1 = 2, j2 = 0)], 14);
EXPECT_EQ(sl[j(j1 = 2, j2 = 1)], 34);
}
TEST_F(ContractionTest, ContractionWorks)
{
MultiArray2dAny ma2(ra);
ma2("alpha","gamma") = ma("alpha","beta","gamma").contract<Range1dAny>(sum<int>(),"beta");
auto i = ma2.beginIndex();
auto i1 = i.template getIndex<0>();
auto i2 = i.template getIndex<1>();
EXPECT_EQ(ma2[i(i1 = 0, i2 = 0)], 33);
EXPECT_EQ(ma2[i(i1 = 0, i2 = 1)], 54);
EXPECT_EQ(ma2[i(i1 = 1, i2 = 0)], 65);
EXPECT_EQ(ma2[i(i1 = 1, i2 = 1)], 82);
EXPECT_EQ(ma2[i(i1 = 2, i2 = 0)], 229);
EXPECT_EQ(ma2[i(i1 = 2, i2 = 1)], -114);
}
} // end namespace