diff --git a/include/alp/reference/matrix.hpp b/include/alp/reference/matrix.hpp index e8b99a902..0419c1f8f 100644 --- a/include/alp/reference/matrix.hpp +++ b/include/alp/reference/matrix.hpp @@ -39,7 +39,7 @@ #include //#include //for help with dealing with pattern matrix input #include -#include +#include #include #include #include @@ -359,28 +359,28 @@ namespace alp { namespace internal { /** Forward declaration */ - template< typename T, typename ImfR, typename ImfC, typename Smf, bool is_original > + template< typename T, typename ImfR, typename ImfC, typename MappingPolynomial, bool is_original > class MatrixContainer; /** Container reference getters used by friend functions of specialized Matrix */ - template< typename T, typename ImfR, typename ImfC, typename Smf, bool is_original > - const Vector< T, reference > & getContainer( const MatrixContainer< T, ImfR, ImfC, Smf, is_original > & A ); + template< typename T, typename ImfR, typename ImfC, typename MappingPolynomial, bool is_original > + const Vector< T, reference > & getContainer( const MatrixContainer< T, ImfR, ImfC, MappingPolynomial, is_original > & A ); - template< typename T, typename ImfR, typename ImfC, typename Smf, bool is_original > - Vector< T, reference > & getContainer( MatrixContainer< T, ImfR, ImfC, Smf, is_original > & A ); + template< typename T, typename ImfR, typename ImfC, typename MappingPolynomial, bool is_original > + Vector< T, reference > & getContainer( MatrixContainer< T, ImfR, ImfC, MappingPolynomial, is_original > & A ); /** Container reference getters. Defer the call to base class friend function */ template< typename T, typename Structure, enum Density density, typename View, typename ImfR, typename ImfC > const Vector< T, reference > & getContainer( const alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference > & A ) { return getContainer( static_cast< const MatrixContainer< T, ImfR, ImfC, - typename alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference >::smf_type, + typename alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference >::mapping_polynomial_type, alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference >::is_original > & >( A ) ); } template< typename T, typename Structure, enum Density density, typename View, typename ImfR, typename ImfC > Vector< T, reference > & getContainer( alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference > & A ) { return getContainer( static_cast< MatrixContainer< T, ImfR, ImfC, - typename alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference >::smf_type, + typename alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference >::mapping_polynomial_type, alp::Matrix< T, Structure, density, View, ImfR, ImfC, reference >::is_original > & >( A ) ); } @@ -480,8 +480,8 @@ namespace alp { * @tparam is_original True if the class is an original container * False if the class is a view of another matrix */ - template< typename T, typename ImfR, typename ImfC, typename Smf, bool is_original > - class MatrixContainer : public MatrixBase< MatrixContainer< T, ImfR, ImfC, Smf, is_original > > { + template< typename T, typename ImfR, typename ImfC, typename MappingPolynomial, bool is_original > + class MatrixContainer : public MatrixBase< MatrixContainer< T, ImfR, ImfC, MappingPolynomial, is_original > > { public: /** Expose static properties */ @@ -495,7 +495,7 @@ namespace alp { typedef ImfC imfc_type; protected: - typedef MatrixContainer< T, ImfR, ImfC, Smf, is_original > self_type; + typedef MatrixContainer< T, ImfR, ImfC, MappingPolynomial, is_original > self_type; friend MatrixBase< self_type >; typedef typename std::conditional< @@ -536,7 +536,7 @@ namespace alp { * \see AMF */ public: // TODO: Temporarily expose AMF publicly until proper getters are implemented - smf::AMF< ImfR, ImfC, Smf > amf; + storage::AMF< ImfR, ImfC, MappingPolynomial > amf; protected: /** @@ -546,7 +546,7 @@ namespace alp { * @return A pair of dimensions. */ std::pair< size_t, size_t > _dims() const noexcept { - return std::make_pair( amf.imf_r.n, amf.imf_c.n ); + return amf.getLogicalDimensions(); } friend const Vector< T, reference > & getContainer( const self_type &A ) { @@ -605,13 +605,13 @@ namespace alp { * TODO: Add the storage scheme a parameter to the constructor * so that allocation can be made accordingly, generalizing the full case. */ - MatrixContainer( smf::AMF< ImfR, ImfC, Smf > amf ) : + MatrixContainer( storage::AMF< ImfR, ImfC, MappingPolynomial > amf ) : // enable only if ImfR and ImfC are imf::Id - container( internal::Vector< T, reference >( 10 /* TODO */ ) ), + container( internal::Vector< T, reference >( amf.getStorageDimensions() ) ), amf( amf ) {} /** View on another container */ - MatrixContainer( Vector< T, reference > &container, smf::AMF< ImfR, ImfC, Smf > amf ) : + MatrixContainer( Vector< T, reference > &container, storage::AMF< ImfR, ImfC, MappingPolynomial > amf ) : container( container ), amf( amf ) {} @@ -733,7 +733,7 @@ namespace alp { */ template< typename T, typename View, typename ImfR, typename ImfC > class Matrix< T, structures::General, Density::Dense, View, ImfR, ImfC, reference > : - public internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, std::is_same< typename View::applied_to, void >::value > { + public internal::MatrixContainer< T, ImfR, ImfC, storage::polynomials::Full_type, std::is_same< typename View::applied_to, void >::value > { private: typedef Matrix< T, structures::General, Density::Dense, View, ImfR, ImfC, reference > self_type; @@ -756,7 +756,7 @@ namespace alp { public: /** Exposes the types and the static properties. */ typedef structures::General structure; - typedef smf::Full_t smf_type; + typedef storage::polynomials::Full_type mapping_polynomial_type; /** * Indicates if a matrix is an original container. * False if it is a view over another matrix or a functor. @@ -767,7 +767,7 @@ namespace alp { * Expose the base type class to enable internal functions to cast * the type of objects of this class to the base class type. */ - typedef internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original > base_type; + typedef internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original > base_type; // A general Structure knows how to define a reference to itself (which is an original reference view) // as well as other static views. @@ -801,11 +801,12 @@ namespace alp { typename = typename std::enable_if< std::is_same< TargetMatrixType, void >::value >::type > Matrix( const size_t rows, const size_t cols, const size_t cap = 0 ) : - internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original >( - smf::AMF< ImfR, ImfC, smf::Full_t >( + internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original >( + storage::AMF< ImfR, ImfC, mapping_polynomial_type >( imf::Id( rows ), imf::Id( cols ), - smf::SMFFactory::Full( cols ) + storage::polynomials::Create< mapping_polynomial_type >( cols ), + rows * cols ) ) { (void)cap; @@ -825,19 +826,9 @@ namespace alp { std::is_same< TargetMatrixType, target_type >::value >::type > Matrix( TargetMatrixType &target_matrix, ImfR imf_r, ImfC imf_c ) : - internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original >( + internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original >( getContainer( target_matrix ), - smf::AMF< ImfR, ImfC, smf::Full_t >( - imf::ComposedFactory::create< typename target_type::imfr_type, ImfR >( - target_matrix.amf.imf_r, - imf_r - ), - imf::ComposedFactory::create< typename target_type::imfc_type, ImfC >( - target_matrix.amf.imf_c, - imf_c - ), - target_matrix.amf.smf - ) + storage::AMFFactory::Create( target_matrix.amf, imf_r, imf_c ) ) {} /** @@ -864,7 +855,7 @@ namespace alp { // Square Matrix template< typename T, typename View, typename ImfR, typename ImfC > class Matrix< T, structures::Square, Density::Dense, View, ImfR, ImfC, reference > : - public internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, std::is_same< typename View::applied_to, void >::value > { + public internal::MatrixContainer< T, ImfR, ImfC, storage::polynomials::Full_type, std::is_same< typename View::applied_to, void >::value > { private: typedef Matrix< T, structures::Square, Density::Dense, View, ImfR, ImfC, reference > self_type; @@ -883,9 +874,9 @@ namespace alp { public: /** Exposes the types and the static properties. */ typedef structures::Square structure; - typedef smf::Full_t smf_type; + typedef storage::polynomials::Full_type mapping_polynomial_type; static constexpr bool is_original = std::is_same< target_type, void >::value; - typedef internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original > base_type; + typedef internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original > base_type; template < view::Views view_tag, bool d=false > struct view_type; @@ -906,11 +897,12 @@ namespace alp { typename = typename std::enable_if< std::is_same< TargetMatrixType, void >::value >::type > Matrix( const size_t dim, const size_t cap = 0 ) : - internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original >( - smf::AMF< ImfR, ImfC, smf::Full_t >( + internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original >( + storage::AMF< ImfR, ImfC, mapping_polynomial_type >( imf::Id( dim ), imf::Id( dim ), - smf::SMFFactory::Full( dim ) + storage::polynomials::Create< mapping_polynomial_type >( dim ), + dim * dim ) ) { (void)cap; @@ -924,19 +916,9 @@ namespace alp { std::is_same< TargetMatrixType, target_type >::value >::type > Matrix( TargetMatrixType &target_matrix, ImfR imf_r, ImfC imf_c ) : - internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original >( + internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original >( getContainer( target_matrix ), - smf::AMF< ImfR, ImfC, smf::Full_t >( - imf::ComposedFactory::create< typename target_type::imfr_type, ImfR >( - target_matrix.amf.imf_r, - imf_r - ), - imf::ComposedFactory::create< typename target_type::imfc_type, ImfC >( - target_matrix.amf.imf_c, - imf_c - ), - target_matrix.amf.smf - ) + storage::AMFFactory::Create( target_matrix.amf, imf_r, imf_c ) ) {} /** @@ -960,7 +942,7 @@ namespace alp { // UpperTriangular Matrix template< typename T, typename View, typename ImfR, typename ImfC > class Matrix< T, structures::UpperTriangular, Density::Dense, View, ImfR, ImfC, reference > : - public internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, std::is_same< typename View::applied_to, void >::value > { + public internal::MatrixContainer< T, ImfR, ImfC, storage::polynomials::Full_type, std::is_same< typename View::applied_to, void >::value > { private: typedef Matrix< T, structures::UpperTriangular, Density::Dense, View, ImfR, ImfC, reference > self_type; @@ -983,9 +965,9 @@ namespace alp { public: /** Exposes the element type and the structure. */ typedef structures::UpperTriangular structure; - typedef smf::Full_t smf_type; + typedef storage::polynomials::Full_type mapping_polynomial_type; static constexpr bool is_original = std::is_same< target_type, void >::value; - typedef internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original > base_type; + typedef internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original > base_type; template < view::Views view_tag, bool d=false > struct view_type; @@ -1006,11 +988,12 @@ namespace alp { typename = typename std::enable_if< std::is_same< TargetMatrixType, void >::value >::type > Matrix( const size_t rows, const size_t cols, const size_t cap = 0 ) : - internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original >( - smf::AMF< ImfR, ImfC, smf::Full_t >( + internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original >( + storage::AMF< ImfR, ImfC, mapping_polynomial_type >( imf::Id( rows ), imf::Id( cols ), - smf::SMFFactory::Full( cols ) + storage::polynomials::Create< mapping_polynomial_type >( cols ), + rows * cols ) ) { (void)cap; @@ -1024,19 +1007,9 @@ namespace alp { std::is_same< TargetMatrixType, target_type >::value >::type > Matrix( TargetMatrixType &target_matrix, ImfR imf_r, ImfC imf_c ) : - internal::MatrixContainer< T, ImfR, ImfC, smf::Full_t, is_original >( + internal::MatrixContainer< T, ImfR, ImfC, mapping_polynomial_type, is_original >( getContainer( target_matrix ), - smf::AMF< ImfR, ImfC, smf::Full_t >( - imf::ComposedFactory::create< typename target_type::imfr_type, ImfR >( - target_matrix.amf.imf_r, - imf_r - ), - imf::ComposedFactory::create< typename target_type::imfc_type, ImfC >( - target_matrix.amf.imf_c, - imf_c - ), - target_matrix.amf.smf - ) + storage::AMFFactory::Create( target_matrix.amf, imf_r, imf_c ) ) {} /** Constructor for a view over another matrix using default IMFs (Identity). @@ -1267,7 +1240,15 @@ namespace alp { template< typename TargetStructure, typename TargetImfR, typename TargetImfC, typename T, typename Structure, enum Density density, typename View, typename ImfR, typename ImfC, enum Backend backend > - alp::Matrix< T, TargetStructure, density, view::Original< alp::Matrix< T, Structure, density, View, ImfR, ImfC, backend > >, TargetImfR, TargetImfC, backend > + alp::Matrix< + T, + TargetStructure, + density, + view::Original< alp::Matrix< T, Structure, density, View, ImfR, ImfC, backend > >, + typename imf::composed_type< TargetImfR, ImfR >::type, + typename imf::composed_type< TargetImfC, ImfC >::type, + backend + > get_view( alp::Matrix< T, Structure, density, View, ImfR, ImfC, backend > &source, TargetImfR imf_r, TargetImfC imf_c ) { @@ -1371,7 +1352,15 @@ namespace alp { template< typename T, typename Structure, enum Density density, typename View, typename ImfR, typename ImfC, enum Backend backend > - Matrix< T, Structure, density, view::Original< Matrix< T, Structure, density, View, ImfR, ImfC, backend > >, imf::Strided, imf::Strided, backend > + Matrix< + T, + Structure, + density, + view::Original< Matrix< T, Structure, density, View, ImfR, ImfC, backend > >, + typename imf::composed_type< imf::Strided, ImfR >::type, + typename imf::composed_type< imf::Strided, ImfC >::type, + backend + > get_view( Matrix< T, Structure, density, View, ImfR, ImfC, backend > &source, const utils::range &rng_r, const utils::range &rng_c ) { @@ -1484,7 +1473,15 @@ namespace alp { typename TargetStructure, typename IndexType, typename IndexStructure, typename IndexView, typename T, typename Structure, enum Density density, typename View, typename ImfR, typename ImfC, enum Backend backend > - alp::Matrix< T, TargetStructure, density, view::Original< alp::Matrix< T, Structure, density, View, ImfR, ImfC, backend > >, imf::Select, imf::Select, backend > + alp::Matrix< + T, + TargetStructure, + density, + view::Original< alp::Matrix< T, Structure, density, View, ImfR, ImfC, backend > >, + typename imf::composed_type< imf::Select, ImfR >::type, + typename imf::composed_type< imf::Select, ImfC >::type, + backend + > get_view( alp::Matrix< T, Structure, density, View, ImfR, ImfC, backend > &source, const Vector< IndexType, IndexStructure, density, IndexView, backend > & sel_r, const Vector< IndexType, IndexStructure, density, IndexView, backend > & sel_c ) { diff --git a/include/alp/smf.hpp b/include/alp/smf.hpp deleted file mode 100644 index ec640d620..000000000 --- a/include/alp/smf.hpp +++ /dev/null @@ -1,244 +0,0 @@ -/* - * Copyright 2021 Huawei Technologies Co., Ltd. - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -/** - * - * @file - * - * This file registers available storage mapping functions (SMFs). - * SMFs are maps between logical and physical storage space. - * - */ - -#ifndef _H_ALP_SMF -#define _H_ALP_SMF - -#include -#include "structures.hpp" - -namespace alp { - - namespace smf { - - namespace polynomials { - - /** - * Implements the polynomial - * ( A*a*x^2 + B*b*y^2 + C*c*x*y + D*d*x + E*e*y + F*f ) / Denominator - * where uppercase coefficients are compile-time constant, - * lowercase coefficients are run-time constant, - * and x and y are variables. - * All coefficients, variables are integers and all operations are integer - * operations. - * - * The purpose of compile-time constant coefficients is to allow compile-time - * optimizations for zero factors. - */ - template< size_t Ax2, size_t Ay2, size_t Axy, size_t Ax, size_t Ay, size_t A0, size_t Denominator > - struct BivariateQuadratic { - - static_assert( Denominator != 0, "Denominator cannot be zero (division by zero)."); - - const size_t ax2, ay2, axy, ax, ay, a0; - - BivariateQuadratic( - const size_t ax2, const size_t ay2, const size_t axy, - const size_t ax, const size_t ay, - const size_t a0 ) : - ax2( ax2 ), ay2( ay2 ), axy( axy ), - ax( ax ), ay( ay ), - a0( a0 ) {} - - size_t f( const size_t x, const size_t y ) const { - return (Ax2 * ax2 * x * x + - Ay2 * ay2 * y * y + - Axy * axy * x * y + - Ax * ax * x + - Ay * ay * y + - A0 * a0) / Denominator; - } - - }; // BivariateQuadratic - - }; - - typedef polynomials::BivariateQuadratic< 0, 0, 0, 0, 0, 0, 1 > None_t; - typedef polynomials::BivariateQuadratic< 0, 0, 0, 1, 1, 0, 1 > Full_t; - typedef polynomials::BivariateQuadratic< 0, 0, 0, 0, 0, 0, 1 > Packed_t; // TODO - typedef polynomials::BivariateQuadratic< 0, 0, 0, 0, 0, 0, 1 > Band_t; // TODO - - /** - * Creates a polynomial object representing desired storage mapping function. - */ - struct SMFFactory { - - /** - * For row-major storage, N represents number of columns. - * For column-major storage, N represents number of rows. - */ - static Full_t Full( size_t N ) { - return Full_t( 0, 0, 0, N, 1, 0 ); - } - - /** - * For triangular and symmetric structures (assuming square matrices) - */ - static Packed_t Packed( size_t N ) { - (void)N; - return Packed_t( /* TODO */ 0, 0, 0, 0, 0, 0 ); - } - - /** - * For Banded structures - */ - static Band_t Band( size_t rows, size_t cols, size_t kl, size_t ku ) { - (void)rows; - (void)cols; - (void)kl; - (void)ku; - return Band_t( /* TODO */ 0, 0, 0, 0, 0, 0 ); - } - - }; // SMFFactory - - /** - * Provides a type of composed Access Mapping Function - * expressed as a BivariateQuadratic polynomial depending - * on the types of the IMFs and the SMF. - */ - template< typename ImfR, typename ImfC, typename Smf > - struct Composition { - typedef polynomials::BivariateQuadratic< 1, 1, 1, 1, 1, 1, 1 > type; - }; - - template<> - struct Composition< imf::Strided, imf::Strided, Full_t > { - typedef polynomials::BivariateQuadratic< 0, 0, 0, 1, 1, 1, 1 > type; - }; - - /** - * Access Mapping Function (AMF) maps a logical matrix coordinates (i, j) - * to a corresponding matrix element's location in the physical container. - * AMF take into account the index mapping function associated to each - * coordinate (rows and columns) and the storage mapping function that - * maps logical coordinates to physical ones. - * - * For certain combinations of IMFs and SMFs it is possible to fuse the - * index computation in a single function call. - */ - template< typename ImfR, typename ImfC, typename Smf > - class AMF { - - public: - - const ImfR imf_r; - const ImfC imf_c; - const Smf smf; - - public: - - AMF( ImfR &&imf_r, ImfC &&imf_c, Smf smf ) : imf_r( imf_r ), imf_c( imf_c ), smf( smf ) {} - - /** - * Returns a storage index based on the coordinates in the - * logical iteration space. - * - * @param[in] i row-coordinate - * @param[in] j column-coordinate - * @param[in] s current process ID - * @param[in] P total number of processes - * - * @return storage index corresponding to the provided logical - * coordinates and parameters s and P. - */ - size_t getStorageIndex( const size_t i, const size_t j, const size_t s, const size_t P ) const { -#ifdef _DEBUG - std::cout << "Calling AMF::getStorageIndex()\n"; -#endif - (void)s; - (void)P; - return smf.f( imf_r.map( i ), imf_c.map( j ) ); - } - - /** - * Returns coordinates in the logical iteration space based on - * the storage index. - * - * @param[in] storageIndex storage index in the physical - * iteration space - * @param[in] s current process ID - * @param[in] P total number of processes - * - * @return a pair of row- and column-coordinates in the - * logical iteration space. - */ - std::pair< size_t, size_t > getCoords( const size_t storageIndex, const size_t s, const size_t P ) const; - }; - - template< typename Smf > - class AMF< imf::Strided, imf::Strided, Smf > { - - public: - - /** For size checks */ - const imf::Strided imf_r; - const imf::Strided imf_c; - const Smf smf; - typedef typename Composition< imf::Strided, imf::Strided, Smf >::type Composition_t; - const Composition_t amf; - - Composition_t fusion( - const imf::Strided &imf_r, - const imf::Strided &imf_c, - const Smf &smf - ) const { - return Composition_t( - smf.ax2 * imf_r.s * imf_r.s, // ax2 ( for x^2 ) - smf.ay2 * imf_c.s * imf_c.s, // ay2 ( for y*2 ) - smf.axy * imf_r.s * imf_c.s, // axy ( for x * y ) - imf_r.s * ( 2 * smf.ax2 * imf_r.b + smf.axy * imf_c.b + smf.ax ), // ax ( for x ) - imf_c.s * ( 2 * smf.ay2 * imf_c.b + smf.axy * imf_r.b + smf.ay ), // ay ( for y ) - smf.ax2 * imf_r.b * imf_r.b + smf.ay2 * imf_c.b * imf_c.b + - smf.axy * imf_r.b * imf_c.b + smf.ax * imf_r.b + smf.ay * imf_c.b + smf.a0 // a0 - ); - } - - public: - - AMF( const imf::Strided &imf_r, const imf::Strided &imf_c, const Smf &smf ) : - imf_r( imf_r ), imf_c( imf_c ), smf( smf ), amf( fusion( imf_r, imf_c, smf ) ) { - } - - size_t getStorageIndex( const size_t i, const size_t j, const size_t s, const size_t P ) const { -#ifdef _DEBUG - std::cout << "Calling AMF::getStorageIndex()\n"; -#endif - // TODO: Maybe these asserts should be pushed to debug mode - // for performance reasons. - (void)s; - (void)P; - assert( i < imf_r.n ); - assert( j < imf_c.n ); - return amf.f( i, j ); - } - - }; // class AMF< imf::Strided, imf::Strided, Smf > - - }; // namespace smf - -} // namespace alp - -#endif // _H_ALP_SMF diff --git a/include/alp/storage.hpp b/include/alp/storage.hpp new file mode 100644 index 000000000..4b954eac6 --- /dev/null +++ b/include/alp/storage.hpp @@ -0,0 +1,336 @@ +/* + * Copyright 2021 Huawei Technologies Co., Ltd. + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + */ + +/** + * + * @file + * + * This file registers mechanisms for coordinate mapping between + * logical and physical iteration spaces. + * + */ + +#ifndef _H_ALP_SMF +#define _H_ALP_SMF + +#include +#include "structures.hpp" + +namespace alp { + + namespace storage { + + /** + * The namespace containts polynomials used to map coordinates + * between logical and physical iteration spaces, + * associated type traits and helper classes. + */ + namespace polynomials { + + /** + * Implements the polynomial + * ( A*a*x^2 + B*b*y^2 + C*c*x*y + D*d*x + E*e*y + F*f ) / Denominator + * where uppercase coefficients are compile-time constant, + * lowercase coefficients are run-time constant, + * and x and y are variables. + * All coefficients and variables are integers and all operations are integer + * operations. + * + * The purpose of compile-time constant coefficients is to allow compile-time + * optimizations for zero terms/monomials. + * + * Denominator allows for implementaiton of polynomials with integer division, + * e.g., n * ( n + 1 ) / 2, + * while avoiding the need for floating point coefficients and operations. + * + * @tparam Ax2 Static coefficient corresponding to x^2 + * @tparam Ay2 Static coefficient corresponding to y^2 + * @tparam Axy Static coefficient corresponding to x*y + * @tparam Ax Static coefficient corresponding to x + * @tparam Ay Static coefficient corresponding to y + * @tparam A0 Static coefficient corresponding to constant term + * @tparam Denominator Static denominator dividing the whole polynomial + */ + template< size_t Ax2, size_t Ay2, size_t Axy, size_t Ax, size_t Ay, size_t A0, size_t Denominator > + struct BivariateQuadratic { + + static_assert( Denominator != 0, "Denominator cannot be zero (division by zero)."); + + const size_t ax2, ay2, axy, ax, ay, a0; + + BivariateQuadratic( + const size_t ax2, const size_t ay2, const size_t axy, + const size_t ax, const size_t ay, + const size_t a0 ) : + ax2( ax2 ), ay2( ay2 ), axy( axy ), + ax( ax ), ay( ay ), + a0( a0 ) {} + + size_t evaluate( const size_t x, const size_t y ) const { + return (Ax2 * ax2 * x * x + + Ay2 * ay2 * y * y + + Axy * axy * x * y + + Ax * ax * x + + Ay * ay * y + + A0 * a0) / Denominator; + } + + }; // BivariateQuadratic + + typedef BivariateQuadratic< 0, 0, 0, 0, 0, 0, 1 > None_type; + typedef BivariateQuadratic< 0, 0, 0, 1, 1, 0, 1 > Full_type; + typedef BivariateQuadratic< 0, 0, 0, 0, 0, 0, 1 > Packed_type; // TODO + typedef BivariateQuadratic< 0, 0, 0, 0, 0, 0, 1 > Band_type; // TODO + + /** + * Polynomial factory method + */ + template< typename PolynomialType, typename... Args > + PolynomialType Create( Args... args ) { + return PolynomialType( args... ); + } + + /** Specialization for Full storage of type i * dim + j */ + template<> + Full_type Create< Full_type >( size_t dim ) { + return Full_type( 0, 0, 0, dim, 1, 0 ); + } + + }; // namespace polynomials + + /** + * Provides a type of composed Access Mapping Function + * expressed as a BivariateQuadratic polynomial depending + * on the types of the IMFs and the Mapping Polynomial. + */ + template< typename ImfR, typename ImfC, typename MappingPolynomial > + struct Composition { + typedef polynomials::BivariateQuadratic< 1, 1, 1, 1, 1, 1, 1 > type; + }; + + template<> + struct Composition< imf::Strided, imf::Strided, polynomials::Full_type > { + typedef polynomials::BivariateQuadratic< 0, 0, 0, 1, 1, 1, 1 > type; + }; + + /** Forward declaration */ + class AMFFactory; + + /** + * Access Mapping Function (AMF) maps logical matrix coordinates (i, j) + * to the corresponding matrix element's location in the physical container. + * + * To calculate the mapping, the AMF first applies logical-to-logical + * mapping provided by one IMF per coordinate (row and column). + * A bivariate polynomial (called mapping polynomial) takes these two + * output coordinates as inputs to calculate the position is physical + * storage of the requested element (logical-to-physical mapping). + * + * For certain combinations of IMFs and mapping polynomial types it is + * possible to fuse the index computation into a single function call. + */ + template< typename ImfR, typename ImfC, typename MappingPolynomial > + class AMF { + + friend AMFFactory; + + private: + + const ImfR imf_r; + const ImfC imf_c; + const MappingPolynomial map_poly; + const size_t storage_dimensions; + + public: + + AMF( ImfR &&imf_r, ImfC &&imf_c, MappingPolynomial map_poly, const size_t storage_dimensions ) : + imf_r( imf_r ), imf_c( imf_c ), map_poly( map_poly ), storage_dimensions( storage_dimensions ) {} + + /** + * Returns dimensions of the logical layout of the associated container. + * + * @return A pair of two values, number of rows and columns, respectively. + */ + std::pair< size_t, size_t> getLogicalDimensions() const { + return std::make_pair( imf_r.n, imf_c.n ); + } + + /** + * Returns dimensions of the physical layout of the associated container. + * + * @return The size of the physical container. + */ + std::size_t getStorageDimensions() const { + return storage_dimensions; + } + + /** + * Returns a storage index based on the coordinates in the + * logical iteration space. + * + * @param[in] i row-coordinate + * @param[in] j column-coordinate + * @param[in] s current process ID + * @param[in] P total number of processes + * + * @return storage index corresponding to the provided logical + * coordinates and parameters s and P. + */ + size_t getStorageIndex( const size_t i, const size_t j, const size_t s, const size_t P ) const { +#ifdef _DEBUG + std::cout << "Calling AMF::getStorageIndex()\n"; +#endif + (void)s; + (void)P; + return map_poly.evaluate( imf_r.map( i ), imf_c.map( j ) ); + } + + /** + * Returns coordinates in the logical iteration space based on + * the storage index. + * + * @param[in] storageIndex storage index in the physical + * iteration space + * @param[in] s current process ID + * @param[in] P total number of processes + * + * @return a pair of row- and column-coordinates in the + * logical iteration space. + */ + std::pair< size_t, size_t > getCoords( const size_t storageIndex, const size_t s, const size_t P ) const; + + }; // class AMF + + /** + * Specialization for strided IMFs + * which allows for fusion of IMFs with the mapping polynomial + */ + template< typename MappingPolynomial > + class AMF< imf::Strided, imf::Strided, MappingPolynomial > { + + friend AMFFactory; + + private: + + /** For size checks */ + const imf::Strided imf_r; + const imf::Strided imf_c; + const MappingPolynomial map_poly; + typedef typename Composition< imf::Strided, imf::Strided, MappingPolynomial >::type Composition_type; + const Composition_type amf; + + const size_t storage_dimensions; + + /** Returns a polynomial representing fusion of the IMFs and the original polynomial */ + Composition_type fusion( + const imf::Strided &imf_r, + const imf::Strided &imf_c, + const MappingPolynomial &map_poly + ) const { + return Composition_type( + map_poly.ax2 * imf_r.s * imf_r.s, // ax2 ( for x^2 ) + map_poly.ay2 * imf_c.s * imf_c.s, // ay2 ( for y*2 ) + map_poly.axy * imf_r.s * imf_c.s, // axy ( for x * y ) + imf_r.s * ( 2 * map_poly.ax2 * imf_r.b + map_poly.axy * imf_c.b + map_poly.ax ), // ax ( for x ) + imf_c.s * ( 2 * map_poly.ay2 * imf_c.b + map_poly.axy * imf_r.b + map_poly.ay ), // ay ( for y ) + map_poly.ax2 * imf_r.b * imf_r.b + map_poly.ay2 * imf_c.b * imf_c.b + + map_poly.axy * imf_r.b * imf_c.b + map_poly.ax * imf_r.b + map_poly.ay * imf_c.b + map_poly.a0 // a0 + ); + } + + public: + + AMF( const imf::Strided &imf_r, const imf::Strided &imf_c, const MappingPolynomial &map_poly, const size_t storage_dimensions ) : + imf_r( imf_r ), imf_c( imf_c ), map_poly( map_poly ), amf( fusion( imf_r, imf_c, map_poly ) ), storage_dimensions( storage_dimensions ) { + } + + std::pair< size_t, size_t> getLogicalDimensions() const { + return std::make_pair( imf_r.n, imf_c.n ); + } + + std::size_t getStorageDimensions() const { + return storage_dimensions; + } + + size_t getStorageIndex( const size_t i, const size_t j, const size_t s, const size_t P ) const { +#ifdef _DEBUG + std::cout << "Calling AMF::getStorageIndex()\n"; +#endif + // TODO: Maybe these asserts should be pushed to debug mode + // for performance reasons. + (void)s; + (void)P; + assert( i < imf_r.n ); + assert( j < imf_c.n ); + return amf.evaluate( i, j ); + } + + std::pair< size_t, size_t > getCoords( const size_t storageIndex, const size_t s, const size_t P ) const; + + }; // class AMF< imf::Strided, imf::Strided, storage > + + /** + * Factory class that creates AMF objects of the appropriate type + * depending on the types of provided IMFs and polynomials. + */ + class AMFFactory { + + public: + + template< typename MappingPolynomial > + static AMF< imf::Id, imf::Id, MappingPolynomial > + Create( imf::Id imf_r, imf::Id imf_c, MappingPolynomial map_poly, const size_t storage_dimensions ) { + return AMF< imf::Id, imf::Id, MappingPolynomial >( imf_r, imf_c, map_poly, storage_dimensions ); + } + + template< + typename OriginalImfR, typename OriginalImfC, typename MappingPolynomial, + typename ViewImfR, typename ViewImfC + > + static AMF< + typename imf::composed_type< ViewImfR, OriginalImfR >::type, + typename imf::composed_type< ViewImfC, OriginalImfC >::type, + MappingPolynomial + > Create( + const AMF< OriginalImfR, OriginalImfC, MappingPolynomial > &original_amf, + ViewImfR view_imf_r, + ViewImfC view_imf_c + ) { + return AMF< + typename imf::composed_type< ViewImfR, OriginalImfR >::type, + typename imf::composed_type< ViewImfC, OriginalImfC >::type, + MappingPolynomial + >( + imf::ComposedFactory::create< ViewImfR, OriginalImfR >( + view_imf_r, + original_amf.imf_r + ), + imf::ComposedFactory::create< ViewImfC, OriginalImfC >( + view_imf_c, + original_amf.imf_c + ), + original_amf.map_poly, + original_amf.storage_dimensions + ); + } + + }; // class AMFFactory + + }; // namespace storage + +} // namespace alp + +#endif // _H_ALP_SMF