From d5a1422f341988eeee7b313d09d9c7da21cf067d Mon Sep 17 00:00:00 2001
From: Christian Zimmermann <christian.zimmermann@ur.de>
Date: Sun, 23 Dec 2018 13:01:13 +0100
Subject: [PATCH] sub range contraction (tested)

---
 src/include/ranges/subrange.h | 17 +++++++++++---
 src/tests/op_unit_test.cc     | 43 +++++++++--------------------------
 2 files changed, 25 insertions(+), 35 deletions(-)

diff --git a/src/include/ranges/subrange.h b/src/include/ranges/subrange.h
index b9b1885..97477ab 100644
--- a/src/include/ranges/subrange.h
+++ b/src/include/ranges/subrange.h
@@ -78,11 +78,10 @@ namespace MultiArrayTools
 	auto ifor(size_t step, Expr ex) const
 	    -> For<SubIndex<Index>,SubExpr<Index,Expr>>;
 
-        /*
 	template <class Expr>
 	auto iforh(size_t step, Expr ex) const
-	    -> For<SubIndex<Index>,Expr,ForType::HIDDEN>;
-        */
+	    -> For<SubIndex<Index>,SubExpr<Index,Expr>,ForType::HIDDEN>;
+
     private:
 	std::shared_ptr<RangeType> mExplicitRangePtr;
 	//const U* mMetaPtr;
@@ -322,6 +321,18 @@ namespace MultiArrayTools
 	       &mExplicitRangePtr->subset(), ex ) );
     }
 
+    template <class Index>
+    template <class Expr>
+    auto SubIndex<Index>::iforh(size_t step, Expr ex) const
+	-> For<SubIndex<Index>,SubExpr<Index,Expr>,ForType::HIDDEN>
+    {
+        return For<SubIndex<Index>,SubExpr<Index,Expr>,ForType::HIDDEN>
+            (this, step, SubExpr<Index,Expr>
+	     ( mFullIndex, reinterpret_cast<std::intptr_t>(this),
+	       &mExplicitRangePtr->subset(), ex ) );
+    }
+
+    
     /************************
      *   SubRangeFactory    *
      ************************/
diff --git a/src/tests/op_unit_test.cc b/src/tests/op_unit_test.cc
index ec59a34..410831e 100644
--- a/src/tests/op_unit_test.cc
+++ b/src/tests/op_unit_test.cc
@@ -546,7 +546,7 @@ namespace {
         SubRangeFactory<SRange> subf(sr2ptr, std::vector<size_t>({0,2}));
         auto subptr = MAT::createExplicit(subf);
         
-        MultiArray<double,SubRange<SRange>,SRange,SRange> res(subptr,sr3ptr,sr1ptr,0.);
+        MultiArray<double,SRange,SRange> res(sr3ptr,sr1ptr,0.);
 	MultiArray<double,SRange,SubRange<SRange>,SRange> res2(sr3ptr,subptr,sr1ptr,0.);
 
         auto i1 = MAT::getIndex( sr1ptr );
@@ -555,42 +555,21 @@ namespace {
         auto si = MAT::getIndex( subptr );
         (*si)(i2);
 
-        res(si,i3,i1) = ma2(i3,i2) - ma1(i1,i2,i3);
+        res(i3,i1) = (ma2(i3,i2) - ma1(i1,i2,i3)).c(si);
 	res2(i3,si,i1) = ma2(i3,i2) - ma1(i1,i2,i3);
 
-        EXPECT_EQ( res.size(), 8 );
-        EXPECT_EQ( res.vdata().size(), 8 );
-        EXPECT_EQ( MAT::rptr<0>( res )->size(), 2 );
+        EXPECT_EQ( res2.size(), 8 );
+        EXPECT_EQ( res2.vdata().size(), 8 );
+        EXPECT_EQ( MAT::rptr<1>( res2 )->size(), 2 );
 
-        EXPECT_EQ( MAT::rptr<0>( res )->isMeta('1'), true );
-        EXPECT_EQ( MAT::rptr<0>( res )->isMeta('3'), true );
+        EXPECT_EQ( MAT::rptr<1>( res2 )->isMeta('1'), true );
+        EXPECT_EQ( MAT::rptr<1>( res2 )->isMeta('3'), true );
         
-        EXPECT_EQ( xround( res.at(mkt('1','a','x')) ), xround(ma2.at(mkt('a','1')) - ma1.at(mkt('x','1','a'))) );
-        EXPECT_EQ( xround( res.at(mkt('1','a','l')) ), xround(ma2.at(mkt('a','1')) - ma1.at(mkt('l','1','a'))) );
-
-        EXPECT_EQ( xround( res.at(mkt('1','b','x')) ), xround(ma2.at(mkt('b','1')) - ma1.at(mkt('x','1','b'))) );
-        EXPECT_EQ( xround( res.at(mkt('1','b','l')) ), xround(ma2.at(mkt('b','1')) - ma1.at(mkt('l','1','b'))) );
-
-        EXPECT_EQ( xround( res.at(mkt('3','a','x')) ), xround(ma2.at(mkt('a','3')) - ma1.at(mkt('x','3','a'))) );
-        EXPECT_EQ( xround( res.at(mkt('3','a','l')) ), xround(ma2.at(mkt('a','3')) - ma1.at(mkt('l','3','a'))) );
-
-        EXPECT_EQ( xround( res.at(mkt('3','b','x')) ), xround(ma2.at(mkt('b','3')) - ma1.at(mkt('x','3','b'))) );
-        EXPECT_EQ( xround( res.at(mkt('3','b','l')) ), xround(ma2.at(mkt('b','3')) - ma1.at(mkt('l','3','b'))) );
-
-        auto resx = res.format( subptr->outRange(), sr3ptr, sr1ptr );
-
-        EXPECT_EQ( xround( resx.at(mkt('1','a','x')) ), xround(ma2.at(mkt('a','1')) - ma1.at(mkt('x','1','a'))) );
-        EXPECT_EQ( xround( resx.at(mkt('1','a','l')) ), xround(ma2.at(mkt('a','1')) - ma1.at(mkt('l','1','a'))) );
-
-        EXPECT_EQ( xround( resx.at(mkt('1','b','x')) ), xround(ma2.at(mkt('b','1')) - ma1.at(mkt('x','1','b'))) );
-        EXPECT_EQ( xround( resx.at(mkt('1','b','l')) ), xround(ma2.at(mkt('b','1')) - ma1.at(mkt('l','1','b'))) );
-
-        EXPECT_EQ( xround( resx.at(mkt('3','a','x')) ), xround(ma2.at(mkt('a','3')) - ma1.at(mkt('x','3','a'))) );
-        EXPECT_EQ( xround( resx.at(mkt('3','a','l')) ), xround(ma2.at(mkt('a','3')) - ma1.at(mkt('l','3','a'))) );
-
-        EXPECT_EQ( xround( resx.at(mkt('3','b','x')) ), xround(ma2.at(mkt('b','3')) - ma1.at(mkt('x','3','b'))) );
-        EXPECT_EQ( xround( resx.at(mkt('3','b','l')) ), xround(ma2.at(mkt('b','3')) - ma1.at(mkt('l','3','b'))) );
+        EXPECT_EQ( xround( res.at(mkt('a','x')) ), xround((ma2.at(mkt('a','1')) - ma1.at(mkt('x','1','a'))) + (ma2.at(mkt('a','3')) - ma1.at(mkt('x','3','a')) ) ));
+	EXPECT_EQ( xround( res.at(mkt('a','l')) ), xround((ma2.at(mkt('a','1')) - ma1.at(mkt('l','1','a'))) + (ma2.at(mkt('a','3')) - ma1.at(mkt('l','3','a')) ) ));
 
+        EXPECT_EQ( xround( res.at(mkt('b','x')) ), xround((ma2.at(mkt('b','1')) - ma1.at(mkt('x','1','b'))) + (ma2.at(mkt('b','3')) - ma1.at(mkt('x','3','b')) ) ) );
+        EXPECT_EQ( xround( res.at(mkt('b','l')) ), xround((ma2.at(mkt('b','1')) - ma1.at(mkt('l','1','b'))) + (ma2.at(mkt('b','3')) - ma1.at(mkt('l','3','b')) ) ) );
 
 	EXPECT_EQ( xround( res2.at(mkt('a','1','x')) ), xround(ma2.at(mkt('a','1')) - ma1.at(mkt('x','1','a'))) );
         EXPECT_EQ( xround( res2.at(mkt('a','1','l')) ), xround(ma2.at(mkt('a','1')) - ma1.at(mkt('l','1','a'))) );