diff --git a/src/multi_array_operation.h b/src/multi_array_operation.h index bf30109..c038870 100644 --- a/src/multi_array_operation.h +++ b/src/multi_array_operation.h @@ -41,7 +41,7 @@ namespace MultiArrayTools template Block makeBlock(const T* vec, size_t stepSize, size_t blockSize); - + template MBlock makeBlock(T* vec, size_t stepSize, size_t blockSize); @@ -92,7 +92,11 @@ namespace MultiArrayTools OperationMaster(MutableMultiArrayBase& ma, const OpClass& second, std::shared_ptr& index); - + + OperationMaster(MutableMultiArrayBase& ma, const OpClass& second, + std::shared_ptr& index, + std::shared_ptr blockIndex); + MBlock& get(); const Block& get() const; @@ -100,8 +104,8 @@ namespace MultiArrayTools const OperationMaster& block() const; protected: - - //void performAssignment(const OperationBase& in); + + void performAssignment(std::intptr_t blockIndexNum); OpClass const& mSecond; MutableMultiArrayBase& mArrayRef; std::shared_ptr mIndex; @@ -156,7 +160,9 @@ namespace MultiArrayTools const MBlock& get() const; MBlock& get(); - + + OperationRoot& set(std::shared_ptr blockIndex); + std::vector block(const std::shared_ptr blockIndex, bool init = false) const; const OperationRoot& block() const; @@ -165,6 +171,7 @@ namespace MultiArrayTools MutableMultiArrayBase& mArrayRef; std::shared_ptr mIndex; mutable MBlock mBlock; + std::shared_ptr mBlockIndex; // predefine to save time }; template @@ -325,6 +332,33 @@ namespace MultiArrayTools block(blockIndex, true); second.block(blockIndex, true); + + performAssignment(blockIndexNum); + } + + template + OperationMaster:: + OperationMaster(MutableMultiArrayBase& ma, const OpClass& second, + std::shared_ptr& index, + std::shared_ptr blockIndex) : + mSecond(second), mArrayRef(ma), mIndex() + { + MultiRangeFactory mrf( index->range() ); + std::shared_ptr > mr = + std::dynamic_pointer_cast >( mrf.create() ); + mIndex = std::make_shared( mr->begin() ); + (*mIndex) = *index; + + std::intptr_t blockIndexNum = blockIndex->getPtrNum(); + second.block(blockIndex, true); + + performAssignment(blockIndexNum); + } + + + template + void OperationMaster::performAssignment(std::intptr_t blockIndexNum) + { //size_t cnt = 0; //std::clock_t cs = clock(); for(*mIndex = 0; mIndex->pos() != mIndex->max(); mIndex->pp(blockIndexNum) ){ @@ -417,7 +451,8 @@ namespace MultiArrayTools OperationRoot(MutableMultiArrayBase& ma, const std::shared_ptr&... indices) : OperationTemplate >(this), - mArrayRef(ma), mIndex( std::make_shared( mArrayRef.range() ) ) + mArrayRef(ma), mIndex( std::make_shared( mArrayRef.range() ) ), + mBlockIndex(nullptr) { (*mIndex)(indices...); } @@ -426,7 +461,12 @@ namespace MultiArrayTools template OperationMaster OperationRoot::operator=(const OpClass& in) { - return OperationMaster(mArrayRef, in, mIndex); + if(mBlockIndex){ + return OperationMaster(mArrayRef, in, mIndex, mBlockIndex); + } + else { + return OperationMaster(mArrayRef, in, mIndex); + } } template @@ -443,6 +483,14 @@ namespace MultiArrayTools return mBlock; } + template + OperationRoot& + OperationRoot::set(std::shared_ptr blockIndex) + { + mBlockIndex = blockIndex; + return *this; + } + template std::vector OperationRoot::block(const std::shared_ptr blockIndex, bool init) const { diff --git a/src/ranges/container_range.h b/src/ranges/container_range.h index 0f7ad5e..888a752 100644 --- a/src/ranges/container_range.h +++ b/src/ranges/container_range.h @@ -77,6 +77,8 @@ namespace MultiArrayTools std::shared_ptr getVPtr(size_t n); size_t getStepSize(size_t n); + + std::vector infoVec() const; std::string id(); void print(size_t offset); @@ -332,6 +334,14 @@ namespace MultiArrayTools return mBlockSizes[n+1]; } + template + std::vector ContainerIndex::infoVec() const + { + std::vector out(sizeof...(Indices)); + RPackNum::buildInfoVec(out, mIPack); + return std::move( out ); + } + template std::string ContainerIndex::id() { diff --git a/src/ranges/index_base.h b/src/ranges/index_base.h index 5439b1f..171f46a 100644 --- a/src/ranges/index_base.h +++ b/src/ranges/index_base.h @@ -11,6 +11,7 @@ #include "range_base.h" #include "index_type.h" #include "vindex_wrapper.h" +#include "index_info.h" namespace MultiArrayTools { @@ -56,6 +57,8 @@ namespace MultiArrayTools { return I::template S_get(THIS()); } */ std::shared_ptr getVPtr(size_t n) const { return THIS().getVPtr(n); } + + std::vector infoVec() const { return THIS().infoVec(); } size_t getStepSize(size_t n) const { return THIS().getStepSize(n); } diff --git a/src/ranges/index_info.cc b/src/ranges/index_info.cc new file mode 100644 index 0000000..f4c4f28 --- /dev/null +++ b/src/ranges/index_info.cc @@ -0,0 +1,33 @@ + +#include "index_info.h" + +namespace MultiArrayTools +{ + + const IndexInfo* IndexInfo::getPtr(size_t inum) const + { + return &mNext[inum]; + } + + std::intptr_t IndexInfo::getPtrNum() const + { + return mPtrNum; + } + + size_t IndexInfo::dim() const + { + return mDim; + } + + size_t IndexInfo::getStepSize(size_t inum) const + { + return mNext[inum]->getStepSzize(); + } + + size_t IndexInfo::getStepSize() const + { + return mStepSize; + } + + +} // end namespace MultiArrayTools diff --git a/src/ranges/index_info.h b/src/ranges/index_info.h new file mode 100644 index 0000000..427aedd --- /dev/null +++ b/src/ranges/index_info.h @@ -0,0 +1,44 @@ +// -*- C++ -*- + +#ifndef __index_info_h__ +#define __index_info_h__ + +#include +#include + +namespace MultiArrayTools +{ + class IndexInfo; + + class IndexInfo + { + public: + IndexInfo() = delete; + + template + IndexInfo(const IndexClass& ind, size_t stepSize = 1); + + const IndexInfo* getPtr(size_t inum) const; + std::intptr_t getPtrNum() const; + size_t dim() const; + size_t getStepSize(size_t inum) const; + size_t getStepSize() const; + + private: + std::vector mNext; + std::intptr_t mPtrNum; + size_t mDim; + size_t mStepSize; + }; + + template + IndexInfo::IndexInfo(const IndexClass& ind, size_t stepSize) : + mNext(ind.infoVec()), + mPtrNum( reinterpret_cast( &ind ) ), + mDim(ind.rangePtr()->dim()), + mStepSize(stepSize) + {} + +} // end namespace MultiArrayTools + +#endif diff --git a/src/ranges/multi_range.h b/src/ranges/multi_range.h index ed376e3..721b4f8 100644 --- a/src/ranges/multi_range.h +++ b/src/ranges/multi_range.h @@ -94,6 +94,8 @@ namespace MultiArrayTools std::shared_ptr getVPtr(size_t n); size_t getStepSize(size_t n); + std::vector infoVec() const; + std::string id(); void print(size_t offset); }; @@ -376,6 +378,14 @@ namespace MultiArrayTools return mBlockSizes[n+1]; } + template + std::vector MultiIndex::infoVec() const + { + std::vector out(sizeof...(Indices)); + RPackNum::buildInfoVec(out, mIPack); + return std::move( out ); + } + template std::string MultiIndex::id() { diff --git a/src/ranges/rpack_num.h b/src/ranges/rpack_num.h index ffe64c9..fa6ad79 100644 --- a/src/ranges/rpack_num.h +++ b/src/ranges/rpack_num.h @@ -4,6 +4,7 @@ #include #include "vindex_wrapper.h" +#include "index_info.h" namespace MultiArrayHelper { @@ -207,6 +208,15 @@ namespace MultiArrayHelper static_assert( Range::defaultable, "not defaultable" ); RPackNum::template checkDefaultable(); } + + template + static void buildInfoVec(std::vector& out, + std::tuple...>& ip) + { + out.emplace_back(*std::get(ip)); + RPackNum::buildInfoVec(out, ip); + } + }; @@ -361,6 +371,13 @@ namespace MultiArrayHelper static_assert( Range::defaultable, "not defaultable" ); } + template + static void buildInfoVec(std::vector& out, + std::tuple...>& ip) + { + out.emplace_back(*std::get(ip)); + } + }; diff --git a/src/ranges/single_range.h b/src/ranges/single_range.h index 4d51ebc..b1ae553 100644 --- a/src/ranges/single_range.h +++ b/src/ranges/single_range.h @@ -52,6 +52,8 @@ namespace MultiArrayTools std::shared_ptr getVPtr(size_t n); size_t getStepSize(size_t n); + + std::vector infoVec() const; std::string id(); void print(size_t offset); @@ -206,6 +208,12 @@ namespace MultiArrayTools return std::shared_ptr(); } + template + std::vector SingleIndex::infoVec() const + { + return std::move( std::vector() ); + } + template size_t SingleIndex::getStepSize(size_t n) { diff --git a/src/tests/op_perf_test.cc b/src/tests/op_perf_test.cc index 3e753ce..a7cc1fb 100644 --- a/src/tests/op_perf_test.cc +++ b/src/tests/op_perf_test.cc @@ -9,6 +9,8 @@ #include #include +#define ONLY_SPIN + namespace MAT = MultiArrayTools; namespace { @@ -62,7 +64,7 @@ namespace { { return std::make_tuple(static_cast( ts )...); } - +#ifndef ONLY_SPIN class OpTest_Performance { public: @@ -119,7 +121,7 @@ namespace { std::vector cv1; std::vector cv2; }; - +#endif class OpTest_Spin { public: @@ -160,11 +162,12 @@ namespace { auto gamma = MAT::getIndex(); auto delta = MAT::getIndex(); auto deltap = MAT::getIndex(); + auto vdeltap = MAT::make_viwb(deltap); auto mix = MAT::mkMIndex( alpha, beta, gamma ); std::clock_t begin = std::clock(); - res1(delta, deltap) = ma(delta, alpha, alpha, beta, beta, gamma, gamma, deltap).c(mix); + res1(delta, deltap).set(vdeltap) = ma(delta, alpha, alpha, beta, beta, gamma, gamma, deltap).c(mix); std::clock_t end = std::clock(); std::cout << "MultiArray time: " << static_cast( end - begin ) / CLOCKS_PER_SEC << std::endl; @@ -212,7 +215,7 @@ namespace { << std::endl; std::cout << "ratio: " << static_cast( end - begin ) / static_cast( end2 - begin2 ) << std::endl; } - +#ifndef ONLY_SPIN void OpTest_Performance::PCheck() { MultiArray ma2(mrptr, cv2); @@ -251,16 +254,18 @@ namespace { //assert( xround( res.at(mkt(700,900)) ) == xround(res2[700*vs1 + 900]) ); } - +#endif } // anonymous namspace int main(int argc, char** argv) { +#ifndef ONLY_SPIN OpTest_Performance pt; - OpTest_Spin st; - - st.contract(); pt.PCheck(); +#endif + OpTest_Spin st; + st.contract(); + return 0; } diff --git a/src/tests/op_unit_test.cc b/src/tests/op_unit_test.cc index 83b6e04..dc1394c 100644 --- a/src/tests/op_unit_test.cc +++ b/src/tests/op_unit_test.cc @@ -215,11 +215,12 @@ namespace { auto gamma = MAT::getIndex(); auto delta = MAT::getIndex(); auto deltap = MAT::getIndex(); + auto vdeltap = MAT::make_viwb(deltap); auto mix = MAT::mkMIndex( alpha, beta, gamma ); std::clock_t begin = std::clock(); - res1(delta, deltap) = ma(delta, alpha, alpha, beta, beta, gamma, gamma, deltap).c(mix); + res1(delta, deltap).set(vdeltap) = ma(delta, alpha, alpha, beta, beta, gamma, gamma, deltap).c(mix); std::clock_t end = std::clock(); std::cout << "MultiArray time: " << static_cast( end - begin ) / CLOCKS_PER_SEC << std::endl;