diff --git a/Grid/serialisation/BaseIO.h b/Grid/serialisation/BaseIO.h index 6adb87d3..771cc628 100644 --- a/Grid/serialisation/BaseIO.h +++ b/Grid/serialisation/BaseIO.h @@ -36,13 +36,19 @@ Author: Guido Cossu #include namespace Grid { + // TODO Support Grid::complex from GPU port + template using Grid_complex = std::complex; + + // Returns original type, except for Grid_complex, where it returns the underlying type + template struct RealType { using type = T; }; + template struct RealType> { using type = T; }; + namespace EigenIO { template struct is_complex : public std::false_type {}; - template struct is_complex> + template struct is_complex> : std::integral_constant::value> {}; // Eigen tensors can be composed of arithmetic scalar and complex types - // TODO Support Grid::comples from GPU port template struct is_scalar : std::integral_constant::value || is_complex::value> {}; @@ -88,11 +94,12 @@ namespace Grid { template struct Traits {}; // C needed for specialisation // This defines the bottom level - i.e. it's a description of the underlying scalar template struct Traits::value, void>::type> { + using scalar_type = T; // Type of the underlying scalar + using scalar_real = typename RealType::type; // real type underlying scalar_type static constexpr unsigned int depth = 0; // How many levels of Grid Tensor there are (TensorLevel) static constexpr unsigned int rank = 0; // The rank of the grid tensor (i.e. how many indices used) static constexpr unsigned int rank_non_trivial = 0; // As per rank, but excludes those of dimension 1 static constexpr unsigned int count = 1; // total number of elements (i.e. product of dimensions) - using scalar_type = T; // Type of the underlying scalar static constexpr std::size_t scalar_size = sizeof(T); // Size of the underlying scalar in bytes static constexpr std::size_t size = scalar_size * count; // total size of elements in bytes static constexpr std::size_t Dimension(unsigned int dim) { return 0; } // Dimension size @@ -114,11 +121,12 @@ namespace Grid { // count = 48 }; template struct Traits> { + using scalar_type = typename Traits::scalar_type; + using scalar_real = typename RealType::type; static constexpr unsigned int depth = 1 + Traits::depth; static constexpr unsigned int rank = 0 + Traits::rank; static constexpr unsigned int rank_non_trivial = 0 + Traits::rank_non_trivial; static constexpr unsigned int count = 1 * Traits::count; - using scalar_type = typename Traits::scalar_type; static constexpr std::size_t scalar_size = Traits::scalar_size; static constexpr std::size_t size = scalar_size * count; static constexpr std::size_t Dimension(unsigned int dim) { @@ -127,11 +135,12 @@ namespace Grid { return Traits::DimensionNT(dim); } }; template struct Traits> { + using scalar_type = typename Traits::scalar_type; + using scalar_real = typename RealType::type; static constexpr unsigned int depth = 1 + Traits::depth; static constexpr unsigned int rank = 1 + Traits::rank; static constexpr unsigned int rank_non_trivial = (N>1 ? 1 : 0) + Traits::rank_non_trivial; static constexpr unsigned int count = N * Traits::count; - using scalar_type = typename Traits::scalar_type; static constexpr std::size_t scalar_size = Traits::scalar_size; static constexpr std::size_t size = scalar_size * count; static constexpr std::size_t Dimension(unsigned int dim) { @@ -141,11 +150,12 @@ namespace Grid { } }; template struct Traits> { + using scalar_type = typename Traits::scalar_type; + using scalar_real = typename RealType::type; static constexpr unsigned int depth = 1 + Traits::depth; static constexpr unsigned int rank = 2 + Traits::rank; static constexpr unsigned int rank_non_trivial = (N>1 ? 2 : 0) + Traits::rank_non_trivial; static constexpr unsigned int count = N * N * Traits::count; - using scalar_type = typename Traits::scalar_type; static constexpr std::size_t scalar_size = Traits::scalar_size; static constexpr std::size_t size = scalar_size * count; static constexpr std::size_t Dimension(unsigned int dim) { @@ -230,16 +240,22 @@ namespace Grid { } } - // Used for sequential initialisations - template constexpr T Flag = 1; - template constexpr std::complex Flag> {1, -1}; - // Returns the type of the real part of an arithmetic type - template struct RealType { using type = T; }; - template struct RealType> { using type = T; }; + // Sequential initialisation of tensors + // Would have preferred to define template variables for this, but that's c++ 17 + template + typename std::enable_if::value && !EigenIO::is_complex::scalar_type>::value, void>::type + SequentialInit( ETensor &ET, typename EigenIO::Traits::scalar_type Inc = 1 ) + { + using Traits = EigenIO::Traits; + using scalar_type = typename Traits::scalar_type; + for_all( ET, [&](scalar_type &c, typename ETensor::Index n, const std::array &Dims ) { + c = Inc * static_cast(n); + } ); + } template - typename std::enable_if::value, void>::type - Sequential_Init( ETensor &ET, typename EigenIO::Traits::scalar_type Inc = Flag::scalar_type> ) + typename std::enable_if::value && EigenIO::is_complex::scalar_type>::value, void>::type + SequentialInit( ETensor &ET, typename EigenIO::Traits::scalar_type Inc={1,-1}) { using Traits = EigenIO::Traits; using scalar_type = typename Traits::scalar_type; @@ -247,7 +263,7 @@ namespace Grid { c = Inc * static_cast::type>(n); } ); } - + // Helper to dump a tensor #ifdef DEBUG #define dump_tensor(args...) dump_tensor_func(args) diff --git a/tests/IO/Test_serialisation.cc b/tests/IO/Test_serialisation.cc index ae61eebe..342fc214 100644 --- a/tests/IO/Test_serialisation.cc +++ b/tests/IO/Test_serialisation.cc @@ -82,7 +82,7 @@ bool b = false; template void ioTest(const std::string &filename, const O &object, const std::string &name) { - std::cout << name << " IO test: writing ..."; + std::cout << "IO test: " << name << " -> " << filename << " ..."; // writer needs to be destroyed so that writing physically happens { W writer(filename); @@ -92,149 +92,105 @@ void ioTest(const std::string &filename, const O &object, const std::string &nam std::cout << " done. reading..."; R reader(filename); - O buf; + std::unique_ptr buf( new O ); // In case object too big for stack - read(reader, "testobject", buf); - bool good = Serializable::CompareMember(object, buf); + read(reader, "testobject", *buf); + bool good = Serializable::CompareMember(object, *buf); if (!good) { std::cout << " failure!" << std::endl; if constexpr (EigenIO::is_tensor::value) - dump_tensor(buf,"???"); + dump_tensor(*buf,"???"); exit(EXIT_FAILURE); } std::cout << " done." << std::endl; } #ifdef HAVE_HDF5 -typedef Eigen::Tensor ShortRank5Tensor; -//typedef int TestScalar; typedef std::complex TestScalar; -typedef Eigen::Tensor TestTensor; -typedef Eigen::TensorFixedSize, Eigen::StorageOptions::RowMajor> TestTensorFixed; -typedef std::vector aTestTensorFixed; -typedef Eigen::TensorFixedSize, Eigen::StorageOptions::RowMajor> LSCTensor; -typedef Eigen::TensorFixedSize, Eigen::StorageOptions::RowMajor> LCMTensor; -// From Test_serialisation.cc +using TensorSingle = Eigen::TensorFixedSize>; +using TensorSimple = Eigen::Tensor, 6>; +typedef Eigen::Tensor TensorRank5UShort; +typedef Eigen::Tensor TensorRank5IntAlt; +typedef Eigen::Tensor TensorRank3; +typedef Eigen::TensorFixedSize, Eigen::StorageOptions::RowMajor> Tensor_9_4_2; +typedef std::vector aTensor_9_4_2; +typedef Eigen::TensorFixedSize> LSCTensor; +#ifdef DEBUG +typedef Eigen::TensorFixedSize,2>,7>,3>, Eigen::Sizes<2,2,11,10,9>, Eigen::StorageOptions::RowMajor> LCMTensor; +#endif + class PerambIOTestClass: Serializable { public: using PerambTensor = Eigen::Tensor; GRID_SERIALIZABLE_CLASS_MEMBERS(PerambIOTestClass - , ShortRank5Tensor, shortRank5Tensor - , PerambTensor, Perambulator + , SpinColourVector, spinColourVector + , SpinColourMatrix, spinColourMatrix , std::vector, DistilParameterNames , std::vector, DistilParameterValues + , PerambTensor, Perambulator , PerambTensor, Perambulator2 - , SpinColourVector, scv - , SpinColourMatrix, scm - //, TestTensor, Critter - //, TestTensorFixed, FixedCritter - //, aTestTensorFixed, aFixedCritter - //, LSCTensor, MyLSCTensor - //, LCMTensor, MyLCMTensor + , TensorRank5UShort, tensorRank5UShort + , TensorRank3, tensorRank3 + , Tensor_9_4_2, tensor_9_4_2 + , aTensor_9_4_2, atensor_9_4_2 + , LSCTensor, MyLSCTensor +#ifdef DEBUG + , LCMTensor, MyLCMTensor +#endif ); - PerambIOTestClass() : Perambulator(2,3,1,4,5,1), Perambulator2(7,1,6,1,5,1), - DistilParameterNames {"alpha", "beta", "gamma", "delta", "epsilon", "what's f?"}, - DistilParameterValues{2,3,1,4,5,1} - , shortRank5Tensor{5,4,3,2,1} - //, Critter(7,3,2)//, aFixedCritter(3) + PerambIOTestClass() + : DistilParameterNames {"alpha", "beta", "gamma", "delta", "epsilon", "zeta"} + , DistilParameterValues{2,3,1,4,5,1} + , Perambulator(2,3,1,4,5,1) + , Perambulator2(7,1,6,1,5,1) + , tensorRank5UShort{5,4,3,2,1} + , tensorRank3(7,3,2) { - Sequential_Init(Perambulator); - Sequential_Init(Perambulator2, {-3.1415927,7}); - Sequential_Init(shortRank5Tensor); + Grid_complex Flag{1,-3.1415927}; + SequentialInit(Perambulator, Flag); + SequentialInit(Perambulator2, Flag); + SequentialInit(tensorRank5UShort); + SequentialInit(tensorRank3, Flag); + SequentialInit(tensor_9_4_2, Flag); + for( auto &t : atensor_9_4_2 ) + SequentialInit(t, Flag); + SequentialInit( MyLSCTensor ); +#ifdef DEBUG + SequentialInit( MyLCMTensor ); +#endif } }; -void EigenHdf5IOTest(void) { - SpinColourVector scv, scv2; - scv2 = scv; - ioTest("iotest_vector.h5", scv, "SpinColourVector"); - SpinColourMatrix scm; - ioTest("iotest_matrix.h5", scm, "SpinColourMatrix"); +#define TensorWriteReadInnerNoInit( T ) \ + ioTest("iotest_"s + std::to_string(++TestNum) + "_" #T ".h5", t, #T); +#define TensorWriteReadInner( T ) SequentialInit( t ); TensorWriteReadInnerNoInit( T ) +#define TensorWriteRead( T ) { T t ; TensorWriteReadInner( T ) } +#define TensorWriteReadV(T, ... ) { T t( __VA_ARGS__ ); TensorWriteReadInner( T ) } +#define TensorWriteReadLarge( T ) { std::unique_ptr p{new T}; T &t{*p}; TensorWriteReadInnerNoInit(T) } - constexpr TestScalar Inc{1,-1}; +void EigenHdf5IOTest(void) +{ + unsigned int TestNum = 0; + using namespace std::string_literals; + TensorWriteRead( TensorSingle ) + TensorWriteReadV( TensorSimple, 1, 1, 1, 1, 1, 1 ) + TensorWriteReadV( TensorRank3, 6, 3, 2 ) + TensorWriteRead ( Tensor_9_4_2 ) + TensorWriteRead ( LSCTensor ) + TensorWriteReadLarge( PerambIOTestClass ) +#ifdef DEBUG + std::cout << "sizeof( LCMTensor ) = " << sizeof( LCMTensor ) / 1024 / 1024 << " MB" << std::endl; + TensorWriteReadLarge ( LCMTensor ) + // Also write > 4GB of complex numbers (I suspect this will fail inside Hdf5) { - using TestTensorSingle = Eigen::TensorFixedSize>; - TestTensorSingle ts; - ts(0) = 7; // lucky - ioTest("iotest_single.h5", ts, "Tensor_single"); + static constexpr size_t Num = 0x11000000; + std::cout << "Stress test: " << Num * sizeof( Grid_complex ) / 1024 / 1024 + << " MB array of complex" << std::endl; + using Stress = std::vector>; + Stress t (Num); + TensorWriteReadInnerNoInit( Stress ); } - - { - using TestTensorSimple = Eigen::Tensor, 6>; - TestTensorSimple ts(1,1,1,1,1,1); - ts(0,0,0,0,0,0) = Inc * 3.1415927; - ioTest("iotest_simple.h5", ts, "Tensor_simple"); - } - - TestTensor t(6,3,2); - TestScalar Val{Inc}; - for( int i = 0 ; i < t.dimension(0) ; i++) - for( int j = 0 ; j < t.dimension(1) ; j++) - for( int k = 0 ; k < t.dimension(2) ; k++) { - t(i,j,k) = Val; - Val += Inc; - } - ioTest("iotest_tensor.h5", t, "eigen_tensor_instance_name"); - //dump_tensor(t, "t"); - - // Now serialise a fixed size tensor - using FixedTensor = Eigen::TensorFixedSize>; - FixedTensor tf; - Val = Inc; - for( int i = 0 ; i < tf.dimension(0) ; i++) - for( int j = 0 ; j < tf.dimension(1) ; j++) - for( int k = 0 ; k < tf.dimension(2) ; k++) { - tf(i,j,k) = Val; - Val += Inc; - } - ioTest("iotest_tensor_fixed.h5", tf, "eigen_tensor_fixed_name"); - //dump_tensor(tf, "tf"); - - PerambIOTestClass o; - //dump_tensor(o.Perambulator, "Perambulator" ); - dump_tensor(o.shortRank5Tensor, "shortRank5Tensor"); - /*for_all( o.FixedCritter, [&](TestScalar &c, float f, const std::array &Dims ){ - c = TestScalar{f,-f}; - //std::cout << " a(" << Dims[0] << "," << Dims[1] << "," << Dims[2] << ")=" << c; - } ); - for( auto &z : o.aFixedCritter ) - for_all( z, [&](TestScalar &c, float f, const std::array &Dims ){ - c = TestScalar{f,-f}; - //std::cout << " a(" << Dims[0] << "," << Dims[1] << "," << Dims[2] << ")=" << c; - } );*/ - ioTest("iotest_object.h5", o, "PerambIOTestClass_object_instance_name"); - //DumpMemoryOrder(o.Perambulator); - - // Tensor of spin colour - LSCTensor l; - Val = 0; - for( int i = 0 ; i < l.dimension(0) ; i++) - for( int j = 0 ; j < l.dimension(1) ; j++) - for( int k = 0 ; k < l.dimension(2) ; k++) - for( int s = 0 ; s < Ns ; s++ ) - for( int c = 0 ; c < Nc ; c++ ) - { - l(i,j,k)()(s)(c) = Val; - Val += Inc; - } - ioTest("iotest_LSCTensor.h5", l, "LSCTensor_object_instance_name"); - //dump_tensor(l, "l"); - - // Tensor of spin colour - LCMTensor l2; - Val = 0; - for( int i = 0 ; i < l2.dimension(0) ; i++) - for( int j = 0 ; j < l2.dimension(1) ; j++) - for( int k = 0 ; k < l2.dimension(2) ; k++) - for( int l = 0 ; l < Ns ; l++ ) - for( int c = 0 ; c < Nc ; c++ ) - for( int c2 = 0 ; c2 < Nc ; c2++ ) - { - l2(i,j,k)(l)()(c,c2) = Val; - Val += Inc; - } - ioTest("iotest_LCMTensor.h5", l2, "LCMTensor_object_instance_name"); - //dump_tensor(l2, "l2"); +#endif } #endif