From 8f94c77f4db253eb35be42e775586dd6d4dda779 Mon Sep 17 00:00:00 2001
From: Christian Zimmermann <christian.zimmermann@ur.de>
Date: Wed, 21 Mar 2018 19:18:57 +0100
Subject: [PATCH] functional multi array test; there are some remaining issues
 (maybe conceptional)...

---
 src/include/arith.h                           |  18 ++-
 src/include/functional_multi_array.h          |  64 ++++++++--
 src/include/multi_array_operation.h           |   8 +-
 src/include/ranges/anonymous_range.h          |   7 ++
 src/include/ranges/range_base.h               |   1 +
 .../ranges/range_types/classic_range.h        |   1 +
 src/include/ranges/range_types/spin_range.h   |   1 +
 src/include/ranges/single_range.h             |  15 ++-
 src/lib/ranges/anonymous_range.cc             |  13 --
 src/tests/op_unit_test.cc                     | 119 +++++++++++++++++-
 10 files changed, 206 insertions(+), 41 deletions(-)

diff --git a/src/include/arith.h b/src/include/arith.h
index b77d9aa..54ed0a1 100644
--- a/src/include/arith.h
+++ b/src/include/arith.h
@@ -8,6 +8,8 @@ namespace MultiArrayHelper
     template <typename T>
     struct plus
     {
+	static constexpr bool FISSTATIC = true;
+	
 	static inline T apply(T a1, T a2)
 	{
 	    return a1 + a2;
@@ -17,6 +19,8 @@ namespace MultiArrayHelper
     template <typename T>
     struct minus
     {
+	static constexpr bool FISSTATIC = true;
+	
 	static inline T apply(T a1, T a2)
 	{
 	    return a1 - a2;
@@ -26,6 +30,8 @@ namespace MultiArrayHelper
     template <typename T>
     struct multiplies
     {
+	static constexpr bool FISSTATIC = true;
+	
 	static inline T apply(T a1, T a2)
 	{
 	    return a1 * a2;
@@ -35,22 +41,26 @@ namespace MultiArrayHelper
     template <typename T>
     struct divides
     {
+	static constexpr bool FISSTATIC = true;
+	
 	static inline T apply(T a1, T a2)
 	{
 	    return a1 / a2;
 	}
     };
-
+    /*
     template <typename T, class Func>
     struct dynamic_function
     {
+	static constexpr bool FISSTATIC = false;
+	
 	template <typename... Us>
-	static inline T apply(const Func& f, Us... args)
+	inline T apply(Us... args)
 	{
 	    return f(args...);
 	}
-    }
-    
+    };
+    */
 } // end namespace MultiArrayHelper
     
 #endif
diff --git a/src/include/functional_multi_array.h b/src/include/functional_multi_array.h
index 6207432..0e15808 100644
--- a/src/include/functional_multi_array.h
+++ b/src/include/functional_multi_array.h
@@ -7,6 +7,11 @@
 
 namespace MultiArrayTools
 {
+
+    template <class Index>
+    auto indexToSlice(const std::shared_ptr<Index>& i)
+	-> Slice<typename Index::MetaType,typename Index::RangeType>;
+
     
     template <typename T, class Function, class... SRanges>
     class FunctionalMultiArray : public MultiArrayBase<T,SRanges...>
@@ -15,9 +20,9 @@ namespace MultiArrayTools
 
 	typedef ContainerRange<T,SRanges...> CRange;
 	typedef MultiArrayBase<T,CRange> MAB;
-	typedef typename MultiArrayBase<T,CRange>::const_iterator const_iterator;
-	typedef typename CRange::IndexType IndexType;
-
+	//typedef typename MultiArrayBase<T,CRange>::const_iterator const_iterator;
+	typedef ContainerIndex<T,typename SRanges::IndexType...> IndexType;
+	//typedef typename CRange::IndexType IndexType;
 
     private:
 	mutable T mVal;
@@ -27,20 +32,21 @@ namespace MultiArrayTools
 	
 	DEFAULT_MEMBERS(FunctionalMultiArray);
 	FunctionalMultiArray(const std::shared_ptr<SRanges>&... ranges, const Function& func);
-
+	FunctionalMultiArray(const std::shared_ptr<SRanges>&... ranges);
+	
 	virtual const T& operator[](const IndexType& i) const override;
-
+	
 	virtual bool isConst() const override;	
 	virtual bool isSlice() const override;
 
 	// EVALUTAION CLASS ??!!!!
 
-	auto operator()(std::shared_ptr<typename SRanges::IndexType>&... inds) const
-	    -> decltype( mkOperation( mFunc, ConstOperationRoot<T,SRanges...>(inds.metaPtr(), inds), ... ) )
+	virtual auto operator()(std::shared_ptr<typename SRanges::IndexType>&... inds) const
+	    -> decltype( mkOperation( mFunc, ConstOperationRoot<T,SRanges>( indexToSlice( inds ), inds) ... ) ) override;
 	
     };
 
-    
+  
 } // namespace MultiArrayTools
 
 /* ========================= *
@@ -49,6 +55,14 @@ namespace MultiArrayTools
 
 namespace MultiArrayTools
 {
+    
+    template <class Index>
+    auto indexToSlice(const std::shared_ptr<Index>& i)
+	-> Slice<typename Index::MetaType, typename Index::RangeType>
+    {
+	return Slice<typename Index::MetaType, typename Index::RangeType>( i->range(), i->metaPtr() );
+    }
+
 
     /****************************
      *   FunctionalMultiArray   *	     
@@ -59,15 +73,39 @@ namespace MultiArrayTools
     FunctionalMultiArray<T,Range,Function>::FunctionalMultiArray(const Range& range) :
 	MultiArrayBase<T,SRanges...>(range), mFunc() {}
     */
+    template <bool FISSTATIC>
+    struct Application
+    {
+	template <typename T, class Function, class Index>
+	static inline T apply(const Function& f, const Index& i)
+	{
+	    return f(i.meta());
+	}
+    };
+
+    template <>
+    struct Application<true>
+    {
+	template <typename T, class Function, class Index>
+	static inline T apply(const Function& f, const Index& i)
+	{
+	    return Function::apply(i.meta());
+	}
+    };
+	
     template <typename T, class Function, class... SRanges>
     FunctionalMultiArray<T,Function,SRanges...>::FunctionalMultiArray(const std::shared_ptr<SRanges>&... ranges,
-								  const Function& func) :
+								      const Function& func) :
 	MultiArrayBase<T,SRanges...>(ranges...), mFunc(func) {}
 
     template <typename T, class Function, class... SRanges>
-    const T& FunctionalMultiArray<T,Function,SRanges...>::operator[](const typename CRange::IndexType& i) const
+    FunctionalMultiArray<T,Function,SRanges...>::FunctionalMultiArray(const std::shared_ptr<SRanges>&... ranges) :
+	MultiArrayBase<T,SRanges...>(ranges...) {}
+
+    template <typename T, class Function, class... SRanges>
+    const T& FunctionalMultiArray<T,Function,SRanges...>::operator[](const IndexType& i) const
     {
-	mVal = mFunc(i.meta());
+	mVal = Application<Function::FISSTATIC>::template apply<T,Function,IndexType>(mFunc, i.meta());
 	return mVal;
     }
 
@@ -86,9 +124,9 @@ namespace MultiArrayTools
     template <typename T, class Function, class... SRanges>
     auto FunctionalMultiArray<T,Function,SRanges...>::
     operator()(std::shared_ptr<typename SRanges::IndexType>&... inds) const
-	-> decltype( mkOperation( mFunc, ConstOperationRoot<T,SRanges...>(inds.metaPtr(), inds), ... ) )
+	-> decltype( mkOperation( mFunc, ConstOperationRoot<T,SRanges>( indexToSlice( inds ), inds) ... ) )
     {
- 	return mkOperation( mFunc, ConstOperationRoot<T,SRanges...>(inds.metaPtr(), inds), ... );
+ 	return mkOperation( mFunc, ConstOperationRoot<T,SRanges>( indexToSlice( inds ), inds ) ... );
     }
 
     
diff --git a/src/include/multi_array_operation.h b/src/include/multi_array_operation.h
index 4cba9e6..6a9f3bc 100644
--- a/src/include/multi_array_operation.h
+++ b/src/include/multi_array_operation.h
@@ -238,10 +238,12 @@ namespace MultiArrayTools
 	typedef OpFunction F;
 	
 	static constexpr size_t SIZE = RootSum<Ops...>::SIZE;
+	static constexpr bool FISSTATIC = OpFunction::FISSTATIC;
 
     private:
 	std::tuple<Ops...> mOps;
-
+	std::shared_ptr<OpFunction> mF; // only if non-static
+	
     public:
 	typedef decltype(PackNum<sizeof...(Ops)-1>::template mkSteps<Ops...>(0, mOps)) ETuple;
 		
@@ -261,9 +263,9 @@ namespace MultiArrayTools
 
     template <class OpFunction, class... Ops>
     auto mkOperation(const OpFunction& f, const Ops&... ops)
-	-> Operation<OpFunction::value_type,OpFunction,Ops...>
+	-> Operation<typename OpFunction::value_type,OpFunction,Ops...>
     {
-	return Operation<OpFunction::value_type,OpFunction,Ops...>(ops...);
+	return Operation<typename OpFunction::value_type,OpFunction,Ops...>(ops...);
     }
     
     template <typename T, class Op, class IndexType>
diff --git a/src/include/ranges/anonymous_range.h b/src/include/ranges/anonymous_range.h
index 8f0a8bc..c22f6dd 100644
--- a/src/include/ranges/anonymous_range.h
+++ b/src/include/ranges/anonymous_range.h
@@ -33,6 +33,13 @@ namespace MultiArrayTools
     
     class AnonymousRange : public RangeInterface<AnonymousIndex>
     {
+    public:
+
+	static constexpr bool defaultable = false;
+        static constexpr size_t ISSTATIC = 0;
+        static constexpr size_t SIZE = -1;
+        static constexpr bool HASMETACONT = false;
+	
 	typedef RangeBase RB;
 	typedef typename RangeInterface<AnonymousIndex>::IndexType IndexType;
 	
diff --git a/src/include/ranges/range_base.h b/src/include/ranges/range_base.h
index d8e3c32..3632d9c 100644
--- a/src/include/ranges/range_base.h
+++ b/src/include/ranges/range_base.h
@@ -78,6 +78,7 @@ namespace MultiArrayTools
     {
     public:
 
+	//typedef typename Index::MetaType MetaType;
 	typedef Index IndexType;
 
 	virtual Index begin() const = 0;
diff --git a/src/include/ranges/range_types/classic_range.h b/src/include/ranges/range_types/classic_range.h
index e236443..d3cf092 100644
--- a/src/include/ranges/range_types/classic_range.h
+++ b/src/include/ranges/range_types/classic_range.h
@@ -29,6 +29,7 @@ namespace MultiArrayTools
 	typedef RangeBase RB;
 	typedef typename RangeInterface<SingleIndex<size_t,SpaceType::NONE> >::IndexType IndexType;
 	typedef SingleRange<size_t,SpaceType::NONE> RangeType;
+	typedef size_t MetaType;
 	
 	virtual size_t size() const override;
 	virtual size_t dim() const override;
diff --git a/src/include/ranges/range_types/spin_range.h b/src/include/ranges/range_types/spin_range.h
index 0a8d118..34f5447 100644
--- a/src/include/ranges/range_types/spin_range.h
+++ b/src/include/ranges/range_types/spin_range.h
@@ -33,6 +33,7 @@ namespace MultiArrayTools
 	typedef RangeBase RB;
 	typedef typename RangeInterface<SingleIndex<size_t,SpaceType::SPIN> >::IndexType IndexType;
 	typedef SingleRange<size_t,SpaceType::SPIN> RangeType;
+	typedef size_t MetaType;
 	
 	virtual size_t size() const override;
 	virtual size_t dim() const override;
diff --git a/src/include/ranges/single_range.h b/src/include/ranges/single_range.h
index 490dd55..8ab53ce 100644
--- a/src/include/ranges/single_range.h
+++ b/src/include/ranges/single_range.h
@@ -74,7 +74,7 @@ namespace MultiArrayTools
 	    -> For<SingleIndex<U,TYPE>,Expr,ForType::HIDDEN>;
 
     private:
-	U* mMetaPtr;
+	const U* mMetaPtr;
     };
 
     template <typename U, SpaceType TYPE>
@@ -97,6 +97,7 @@ namespace MultiArrayTools
 	typedef RangeBase RB;
 	typedef SingleIndex<U,TYPE> IndexType;
 	typedef SingleRange RangeType;
+	typedef U MetaType;
 	//typedef typename RangeInterface<SingleIndex<U,TYPE> >::IndexType IndexType;
 	
 	virtual size_t size() const override;
@@ -140,18 +141,19 @@ namespace MultiArrayTools
     template <bool HASMETACONT>
     struct MetaPtrSetter
     {
-	template <typename U, class Range>
-	U* set(Range* r)
+	template <class Range>
+	static const typename Range::MetaType* set(Range* r)
 	{
 	    return &r->get(0);
 	}
+
     };
 
     template <>
     struct MetaPtrSetter<false>
     {
-	template <typename U, class Range>
-	U* set(Range* r)
+	template <class Range>
+	static const typename Range::MetaType* set(Range* r)
 	{
 	    return nullptr;
 	}
@@ -160,7 +162,8 @@ namespace MultiArrayTools
     template <typename U, SpaceType TYPE>
     SingleIndex<U,TYPE>::SingleIndex(const std::shared_ptr<SingleRange<U,TYPE> >& range) :
 	IndexInterface<SingleIndex<U,TYPE>,U>(range, 0),
-	mMetaPtr(MetaPtrSetter<SingleIndex<U,TYPE>::HASMETACONT>::set(IB::mRangePtr)) {}
+	mMetaPtr(MetaPtrSetter<SingleIndex<U,TYPE>::RangeType::HASMETACONT>::set
+		 ( dynamic_cast<RangeType*>(IB::mRangePtr.get() ) ) ) {}
 
     template <typename U, SpaceType TYPE>
     IndexType SingleIndex<U,TYPE>::type() const
diff --git a/src/lib/ranges/anonymous_range.cc b/src/lib/ranges/anonymous_range.cc
index 0d3bdd5..a302dc1 100644
--- a/src/lib/ranges/anonymous_range.cc
+++ b/src/lib/ranges/anonymous_range.cc
@@ -54,19 +54,6 @@ namespace MultiArrayTools
 	i = size();
 	return i;
     }
-
-    // put this in the interface class !!!
-    /*
-    std::shared_ptr<VIWB> AnonymousRange::index() const
-    {
-	//typedef IndexWrapper<IndexType> IW;
-	return std::shared_ptr<VIWB>();
-		    std::make_shared<IW>
-	    (std::make_shared<IndexType>
-	     ( std::dynamic_pointer_cast<AnonymousRange>
-	     ( std::shared_ptr<RangeBase>( RB::mThis ) ) ) );  //!!!!
-    }
-    */
     
     /*****************
      *   Functions   *
diff --git a/src/tests/op_unit_test.cc b/src/tests/op_unit_test.cc
index 7c9a99c..2b6f66f 100644
--- a/src/tests/op_unit_test.cc
+++ b/src/tests/op_unit_test.cc
@@ -62,6 +62,31 @@ namespace {
     {
 	return std::make_tuple(static_cast<size_t>( ts )...);
     }
+
+    template <typename T>
+    struct Pow
+    {
+	static constexpr bool FISSTATIC = true;
+	typedef T value_type;
+	
+	static inline T apply(T b, T e)
+	{
+	    return pow(b, e);
+	}
+
+	template <class... Ops>
+	static auto mk(const Ops&... ops)
+	    -> Operation<T,Pow<T>,Ops...>
+	{
+	    return Operation<T,Pow<T>,Ops...>(ops...);
+	}
+
+	static inline T apply(const std::tuple<double, double>& arg)
+	{
+	    return pow(std::get<0>(arg), std::get<1>(arg));
+	}
+    };
+
     
     class OpTest_1Dim : public ::testing::Test
     {
@@ -208,6 +233,50 @@ namespace {
 	std::shared_ptr<SR> sr;
     };
 
+    class MetaOp_Test : public ::testing::Test
+    {
+    protected:
+
+	typedef SingleRangeFactory<double,SpaceType::ANY> SRF;
+	typedef SRF::oType SR;
+
+	const size_t s1 = 10;
+	const size_t s2 = 17;
+	
+	MetaOp_Test()
+	{
+	    mv1 = { 2.476, 9.665, 1.289, 2.89, 77.04, -11.09, 100.4, 2.0, -26.5, -0.001 };
+	    mv2 = { 44.56, 23.097, -117.3, -0.0765, 3.445, 0.02389, -4.0112, 10.567, 8.99, -177.2, 475.3,
+		    11.111, 13.108, -35.6, 64.32, 2.44, -12.};
+
+	    assert(mv1.size() == s1); // just to prevent typos...
+	    assert(mv1.size() == s2);
+	    
+	    swapFactory<SRF>(rfbptr, mv1);
+	    sr1ptr = std::dynamic_pointer_cast<SR>(rfbptr->create());
+
+	    swapFactory<SRF>(rfbptr, mv2);
+	    sr2ptr = std::dynamic_pointer_cast<SR>(rfbptr->create());
+	}
+
+	std::shared_ptr<RangeFactoryBase> rfbptr;
+	std::shared_ptr<SR> sr1ptr;
+	std::shared_ptr<SR> sr2ptr;
+
+	std::vector<double> mv1;
+	std::vector<double> mv2;
+    };
+
+    
+    TEST_F(MetaOp_Test, SimpleCall)
+    {
+	FunctionalMultiArray<double,Pow<double>,SR,SR> fma(sr1ptr, sr2ptr);
+
+	auto i = fma.begin();
+	
+	EXPECT_EQ( xround( i.at( mkt(9.665, -0.0765) ) ), xround( pow(9.665, -0.0765) ) );
+    }
+    
     TEST_F(OpTest_Spin, Contract)
     {
 	MultiArray<double,SR,SR,SR,SR,SR,SR,SR,SR> ma(sr, sr, sr, sr, sr, sr, sr, sr, data);
@@ -392,6 +461,52 @@ namespace {
 	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'B')) ), xround(2.911 + 4.790 - 2.210)  );
     }
 
+    template <typename T>
+    struct Monopole
+    {
+	static constexpr bool FISSTATIC = true;
+
+	static inline T apply(T f0, T x, T n)
+	{
+	    return f0 / ( 1 + x / n );
+	}
+
+	template <class... Ops>
+	static auto mk(const Ops&... ops)
+	    -> Operation<T,Monopole<T>,Ops...>
+	{
+	    return Operation<T,Monopole<T>,Ops...>(ops...);
+	}
+    };
+    
+    TEST_F(OpTest_MDim, ExecFOp)
+    {
+	MultiArray<double,MRange,SRange> res(mr1ptr,sr4ptr);
+	MultiArray<double,MRange> ma1(mr1ptr, v3);
+	MultiArray<double,SRange> ma2(sr4ptr, v2);
+	MultiArray<double,SRange> ma3(sr4ptr, v4);
+	
+	auto i1 = MAT::getIndex( mr1ptr );
+	auto i2 = MAT::getIndex( sr4ptr );
+	
+	res(i1,i2) = Monopole<double>::mk( ma1(i1) , ma2(i2) , ma3(i2) );
+
+	EXPECT_EQ( xround( res.at(mkt(mkt('1','a'),'A')) ), xround(0.353 / ( 1 + 8.870 / 1.470) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('1','a'),'B')) ), xround(0.353 / ( 1 + 4.790 / 2.210) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('1','b'),'A')) ), xround(4.005 / ( 1 + 8.870 / 1.470) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('1','b'),'B')) ), xround(4.005 / ( 1 + 4.790 / 2.210) ) );
+
+	EXPECT_EQ( xround( res.at(mkt(mkt('2','a'),'A')) ), xround(1.070 / ( 1 + 8.870 / 1.470) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('2','a'),'B')) ), xround(1.070 / ( 1 + 4.790 / 2.210) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('2','b'),'A')) ), xround(2.310 / ( 1 + 8.870 / 1.470) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('2','b'),'B')) ), xround(2.310 / ( 1 + 4.790 / 2.210) ) );
+
+	EXPECT_EQ( xround( res.at(mkt(mkt('3','a'),'A')) ), xround(9.243 / ( 1 + 8.870 / 1.470) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('3','a'),'B')) ), xround(9.243 / ( 1 + 4.790 / 2.210) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'A')) ), xround(2.911 / ( 1 + 8.870 / 1.470) ) );
+	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'B')) ), xround(2.911 / ( 1 + 4.790 / 2.210) ) );
+    }
+
     TEST_F(OpTest_MDim, ExecOp3)
     {
 	MultiArray<double,MRange,SRange> res(mr1ptr,sr4ptr);
@@ -422,7 +537,7 @@ namespace {
 	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'A')) ), xround(2.911 + 0.373 + 1.470) );
 	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'B')) ), xround(2.911 + 0.373 + 2.210) );
     }
-
+    /*
     TEST_F(OpTest_MDim, ExecAnonOp1)
     {
 	MultiArray<double,MRange,SRange> res(mr1ptr,sr4ptr);
@@ -453,7 +568,7 @@ namespace {
 	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'A')) ), xround(2.911 + 0.373 + 1.470) );
 	EXPECT_EQ( xround( res.at(mkt(mkt('3','b'),'B')) ), xround(2.911 + 0.373 + 2.210) );
     }
-
+    */
     
 } // anonymous namspace