From 35f05e83cc3584c33180e148665f74a33a1ad256 Mon Sep 17 00:00:00 2001 From: gray95 Date: Thu, 30 Oct 2025 12:35:36 +0000 Subject: [PATCH 01/67] add element to ildg_format and put ildg_lfn in it's own message after ildg_binary_data --- Grid/parallelIO/IldgIO.h | 19 ++++++++++++++----- Grid/parallelIO/IldgIOtypes.h | 1 + 2 files changed, 15 insertions(+), 5 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 688d1ac9b8..a0bdccc99e 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -484,7 +484,7 @@ class GridLimeWriter : public BinaryIO err=limeWriterCloseRecord(LimeW); assert(err>=0); } //////////////////////////////////////// - // Write checksum element, propagaing forward from the BinaryIO + // Write checksum element, propagating forward from the BinaryIO // Always pair a checksum with a binary object, and close message //////////////////////////////////////// scidacChecksum checksum; @@ -620,7 +620,7 @@ class IldgWriter : public ScidacWriter { { uint64_t PayloadSize = LFN.size(); int err; - createLimeRecordHeader(ILDG_DATA_LFN, 0 , 0, PayloadSize); + createLimeRecordHeader(ILDG_DATA_LFN, 1 , 1, PayloadSize); err=limeWriteRecordData(const_cast(LFN.c_str()), &PayloadSize,LimeW); assert(err>=0); err=limeWriterCloseRecord(LimeW); assert(err>=0); } @@ -631,7 +631,7 @@ class IldgWriter : public ScidacWriter { // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// template - void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) + void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description, bool reducedStorage = false) { GridBase * grid = Umu.Grid(); typedef Lattice GaugeField; @@ -661,17 +661,25 @@ class IldgWriter : public ScidacWriter { ////////////////////////////////////////////////////// // Fill ILDG header data struct + // Needs to be able to figure out what type of field + // is being written and change ilgfmt.field accordngly ////////////////////////////////////////////////////// ildgFormat ildgfmt ; const std::string stNC = std::to_string( Nc ) ; ildgfmt.field = std::string("su"+stNC+"gauge"); + if ( reducedStorage ) { + ildgfmt.rows = Nc-1; + } else { + ildgfmt.rows = Nc; + } + if ( format == std::string("IEEE32BIG") ) { ildgfmt.precision = 32; } else { ildgfmt.precision = 64; } - ildgfmt.version = 1.0; + ildgfmt.version = 1.0; //update to 1.2 ildgfmt.lx = header.dimension[0]; ildgfmt.ly = header.dimension[1]; ildgfmt.lz = header.dimension[2]; @@ -704,8 +712,9 @@ class IldgWriter : public ScidacWriter { writeLimeObject(1,0,_scidacRecord,_scidacRecord.SerialisableClassName(),std::string(SCIDAC_PRIVATE_RECORD_XML)); writeLimeObject(0,0,info,info.SerialisableClassName(),std::string(SCIDAC_RECORD_XML)); writeLimeObject(0,0,ildgfmt,std::string("ildgFormat") ,std::string(ILDG_FORMAT)); // rec - writeLimeIldgLFN(header.ildg_lfn); // rec + //writeLimeIldgLFN(header.ildg_lfn); // rec writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum + writeLimeIldgLFN(header.ildg_lfn); // rec // limeDestroyWriter(LimeW); } }; diff --git a/Grid/parallelIO/IldgIOtypes.h b/Grid/parallelIO/IldgIOtypes.h index ddc0969c60..4fae9d9dac 100644 --- a/Grid/parallelIO/IldgIOtypes.h +++ b/Grid/parallelIO/IldgIOtypes.h @@ -150,6 +150,7 @@ struct ildgFormat : Serializable { double, version, std::string, field, int, precision, + int, rows, int, lx, int, ly, int, lz, From a73b6bce505cecd616538d013bd72a43965b9f6a Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 11 Nov 2025 17:31:11 +0000 Subject: [PATCH 02/67] readConfiguration is backwards compatible with ildg<1.2 lattices. writeConfiguration now has template args specifying the gauge group and whether to write to reduced storage format. Started working on a reduced storage test. --- Grid/parallelIO/IldgIO.h | 88 +++++++++++++++++++------ tests/IO/Test_ildg_reducedStorage.cc | 98 ++++++++++++++++++++++++++++ 2 files changed, 167 insertions(+), 19 deletions(-) create mode 100644 tests/IO/Test_ildg_reducedStorage.cc diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index a0bdccc99e..ed23e8058c 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -38,6 +38,10 @@ directory #include #include +//#include +//#include +//#include + //C-Lime is a must have for this functionality extern "C" { #include "lime.h" @@ -630,8 +634,8 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template - void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description, bool reducedStorage = false) + template + void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) { GridBase * grid = Umu.Grid(); typedef Lattice GaugeField; @@ -665,21 +669,38 @@ class IldgWriter : public ScidacWriter { // is being written and change ilgfmt.field accordngly ////////////////////////////////////////////////////// ildgFormat ildgfmt ; - const std::string stNC = std::to_string( Nc ) ; - ildgfmt.field = std::string("su"+stNC+"gauge"); + const std::string stNC = std::to_string(Nc) ; + + // check which gauge group is being written + is_su is_su; + is_sp is_sp; + + if constexpr ( is_su.value ) { + std::cout << GridLogMessage << "is_su returned TRUE" << std::endl; + ildgfmt.field = std::string("su"+stNC+"gauge"); + } else if constexpr ( is_sp.value ) { + std::cout << GridLogMessage << "is_sp returned TRUE" << std::endl; + ildgfmt.field = std::string("sp"+stNC+"gauge"); + } else { + static_assert(is_su.value || is_sp.value, "Can't tell whether SU or Sp is being written!"); + } + + // populate 'rows' element accordingly + if constexpr( reducedStorage && is_su.value ) { + ildgfmt.rows = Nc-1 ; + } else if constexpr( reducedStorage && is_sp.value ) { + ildgfmt.rows = Nc/2 ; + } else { + ildgfmt.rows = Nc ; + } - if ( reducedStorage ) { - ildgfmt.rows = Nc-1; - } else { - ildgfmt.rows = Nc; - } - if ( format == std::string("IEEE32BIG") ) { ildgfmt.precision = 32; } else { ildgfmt.precision = 64; } - ildgfmt.version = 1.0; //update to 1.2 + + ildgfmt.version = 1.2; ildgfmt.lx = header.dimension[0]; ildgfmt.ly = header.dimension[1]; ildgfmt.lz = header.dimension[2]; @@ -753,12 +774,12 @@ class IldgReader : public GridLimeReader { FieldNormMetaData FieldNormMetaData_; // track what we read from file - int found_ildgFormat =0; - int found_ildgLFN =0; - int found_scidacChecksum=0; - int found_usqcdInfo =0; - int found_ildgBinary =0; - int found_FieldMetaData =0; + int found_ildgFormat =0; + int found_ildgLFN =0; + int found_scidacChecksum =0; + int found_usqcdInfo =0; + int found_ildgBinary =0; + int found_FieldMetaData =0; int found_FieldNormMetaData =0; uint32_t nersc_csum; uint32_t scidac_csuma; @@ -788,20 +809,49 @@ class IldgReader : public GridLimeReader { // Copy out the string std::vector xmlc(nbytes+1,'\0'); limeReaderReadData((void *)&xmlc[0], &nbytes, LimeR); - // std::cout << GridLogMessage<< "Non binary record :" < +Author: Peter Boyle +Author: paboyle + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License along + with this program; if not, write to the Free Software Foundation, Inc., + 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + + See the full license in the file "LICENSE" in the top level distribution directory + *************************************************************************************/ + /* END LEGAL */ +#include + +using namespace std; +using namespace Grid; + ; + + +int main (int argc, char ** argv) +{ +#ifdef HAVE_LIME + Grid_init(&argc,&argv); + + std::cout <({45,12,81,9})); + sRNGa.SeedFixedIntegers(std::vector({45,12,81,9})); + std::cout < U(4,&Fine); + + SU::HotConfiguration(pRNGa,Umu); + + + FieldMetaData header; + + std::cout < + template void writeLimeLatticeBinaryObject(Lattice &field,std::string record_name,int control=BINARYIO_LEXICOGRAPHIC) { //////////////////////////////////////////////////////////////////// @@ -442,7 +442,7 @@ class GridLimeWriter : public BinaryIO //////////////////////////////////////////// // Create record header //////////////////////////////////////////// - typedef typename vobj::scalar_object sobj; + typedef typename GaugeUnMunger::out_type sobj; int err; uint32_t nersc_csum,scidac_csuma,scidac_csumb; uint64_t PayloadSize = sizeof(sobj) * grid->_gsites; @@ -451,9 +451,9 @@ class GridLimeWriter : public BinaryIO fflush(File); } - // std::cout << "W sizeof(sobj)" <_gsites<_gsites<(); - BinarySimpleMunger munge; - BinaryIO::writeLatticeObject(field, filename, munge, offset1, format,nersc_csum,scidac_csuma,scidac_csumb,control); + std::string format = getFormatString(); + + GaugeUnMunger unmunger; + + BinaryIO::writeLatticeObject(field,filename,unmunger,offset1,format,nersc_csum,scidac_csuma,scidac_csumb,control); /////////////////////////////////////////// // Wind forward and close the record @@ -634,13 +636,13 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template - void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) + template + void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) { GridBase * grid = Umu.Grid(); - typedef Lattice GaugeField; - typedef vLorentzColourMatrixD vobj; - typedef typename vobj::scalar_object sobj; + //typedef Lattice GaugeField; + //typedef vLorentzColourMatrixD vobj; + //typedef typename vobj::scalar_object sobj; //////////////////////////////////////// // fill the Grid header @@ -653,13 +655,13 @@ class IldgWriter : public ScidacWriter { stats Stats; Stats(Umu,header); - + + //header.floating_point = std::string("IEEE32BIG"); std::string format = header.floating_point; header.ensemble_id = description; header.ensemble_label = description; header.sequence_number = sequence; header.ildg_lfn = LFN; - assert ( (format == std::string("IEEE32BIG")) ||(format == std::string("IEEE64BIG")) ); @@ -671,29 +673,30 @@ class IldgWriter : public ScidacWriter { ildgFormat ildgfmt ; const std::string stNC = std::to_string(Nc) ; - // check which gauge group is being written - is_su is_su; - is_sp is_sp; + // find out which gauge group is going to be written + is_su su; + is_sp sp; - if constexpr ( is_su.value ) { + if constexpr ( su.value ) { std::cout << GridLogMessage << "is_su returned TRUE" << std::endl; ildgfmt.field = std::string("su"+stNC+"gauge"); - } else if constexpr ( is_sp.value ) { + } else if constexpr ( sp.value ) { std::cout << GridLogMessage << "is_sp returned TRUE" << std::endl; ildgfmt.field = std::string("sp"+stNC+"gauge"); } else { - static_assert(is_su.value || is_sp.value, "Can't tell whether SU or Sp is being written!"); + static_assert(su.value || sp.value, "Can't tell whether SU or Sp is being written!"); } // populate 'rows' element accordingly - if constexpr( reducedStorage && is_su.value ) { + if constexpr( reducedStorage && su.value ) { ildgfmt.rows = Nc-1 ; - } else if constexpr( reducedStorage && is_sp.value ) { + } else if constexpr( reducedStorage && sp.value ) { ildgfmt.rows = Nc/2 ; } else { - ildgfmt.rows = Nc ; + ildgfmt.rows = Nc; } - + + if ( format == std::string("IEEE32BIG") ) { ildgfmt.precision = 32; } else { @@ -734,7 +737,7 @@ class IldgWriter : public ScidacWriter { writeLimeObject(0,0,info,info.SerialisableClassName(),std::string(SCIDAC_RECORD_XML)); writeLimeObject(0,0,ildgfmt,std::string("ildgFormat") ,std::string(ILDG_FORMAT)); // rec //writeLimeIldgLFN(header.ildg_lfn); // rec - writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum + writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum writeLimeIldgLFN(header.ildg_lfn); // rec // limeDestroyWriter(LimeW); } @@ -757,8 +760,15 @@ class IldgReader : public GridLimeReader { typedef typename GaugeField::vector_object vobj; typedef typename vobj::scalar_object sobj; - typedef LorentzColourMatrixF fobj; - typedef LorentzColourMatrixD dobj; + // need all the types we could possibly read from, + // including the intermediate data types for reduced + // format lattices. + typedef LorentzColourMatrixF fobj; + typedef LorentzColourMatrixD dobj; + typedef LorentzColour2x3F fobjsuR; + typedef LorentzColour2x3D dobjsuR; + typedef LorentzColourNx2NF fobjspR; + typedef LorentzColourNx2ND dobjspR; GridBase *grid = Umu.Grid(); @@ -785,6 +795,11 @@ class IldgReader : public GridLimeReader { uint32_t scidac_csuma; uint32_t scidac_csumb; + // do we need to reconstruct a reduced format lattice? + bool reducedStorage; + bool is_su; + bool is_sp; + // Binary format std::string format; @@ -804,7 +819,7 @@ class IldgReader : public GridLimeReader { ////////////////////////////////////////////////////////////////// // If not BINARY_DATA read a string and parse ////////////////////////////////////////////////////////////////// - if ( strncmp(limeReaderType(LimeR), ILDG_BINARY_DATA,strlen(ILDG_BINARY_DATA) ) ) { + if ( strncmp(limeReaderType(LimeR), ILDG_BINARY_DATA,strlen(ILDG_BINARY_DATA) ) ) { // Copy out the string std::vector xmlc(nbytes+1,'\0'); @@ -826,31 +841,29 @@ class IldgReader : public GridLimeReader { if(doc.child("ildgFormat").child("rows")) { std::string rows = doc.child("ildgFormat").child("rows").child_value(); std::cout << GridLogMessage << "rows element found in ildg-format and is " << rows << " so this lattice could be ildg 1.2 compliant!" << std::endl; - } else { + } else { std::cout << GridLogMessage << "rows element is not present in ildg-format - let's add it to make cfg ildg 1.2 compliant when written out." << std::endl; pugi::xml_node ildgfmt = doc.child("ildgFormat"); const std::string stNC = std::to_string(Nc); ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); - //pugi::xml_node row = ildgfmt.insert_child_after("rows", ildgfmt.child("field")); - //row.append_child(pugi::node_pcdata).set_value(stNC.c_str()); // write result back into xmlstring std::ostringstream ss; doc.save(ss); xmlstring = ss.str(); - - //std::cout << GridLogMessage << "[inside if] size of xmlstring" << xmlstring.size() << std::endl; - //std::cout << GridLogMessage << "[inside if] size of xmlc" << xmlc.size() << std::endl; - //xmlc.resize(xmlstring.size()); - } + } XmlReader RD(xmlstring, true, ""); read(RD,"ildgFormat",ildgFormat_); - if ( ildgFormat_.precision == 64 ) format = std::string("IEEE64BIG"); - if ( ildgFormat_.precision == 32 ) format = std::string("IEEE32BIG"); + if ( ildgFormat_.precision == 64 ) { format = std::string("IEEE64BIG"); } + if ( ildgFormat_.precision == 32 ) { format = std::string("IEEE32BIG"); } + // do we need to reconstruct the full field with a munger? std::cout << GridLogMessage << "ildgFormat_.rows is " << ildgFormat_.rows << std::endl; + reducedStorage = (ildgFormat_.rows < Nc) ? true : false; + is_su = (!strncmp(ildgFormat_.field.c_str(),"su",2)) ? true : false; + is_sp = (!strncmp(ildgFormat_.field.c_str(),"sp",2)) ? true : false; assert( ildgFormat_.lx == dims[0]); assert( ildgFormat_.ly == dims[1]); @@ -907,14 +920,36 @@ class IldgReader : public GridLimeReader { // Binary data ///////////////////////////////// // std::cout << GridLogMessage << "ILDG Binary record found : " ILDG_BINARY_DATA << std::endl; + uint64_t offset= ftello(File); - if ( format == std::string("IEEE64BIG") ) { - GaugeSimpleMunger munge; - BinaryIO::readLatticeObject< vobj, dobj >(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else { - GaugeSimpleMunger munge; - BinaryIO::readLatticeObject< vobj, fobj >(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } + // read lattice into Umu + // depends on precision, + // reducedStorage, + // and whether reading su or sp fields + + if ( format == std::string("IEEE64BIG") ) { + if ( is_su && reducedStorage ) { + Gauge3x2munger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else if ( is_sp && reducedStorage ) { + GaugeSpmunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else { + GaugeSimpleMunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } + } else { + if ( is_su && reducedStorage ) { + Gauge3x2munger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else if ( is_sp && reducedStorage ) { + GaugeSpmunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else { + GaugeSimpleMunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } + } found_ildgBinary = 1; } diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 6b9d870855..d5071bee5d 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -209,7 +209,7 @@ inline void reconstruct3(LorentzColourMatrix & cm) cm(mu)()(1,0) = -adj(cm(mu)()(0,y)) ; cm(mu)()(1,1) = adj(cm(mu)()(0,x)) ; #else - const int x=0 , y=1 , z=2 ; // a little disinenuous labelling + const int x=0 , y=1 , z=2 ; // a little disingenuous labelling cm(mu)()(2,x) = adj(cm(mu)()(0,y)*cm(mu)()(1,z)-cm(mu)()(0,z)*cm(mu)()(1,y)); //x= yz-zy cm(mu)()(2,y) = adj(cm(mu)()(0,z)*cm(mu)()(1,x)-cm(mu)()(0,x)*cm(mu)()(1,z)); //y= zx-xz cm(mu)()(2,z) = adj(cm(mu)()(0,x)*cm(mu)()(1,y)-cm(mu)()(0,y)*cm(mu)()(1,x)); //z= xy-yx @@ -217,15 +217,32 @@ inline void reconstruct3(LorentzColourMatrix & cm) } } +inline void reconstructSp(LorentzColourMatrix & cm) { + assert( Nc%2==0 ); + int N = Nc/2; + for(int mu=0;mu using iLorentzColour2x3 = iVector, Nc-1>, Nd >; +template using iLorentzColourNx2N = iVector, Nc/2>, Nd >; typedef iLorentzColour2x3 LorentzColour2x3; typedef iLorentzColour2x3 LorentzColour2x3F; typedef iLorentzColour2x3 LorentzColour2x3D; +typedef iLorentzColourNx2N LorentzColourNx2N; +typedef iLorentzColourNx2N LorentzColourNx2NF; +typedef iLorentzColourNx2N LorentzColourNx2ND; + ///////////////////////////////////////////////////////////////////////////////// // Simple classes for precision conversion ///////////////////////////////////////////////////////////////////////////////// @@ -235,7 +252,7 @@ struct BinarySimpleUnmunger { typedef typename getPrecision::real_scalar_type sobj_stype; void operator()(sobj &in, fobj &out) { - // take word by word and transform accoding to the status + // take word by word and transform according to the status fobj_stype *out_buffer = (fobj_stype *)&out; sobj_stype *in_buffer = (sobj_stype *)∈ size_t fobj_words = sizeof(out) / sizeof(fobj_stype); @@ -273,7 +290,7 @@ struct GaugeSimpleMunger{ void operator()(fobj &in, sobj &out) { for (int mu = 0; mu < Nd; mu++) { for (int i = 0; i < Nc; i++) { - for (int j = 0; j < Nc; j++) { +for (int j = 0; j < Nc; j++) { out(mu)()(i, j) = in(mu)()(i, j); }} } @@ -341,5 +358,92 @@ struct Gauge3x2unmunger{ } }; +template +struct GaugeSpmunger{ + void operator() (fobj &in,sobj &out){ + for(int mu=0;mu +struct GaugeSpunmunger{ + void operator() (sobj &in,fobj &out){ + for(int mu=0;mu +struct GaugeUnMunger; + +// no group reduction +template +struct GaugeUnMunger +{ + using in_type = typename vobj::scalar_object; // LorentzColourMatrixF/D + using out_type = typename std::conditional::type; + BinarySimpleUnmunger unmunger; + + void operator() (in_type &in, out_type &out){ + unmunger(in,out); + } + +}; + +//group reduction for SU +template +struct GaugeUnMunger +{ + using in_type = typename vobj::scalar_object; // LorentzColourMatrixF/D + using tmp_type = typename std::conditional::type; + using out_type = typename std::conditional::type; + + BinarySimpleUnmunger binary_unmunger; + Gauge3x2unmunger gauge_unmunger; + + void operator() (in_type &in, out_type &out){ + tmp_type tmp; + binary_unmunger(in, tmp); + gauge_unmunger(tmp,out); + } + +}; + +//group reduction for Sp +template +struct GaugeUnMunger +{ + using in_type = typename vobj::scalar_object; // LorentzColourMatrixF/D + using tmp_type = typename std::conditional::type; + using out_type = typename std::conditional::type; + + BinarySimpleUnmunger binary_unmunger; + GaugeSpunmunger gauge_unmunger; + + void operator() (in_type &in, out_type &out){ + tmp_type tmp; + binary_unmunger(in, tmp); + gauge_unmunger(tmp, out); + } + +}; + NAMESPACE_END(Grid); From d5e51c8f1e164e0c6667e047c9a7833693747013 Mon Sep 17 00:00:00 2001 From: gray95 Date: Thu, 27 Nov 2025 20:08:48 +0000 Subject: [PATCH 04/67] writeConfiguration correctly writes out precision into ildg header. now the default is to write to double (with no group reduction). added Test_ildg_reducedfmt_io for testing writing/reading to reduced group format and single/double precision. --- Grid/parallelIO/IldgIO.h | 50 ++++++++--------- Grid/parallelIO/MetaData.h | 30 +++++----- ...dStorage.cc => Test_ildg_reducedfmt_io.cc} | 56 ++++++++++++++++--- 3 files changed, 87 insertions(+), 49 deletions(-) rename tests/IO/{Test_ildg_reducedStorage.cc => Test_ildg_reducedfmt_io.cc} (52%) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 57c83c34c1..5aaf7b4ecc 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -417,7 +417,7 @@ class GridLimeWriter : public BinaryIO // This routine is Collectively called by all nodes // in communicator used by the field.Grid() //////////////////////////////////////////////////// - template + template void writeLimeLatticeBinaryObject(Lattice &field,std::string record_name,int control=BINARYIO_LEXICOGRAPHIC) { //////////////////////////////////////////////////////////////////// @@ -442,7 +442,7 @@ class GridLimeWriter : public BinaryIO //////////////////////////////////////////// // Create record header //////////////////////////////////////////// - typedef typename GaugeUnMunger::out_type sobj; + typedef typename GaugeUnMunger::out_type sobj; int err; uint32_t nersc_csum,scidac_csuma,scidac_csumb; uint64_t PayloadSize = sizeof(sobj) * grid->_gsites; @@ -469,7 +469,7 @@ class GridLimeWriter : public BinaryIO /////////////////////////////////////////// std::string format = getFormatString(); - GaugeUnMunger unmunger; + GaugeUnMunger unmunger; BinaryIO::writeLatticeObject(field,filename,unmunger,offset1,format,nersc_csum,scidac_csuma,scidac_csumb,control); @@ -636,13 +636,11 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template + template void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) { + GridBase * grid = Umu.Grid(); - //typedef Lattice GaugeField; - //typedef vLorentzColourMatrixD vobj; - //typedef typename vobj::scalar_object sobj; //////////////////////////////////////// // fill the Grid header @@ -656,22 +654,16 @@ class IldgWriter : public ScidacWriter { stats Stats; Stats(Umu,header); - //header.floating_point = std::string("IEEE32BIG"); - std::string format = header.floating_point; header.ensemble_id = description; header.ensemble_label = description; header.sequence_number = sequence; header.ildg_lfn = LFN; - assert ( (format == std::string("IEEE32BIG")) - ||(format == std::string("IEEE64BIG")) ); - ////////////////////////////////////////////////////// // Fill ILDG header data struct - // Needs to be able to figure out what type of field - // is being written and change ilgfmt.field accordngly ////////////////////////////////////////////////////// ildgFormat ildgfmt ; const std::string stNC = std::to_string(Nc) ; + std::string format; // find out which gauge group is going to be written is_su su; @@ -695,13 +687,19 @@ class IldgWriter : public ScidacWriter { } else { ildgfmt.rows = Nc; } - - - if ( format == std::string("IEEE32BIG") ) { - ildgfmt.precision = 32; - } else { - ildgfmt.precision = 64; - } + // set fmt string for single/double precision + if constexpr( store_as_float) { + format = std::string("IEEE32BIG"); + header.floating_point = std::string("IEEE32BIG"); + ildgfmt.precision = 32; + } else { + format = std::string("IEEE64BIG"); + header.floating_point = std::string("IEEE64BIG"); + ildgfmt.precision = 64; + } + + std::cout << GridLogMessage << "format is " << format << std::endl; + std::cout << GridLogMessage << "header.floating_point is " << header.floating_point << std::endl; ildgfmt.version = 1.2; ildgfmt.lx = header.dimension[0]; @@ -737,7 +735,7 @@ class IldgWriter : public ScidacWriter { writeLimeObject(0,0,info,info.SerialisableClassName(),std::string(SCIDAC_RECORD_XML)); writeLimeObject(0,0,ildgfmt,std::string("ildgFormat") ,std::string(ILDG_FORMAT)); // rec //writeLimeIldgLFN(header.ildg_lfn); // rec - writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum + writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum writeLimeIldgLFN(header.ildg_lfn); // rec // limeDestroyWriter(LimeW); } @@ -753,11 +751,11 @@ class IldgReader : public GridLimeReader { // Else use ILDG MetaData object if present. // Else use SciDAC MetaData object if present. //////////////////////////////////////////////////////////////// - template - void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_) { + template + void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_) { - typedef Lattice GaugeField; - typedef typename GaugeField::vector_object vobj; + //typedef Lattice GaugeField; + //typedef typename GaugeField::vector_object vobj; typedef typename vobj::scalar_object sobj; // need all the types we could possibly read from, diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index d5071bee5d..6eaf160cd9 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -390,15 +390,15 @@ struct GaugeSpunmunger{ // the munger (for reading) can only be chosen // at runtime since you don't know beforehand // which format the fields were written in. -template +template struct GaugeUnMunger; // no group reduction -template -struct GaugeUnMunger +template +struct GaugeUnMunger { - using in_type = typename vobj::scalar_object; // LorentzColourMatrixF/D - using out_type = typename std::conditional::type; + using in_type = typename vobj::scalar_object; + using out_type = typename std::conditional::type; BinarySimpleUnmunger unmunger; void operator() (in_type &in, out_type &out){ @@ -408,12 +408,12 @@ struct GaugeUnMunger }; //group reduction for SU -template -struct GaugeUnMunger +template +struct GaugeUnMunger { - using in_type = typename vobj::scalar_object; // LorentzColourMatrixF/D - using tmp_type = typename std::conditional::type; - using out_type = typename std::conditional::type; + using in_type = typename vobj::scalar_object; + using tmp_type = typename std::conditional::type; + using out_type = typename std::conditional::type; BinarySimpleUnmunger binary_unmunger; Gauge3x2unmunger gauge_unmunger; @@ -427,12 +427,12 @@ struct GaugeUnMunger }; //group reduction for Sp -template -struct GaugeUnMunger +template +struct GaugeUnMunger { - using in_type = typename vobj::scalar_object; // LorentzColourMatrixF/D - using tmp_type = typename std::conditional::type; - using out_type = typename std::conditional::type; + using in_type = typename vobj::scalar_object; + using tmp_type = typename std::conditional::type; + using out_type = typename std::conditional::type; BinarySimpleUnmunger binary_unmunger; GaugeSpunmunger gauge_unmunger; diff --git a/tests/IO/Test_ildg_reducedStorage.cc b/tests/IO/Test_ildg_reducedfmt_io.cc similarity index 52% rename from tests/IO/Test_ildg_reducedStorage.cc rename to tests/IO/Test_ildg_reducedfmt_io.cc index 2c3c8d1f81..6afd3742d1 100644 --- a/tests/IO/Test_ildg_reducedStorage.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -43,8 +43,10 @@ int main (int argc, char ** argv) auto simd_layout = GridDefaultSimd(4,vComplex::Nsimd()); auto mpi_layout = GridDefaultMpi(); + //std::vector latt_size ({48,48,48,96}); + //std::vector latt_size ({32,32,32,32}); Coordinate latt_size ({8,8,8,16}); - Coordinate clatt_size ({4,4,4,8}); + Coordinate clatt_size ({2,2,2,4}); int orthodir=3; int orthosz =latt_size[orthodir]; @@ -65,22 +67,23 @@ int main (int argc, char ** argv) LatticeGaugeField Umu_diff(&Fine); LatticeGaugeField Umu_saved(&Fine); - std::vector U(4,&Fine); + //std::vector U(4,&Fine); - SU::HotConfiguration(pRNGa,Umu); - + //SU::HotConfiguration(pRNGa,Umu); + Sp::HotConfiguration(pRNGa,Umu); + //const bool reducedStorage = true; + //const bool store_as_float = false; //write to disk in double or single precision FieldMetaData header; std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); - Umu_saved = Umu; std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter1.close(); + Umu_saved = Umu; + std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter2.close(); + Umu_saved = Umu; + std::cout <_gsites<_gsites<(); - GaugeUnMunger unmunger; + GaugeUnMunger unmunger; BinaryIO::writeLatticeObject(field,filename,unmunger,offset1,format,nersc_csum,scidac_csuma,scidac_csumb,control); @@ -636,7 +636,7 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template + template void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) { @@ -680,22 +680,24 @@ class IldgWriter : public ScidacWriter { } // populate 'rows' element accordingly - if constexpr( reducedStorage && su.value ) { + if constexpr( reducedGrpStorage && su.value ) { ildgfmt.rows = Nc-1 ; - } else if constexpr( reducedStorage && sp.value ) { + } else if constexpr( reducedGrpStorage && sp.value ) { ildgfmt.rows = Nc/2 ; } else { ildgfmt.rows = Nc; } // set fmt string for single/double precision - if constexpr( store_as_float) { + if constexpr( fp_fmt == FP_FMT::IEEE32BIG ) { format = std::string("IEEE32BIG"); header.floating_point = std::string("IEEE32BIG"); ildgfmt.precision = 32; - } else { + } else if constexpr ( fp_fmt == FP_FMT::IEEE64BIG ) { format = std::string("IEEE64BIG"); header.floating_point = std::string("IEEE64BIG"); ildgfmt.precision = 64; + } else { + static_assert(1, "Grid does NOT recognise FP_FMT enumerator"); } std::cout << GridLogMessage << "format is " << format << std::endl; @@ -734,8 +736,7 @@ class IldgWriter : public ScidacWriter { writeLimeObject(1,0,_scidacRecord,_scidacRecord.SerialisableClassName(),std::string(SCIDAC_PRIVATE_RECORD_XML)); writeLimeObject(0,0,info,info.SerialisableClassName(),std::string(SCIDAC_RECORD_XML)); writeLimeObject(0,0,ildgfmt,std::string("ildgFormat") ,std::string(ILDG_FORMAT)); // rec - //writeLimeIldgLFN(header.ildg_lfn); // rec - writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum + writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum writeLimeIldgLFN(header.ildg_lfn); // rec // limeDestroyWriter(LimeW); } @@ -794,7 +795,8 @@ class IldgReader : public GridLimeReader { uint32_t scidac_csumb; // do we need to reconstruct a reduced format lattice? - bool reducedStorage; + FP_FMT fp_fmt; + bool reducedGrpStorage; bool is_su; bool is_sp; @@ -854,15 +856,18 @@ class IldgReader : public GridLimeReader { XmlReader RD(xmlstring, true, ""); read(RD,"ildgFormat",ildgFormat_); - if ( ildgFormat_.precision == 64 ) { format = std::string("IEEE64BIG"); } - if ( ildgFormat_.precision == 32 ) { format = std::string("IEEE32BIG"); } + if ( ildgFormat_.precision == 64 ) { format = std::string("IEEE64BIG"); fp_fmt = FP_FMT::IEEE64BIG; } + if ( ildgFormat_.precision == 32 ) { format = std::string("IEEE32BIG"); fp_fmt = FP_FMT::IEEE32BIG;} // do we need to reconstruct the full field with a munger? std::cout << GridLogMessage << "ildgFormat_.rows is " << ildgFormat_.rows << std::endl; - reducedStorage = (ildgFormat_.rows < Nc) ? true : false; + reducedGrpStorage = (ildgFormat_.rows < Nc) ? true : false; is_su = (!strncmp(ildgFormat_.field.c_str(),"su",2)) ? true : false; is_sp = (!strncmp(ildgFormat_.field.c_str(),"sp",2)) ? true : false; + // enum logic in here + //if ( !reducedGrpStorage && + assert( ildgFormat_.lx == dims[0]); assert( ildgFormat_.ly == dims[1]); assert( ildgFormat_.lz == dims[2]); @@ -922,32 +927,32 @@ class IldgReader : public GridLimeReader { uint64_t offset= ftello(File); // read lattice into Umu // depends on precision, - // reducedStorage, + // reducedGrpStorage, // and whether reading su or sp fields - if ( format == std::string("IEEE64BIG") ) { - if ( is_su && reducedStorage ) { + if ( fp_fmt == FP_FMT::IEEE64BIG ) { + if ( is_su && reducedGrpStorage ) { Gauge3x2munger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && reducedStorage ) { + } else if ( is_sp && reducedGrpStorage ) { GaugeSpmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else { GaugeSimpleMunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } - } else { - if ( is_su && reducedStorage ) { - Gauge3x2munger munge; + } else if ( fp_fmt == FP_FMT::IEEE32BIG) { + if ( is_su && reducedGrpStorage ) { + Gauge3x2munger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && reducedStorage ) { - GaugeSpmunger munge; + } else if ( is_sp && reducedGrpStorage ) { + GaugeSpmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else { - GaugeSimpleMunger munge; + GaugeSimpleMunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } - } + } else { assert("FP_FMT not found"); } found_ildgBinary = 1; } diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 6eaf160cd9..ed001274e3 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -239,7 +239,7 @@ typedef iLorentzColour2x3 LorentzColour2x3; typedef iLorentzColour2x3 LorentzColour2x3F; typedef iLorentzColour2x3 LorentzColour2x3D; -typedef iLorentzColourNx2N LorentzColourNx2N; +//typedef iLorentzColourNx2N LorentzColourNx2N; // better to always specify F or D? typedef iLorentzColourNx2N LorentzColourNx2NF; typedef iLorentzColourNx2N LorentzColourNx2ND; @@ -385,11 +385,9 @@ struct GaugeSpunmunger{ } }; +/* // this struct is used to choose the appropriate // unmunger (for writing) at compile time. -// the munger (for reading) can only be chosen -// at runtime since you don't know beforehand -// which format the fields were written in. template struct GaugeUnMunger; @@ -444,6 +442,68 @@ struct GaugeUnMunger } }; +*/ + +// this enum should cover all the cases for reading and writing +enum struct FP_FMT { IEEE64BIG, IEEE32BIG }; + +// this struct is used to choose the appropriate +// unmunger (for writing) at compile time. +template +struct GaugeUnMunger; + +// no group reduction +template +struct GaugeUnMunger +{ + using in_type = typename vobj::scalar_object; + using out_type = typename std::tuple_element_t(fmt),std::tuple>; + BinarySimpleUnmunger unmunger; + + void operator() (in_type &in, out_type &out){ + unmunger(in,out); + } + +}; + +//group reduction for SU +template +struct GaugeUnMunger +{ + using in_type = typename vobj::scalar_object; + using tmp_type = typename std::tuple_element_t(fmt),std::tuple>; + using out_type = typename std::tuple_element_t(fmt),std::tuple>; + + BinarySimpleUnmunger binary_unmunger; + Gauge3x2unmunger gauge_unmunger; + + void operator() (in_type &in, out_type &out){ + tmp_type tmp; + binary_unmunger(in, tmp); + gauge_unmunger(tmp,out); + } + +}; + +//group reduction for Sp +template +struct GaugeUnMunger +{ + using in_type = typename vobj::scalar_object; + using tmp_type = typename std::tuple_element_t(fmt),std::tuple>; + using out_type = typename std::tuple_element_t(fmt),std::tuple>; + + BinarySimpleUnmunger binary_unmunger; + GaugeSpunmunger gauge_unmunger; + + void operator() (in_type &in, out_type &out){ + tmp_type tmp; + binary_unmunger(in, tmp); + gauge_unmunger(tmp, out); + } + +}; + NAMESPACE_END(Grid); From d0a7f3ac69762d29900011eb65e70fd31129b7a5 Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 3 Dec 2025 14:34:05 +0000 Subject: [PATCH 06/67] cleaned up some comments and set test_ildg_reducedfmt to use SU by default. --- Grid/parallelIO/IldgIO.h | 72 +++++++++++------------ Grid/parallelIO/MetaData.h | 90 +++++------------------------ tests/IO/Test_ildg_reducedfmt_io.cc | 26 +++++---- 3 files changed, 65 insertions(+), 123 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 22dc5f7d69..bcc38ddbc1 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -417,7 +417,7 @@ class GridLimeWriter : public BinaryIO // This routine is Collectively called by all nodes // in communicator used by the field.Grid() //////////////////////////////////////////////////// - template + template void writeLimeLatticeBinaryObject(Lattice &field,std::string record_name,int control=BINARYIO_LEXICOGRAPHIC) { //////////////////////////////////////////////////////////////////// @@ -442,9 +442,10 @@ class GridLimeWriter : public BinaryIO //////////////////////////////////////////// // Create record header //////////////////////////////////////////// - typedef typename GaugeUnMunger::out_type sobj; int err; uint32_t nersc_csum,scidac_csuma,scidac_csumb; + + typedef typename GaugeUnMunger::out_type sobj; uint64_t PayloadSize = sizeof(sobj) * grid->_gsites; if ( boss_node ) { createLimeRecordHeader(record_name, 0, 0, PayloadSize); @@ -469,7 +470,7 @@ class GridLimeWriter : public BinaryIO /////////////////////////////////////////// std::string format = getFormatString(); - GaugeUnMunger unmunger; + GaugeUnMunger unmunger; BinaryIO::writeLatticeObject(field,filename,unmunger,offset1,format,nersc_csum,scidac_csuma,scidac_csumb,control); @@ -636,7 +637,7 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template + template void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) { @@ -663,7 +664,7 @@ class IldgWriter : public ScidacWriter { ////////////////////////////////////////////////////// ildgFormat ildgfmt ; const std::string stNC = std::to_string(Nc) ; - std::string format; + //std::string format; // find out which gauge group is going to be written is_su su; @@ -680,29 +681,29 @@ class IldgWriter : public ScidacWriter { } // populate 'rows' element accordingly - if constexpr( reducedGrpStorage && su.value ) { + if constexpr( matrix_fmt==MATRIX_FMT::REDUCED && su.value ) { ildgfmt.rows = Nc-1 ; - } else if constexpr( reducedGrpStorage && sp.value ) { + } else if constexpr( matrix_fmt==MATRIX_FMT::REDUCED && sp.value ) { ildgfmt.rows = Nc/2 ; - } else { + } else if constexpr( matrix_fmt==MATRIX_FMT::FULL ) { ildgfmt.rows = Nc; + } else { + static_assert(1, "MATRIX_FMT might be incorrectly specified"); } + // set fmt string for single/double precision if constexpr( fp_fmt == FP_FMT::IEEE32BIG ) { - format = std::string("IEEE32BIG"); + //format = std::string("IEEE32BIG"); header.floating_point = std::string("IEEE32BIG"); ildgfmt.precision = 32; } else if constexpr ( fp_fmt == FP_FMT::IEEE64BIG ) { - format = std::string("IEEE64BIG"); + //format = std::string("IEEE64BIG"); header.floating_point = std::string("IEEE64BIG"); ildgfmt.precision = 64; } else { static_assert(1, "Grid does NOT recognise FP_FMT enumerator"); } - std::cout << GridLogMessage << "format is " << format << std::endl; - std::cout << GridLogMessage << "header.floating_point is " << header.floating_point << std::endl; - ildgfmt.version = 1.2; ildgfmt.lx = header.dimension[0]; ildgfmt.ly = header.dimension[1]; @@ -736,7 +737,7 @@ class IldgWriter : public ScidacWriter { writeLimeObject(1,0,_scidacRecord,_scidacRecord.SerialisableClassName(),std::string(SCIDAC_PRIVATE_RECORD_XML)); writeLimeObject(0,0,info,info.SerialisableClassName(),std::string(SCIDAC_RECORD_XML)); writeLimeObject(0,0,ildgfmt,std::string("ildgFormat") ,std::string(ILDG_FORMAT)); // rec - writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum + writeLimeLatticeBinaryObject(Umu,std::string(ILDG_BINARY_DATA)); // Closes message with checksum writeLimeIldgLFN(header.ildg_lfn); // rec // limeDestroyWriter(LimeW); } @@ -755,13 +756,10 @@ class IldgReader : public GridLimeReader { template void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_) { - //typedef Lattice GaugeField; - //typedef typename GaugeField::vector_object vobj; typedef typename vobj::scalar_object sobj; // need all the types we could possibly read from, - // including the intermediate data types for reduced - // format lattices. + // including the intermediate data types for reduced format lattices. typedef LorentzColourMatrixF fobj; typedef LorentzColourMatrixD dobj; typedef LorentzColour2x3F fobjsuR; @@ -794,9 +792,9 @@ class IldgReader : public GridLimeReader { uint32_t scidac_csuma; uint32_t scidac_csumb; - // do we need to reconstruct a reduced format lattice? + // do we need to reconstruct the lattice from a reduced format? FP_FMT fp_fmt; - bool reducedGrpStorage; + MATRIX_FMT matrix_fmt; bool is_su; bool is_sp; @@ -834,15 +832,15 @@ class IldgReader : public GridLimeReader { if ( !strncmp(limeReaderType(LimeR), ILDG_FORMAT,strlen(ILDG_FORMAT)) ) { - // if there is no rows element in ildg-format insert such an element with value of Nc + // if no rows element in ildg-format insert such an element with value of Nc pugi::xml_document doc; doc.load_string(xmlstring.c_str()); if(doc.child("ildgFormat").child("rows")) { std::string rows = doc.child("ildgFormat").child("rows").child_value(); - std::cout << GridLogMessage << "rows element found in ildg-format and is " << rows << " so this lattice could be ildg 1.2 compliant!" << std::endl; + std::cout << GridLogMessage << "rows element found in ildg-format and is " << rows << " so this lattice might be ildg 1.2 compliant." << std::endl; } else { - std::cout << GridLogMessage << "rows element is not present in ildg-format - let's add it to make cfg ildg 1.2 compliant when written out." << std::endl; + std::cout << GridLogMessage << "rows element is not present in ildg-format - adding it after ." << std::endl; pugi::xml_node ildgfmt = doc.child("ildgFormat"); const std::string stNC = std::to_string(Nc); ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); @@ -856,18 +854,16 @@ class IldgReader : public GridLimeReader { XmlReader RD(xmlstring, true, ""); read(RD,"ildgFormat",ildgFormat_); - if ( ildgFormat_.precision == 64 ) { format = std::string("IEEE64BIG"); fp_fmt = FP_FMT::IEEE64BIG; } - if ( ildgFormat_.precision == 32 ) { format = std::string("IEEE32BIG"); fp_fmt = FP_FMT::IEEE32BIG;} + // what if data is in another format ieee64[?] + if ( ildgFormat_.precision == 64 ) { fp_fmt = FP_FMT::IEEE64BIG; } + if ( ildgFormat_.precision == 32 ) { fp_fmt = FP_FMT::IEEE32BIG; } // do we need to reconstruct the full field with a munger? std::cout << GridLogMessage << "ildgFormat_.rows is " << ildgFormat_.rows << std::endl; - reducedGrpStorage = (ildgFormat_.rows < Nc) ? true : false; + matrix_fmt = (ildgFormat_.rows < Nc) ? MATRIX_FMT::REDUCED : MATRIX_FMT::FULL; is_su = (!strncmp(ildgFormat_.field.c_str(),"su",2)) ? true : false; is_sp = (!strncmp(ildgFormat_.field.c_str(),"sp",2)) ? true : false; - // enum logic in here - //if ( !reducedGrpStorage && - assert( ildgFormat_.lx == dims[0]); assert( ildgFormat_.ly == dims[1]); assert( ildgFormat_.lz == dims[2]); @@ -925,16 +921,13 @@ class IldgReader : public GridLimeReader { // std::cout << GridLogMessage << "ILDG Binary record found : " ILDG_BINARY_DATA << std::endl; uint64_t offset= ftello(File); - // read lattice into Umu - // depends on precision, - // reducedGrpStorage, - // and whether reading su or sp fields - + if ( fp_fmt == FP_FMT::IEEE64BIG ) { - if ( is_su && reducedGrpStorage ) { + format = std::string("IEEE64BIG"); + if ( is_su && matrix_fmt==MATRIX_FMT::REDUCED ) { Gauge3x2munger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && reducedGrpStorage ) { + } else if ( is_sp && matrix_fmt==MATRIX_FMT::REDUCED ) { GaugeSpmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else { @@ -942,17 +935,18 @@ class IldgReader : public GridLimeReader { BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } } else if ( fp_fmt == FP_FMT::IEEE32BIG) { - if ( is_su && reducedGrpStorage ) { + format = std::string("IEEE32BIG"); + if ( is_su && matrix_fmt==MATRIX_FMT::REDUCED ) { Gauge3x2munger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && reducedGrpStorage ) { + } else if ( is_sp && matrix_fmt==MATRIX_FMT::REDUCED ) { GaugeSpmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else { GaugeSimpleMunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } - } else { assert("FP_FMT not found"); } + } else { assert("Unable to determine which readLatticeObject function template to instantiate."); } found_ildgBinary = 1; } diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index ed001274e3..755cd8ddbf 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -217,7 +217,8 @@ inline void reconstruct3(LorentzColourMatrix & cm) } } -inline void reconstructSp(LorentzColourMatrix & cm) { +inline void reconstructSp(LorentzColourMatrix & cm) +{ assert( Nc%2==0 ); int N = Nc/2; for(int mu=0;mu LorentzColour2x3; typedef iLorentzColour2x3 LorentzColour2x3F; typedef iLorentzColour2x3 LorentzColour2x3D; -//typedef iLorentzColourNx2N LorentzColourNx2N; // better to always specify F or D? typedef iLorentzColourNx2N LorentzColourNx2NF; typedef iLorentzColourNx2N LorentzColourNx2ND; @@ -385,79 +385,21 @@ struct GaugeSpunmunger{ } }; -/* -// this struct is used to choose the appropriate -// unmunger (for writing) at compile time. -template -struct GaugeUnMunger; - -// no group reduction -template -struct GaugeUnMunger -{ - using in_type = typename vobj::scalar_object; - using out_type = typename std::conditional::type; - BinarySimpleUnmunger unmunger; - - void operator() (in_type &in, out_type &out){ - unmunger(in,out); - } - -}; - -//group reduction for SU -template -struct GaugeUnMunger -{ - using in_type = typename vobj::scalar_object; - using tmp_type = typename std::conditional::type; - using out_type = typename std::conditional::type; - - BinarySimpleUnmunger binary_unmunger; - Gauge3x2unmunger gauge_unmunger; - - void operator() (in_type &in, out_type &out){ - tmp_type tmp; - binary_unmunger(in, tmp); - gauge_unmunger(tmp,out); - } - -}; - -//group reduction for Sp -template -struct GaugeUnMunger -{ - using in_type = typename vobj::scalar_object; - using tmp_type = typename std::conditional::type; - using out_type = typename std::conditional::type; - - BinarySimpleUnmunger binary_unmunger; - GaugeSpunmunger gauge_unmunger; - - void operator() (in_type &in, out_type &out){ - tmp_type tmp; - binary_unmunger(in, tmp); - gauge_unmunger(tmp, out); - } - -}; -*/ - -// this enum should cover all the cases for reading and writing +// add other options? eg. IEEE64LITTLE ? enum struct FP_FMT { IEEE64BIG, IEEE32BIG }; +enum struct MATRIX_FMT { FULL, REDUCED }; // this struct is used to choose the appropriate // unmunger (for writing) at compile time. -template +template struct GaugeUnMunger; // no group reduction -template -struct GaugeUnMunger +template +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; - using out_type = typename std::tuple_element_t(fmt),std::tuple>; + using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; BinarySimpleUnmunger unmunger; void operator() (in_type &in, out_type &out){ @@ -467,12 +409,12 @@ struct GaugeUnMunger }; //group reduction for SU -template -struct GaugeUnMunger +template +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; - using tmp_type = typename std::tuple_element_t(fmt),std::tuple>; - using out_type = typename std::tuple_element_t(fmt),std::tuple>; + using tmp_type = typename std::tuple_element_t(fp_fmt),std::tuple>; + using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; BinarySimpleUnmunger binary_unmunger; Gauge3x2unmunger gauge_unmunger; @@ -486,12 +428,12 @@ struct GaugeUnMunger }; //group reduction for Sp -template -struct GaugeUnMunger +template +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; - using tmp_type = typename std::tuple_element_t(fmt),std::tuple>; - using out_type = typename std::tuple_element_t(fmt),std::tuple>; + using tmp_type = typename std::tuple_element_t(fp_fmt),std::tuple>; + using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; BinarySimpleUnmunger binary_unmunger; GaugeSpunmunger gauge_unmunger; diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index 6afd3742d1..f7a9a7c89a 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -69,20 +69,26 @@ int main (int argc, char ** argv) //std::vector U(4,&Fine); - //SU::HotConfiguration(pRNGa,Umu); - Sp::HotConfiguration(pRNGa,Umu); - //const bool reducedStorage = true; - //const bool store_as_float = false; //write to disk in double or single precision + SU::HotConfiguration(pRNGa,Umu); + using grpName = GroupName::SU; + + //Sp::HotConfiguration(pRNGa,Umu); + //using grpName = GroupName::Sp; // ideally want to automate this + + FP_FMT const fmt64 = FP_FMT::IEEE64BIG; + FP_FMT const fmt32 = FP_FMT::IEEE32BIG; + MATRIX_FMT const noGrpRdc = MATRIX_FMT::FULL; + MATRIX_FMT const GrpRdc = MATRIX_FMT::REDUCED; FieldMetaData header; std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); Umu_saved = Umu; std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter1.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter1.close(); Umu_saved = Umu; std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter2.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter2.close(); Umu_saved = Umu; std::cout < Date: Mon, 12 Jan 2026 19:22:07 +0000 Subject: [PATCH 07/67] started updating manual.rst and commit some trailing changes from before xmas --- Grid/parallelIO/IldgIO.h | 96 +++++++++++++---------------- documentation/manual.rst | 15 +++-- tests/IO/Test_ildg_reducedfmt_io.cc | 41 +++++++----- 3 files changed, 78 insertions(+), 74 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index bcc38ddbc1..f755870125 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -38,10 +38,6 @@ directory #include #include -//#include -//#include -//#include - //C-Lime is a must have for this functionality extern "C" { #include "lime.h" @@ -349,7 +345,7 @@ class GridLimeWriter : public BinaryIO boss_node = isboss; } void open(const std::string &_filename) { - filename = _filename; + filename= _filename; if ( boss_node ) { File = fopen(filename.c_str(), "w"); LimeW = limeCreateWriter(File); assert(LimeW != NULL ); @@ -417,7 +413,7 @@ class GridLimeWriter : public BinaryIO // This routine is Collectively called by all nodes // in communicator used by the field.Grid() //////////////////////////////////////////////////// - template + template void writeLimeLatticeBinaryObject(Lattice &field,std::string record_name,int control=BINARYIO_LEXICOGRAPHIC) { //////////////////////////////////////////////////////////////////// @@ -444,7 +440,6 @@ class GridLimeWriter : public BinaryIO //////////////////////////////////////////// int err; uint32_t nersc_csum,scidac_csuma,scidac_csumb; - typedef typename GaugeUnMunger::out_type sobj; uint64_t PayloadSize = sizeof(sobj) * grid->_gsites; if ( boss_node ) { @@ -452,9 +447,9 @@ class GridLimeWriter : public BinaryIO fflush(File); } - // std::cout << "W sizeof(sobj)" <_gsites<_gsites< unmunger; - BinaryIO::writeLatticeObject(field,filename,unmunger,offset1,format,nersc_csum,scidac_csuma,scidac_csumb,control); + BinaryIO::writeLatticeObject(field, filename, unmunger, offset1, format, nersc_csum, scidac_csuma, scidac_csumb, control); /////////////////////////////////////////// // Wind forward and close the record @@ -637,10 +632,9 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template - void writeConfiguration(Lattice &Umu,int sequence,std::string LFN,std::string description) + template + void writeConfiguration(Lattice &Umu, int sequence, std::string LFN, std::string description) { - GridBase * grid = Umu.Grid(); //////////////////////////////////////// @@ -654,7 +648,7 @@ class IldgWriter : public ScidacWriter { stats Stats; Stats(Umu,header); - + header.ensemble_id = description; header.ensemble_label = description; header.sequence_number = sequence; @@ -663,8 +657,7 @@ class IldgWriter : public ScidacWriter { // Fill ILDG header data struct ////////////////////////////////////////////////////// ildgFormat ildgfmt ; - const std::string stNC = std::to_string(Nc) ; - //std::string format; + const std::string stNC = std::to_string( Nc ) ; // find out which gauge group is going to be written is_su su; @@ -681,27 +674,25 @@ class IldgWriter : public ScidacWriter { } // populate 'rows' element accordingly - if constexpr( matrix_fmt==MATRIX_FMT::REDUCED && su.value ) { + if constexpr( matrix_fmt==MatrixFormat::REDUCED && su.value ) { ildgfmt.rows = Nc-1 ; - } else if constexpr( matrix_fmt==MATRIX_FMT::REDUCED && sp.value ) { + } else if constexpr( matrix_fmt==MatrixFormat::REDUCED && sp.value ) { ildgfmt.rows = Nc/2 ; - } else if constexpr( matrix_fmt==MATRIX_FMT::FULL ) { + } else if constexpr( matrix_fmt==MatrixFormat::FULL ) { ildgfmt.rows = Nc; } else { - static_assert(1, "MATRIX_FMT might be incorrectly specified"); + static_assert(1, "Unknown MatrixFormat specified"); } // set fmt string for single/double precision - if constexpr( fp_fmt == FP_FMT::IEEE32BIG ) { - //format = std::string("IEEE32BIG"); + if constexpr( fp_fmt == FloatingPointFormat::IEEE32BIG ) { header.floating_point = std::string("IEEE32BIG"); ildgfmt.precision = 32; - } else if constexpr ( fp_fmt == FP_FMT::IEEE64BIG ) { - //format = std::string("IEEE64BIG"); + } else if constexpr ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { header.floating_point = std::string("IEEE64BIG"); ildgfmt.precision = 64; } else { - static_assert(1, "Grid does NOT recognise FP_FMT enumerator"); + static_assert(1, "Unknown FloatingPointFormat specified"); } ildgfmt.version = 1.2; @@ -762,10 +753,10 @@ class IldgReader : public GridLimeReader { // including the intermediate data types for reduced format lattices. typedef LorentzColourMatrixF fobj; typedef LorentzColourMatrixD dobj; - typedef LorentzColour2x3F fobjsuR; - typedef LorentzColour2x3D dobjsuR; - typedef LorentzColourNx2NF fobjspR; - typedef LorentzColourNx2ND dobjspR; + typedef LorentzColour2x3F fobjsuR; + typedef LorentzColour2x3D dobjsuR; + typedef LorentzColourNx2NF fobjspR; + typedef LorentzColourNx2ND dobjspR; GridBase *grid = Umu.Grid(); @@ -781,23 +772,24 @@ class IldgReader : public GridLimeReader { FieldNormMetaData FieldNormMetaData_; // track what we read from file - int found_ildgFormat =0; - int found_ildgLFN =0; - int found_scidacChecksum =0; - int found_usqcdInfo =0; - int found_ildgBinary =0; - int found_FieldMetaData =0; + int found_ildgFormat =0; + int found_ildgLFN =0; + int found_scidacChecksum=0; + int found_usqcdInfo =0; + int found_ildgBinary =0; + int found_FieldMetaData =0; int found_FieldNormMetaData =0; uint32_t nersc_csum; uint32_t scidac_csuma; uint32_t scidac_csumb; - // do we need to reconstruct the lattice from a reduced format? - FP_FMT fp_fmt; - MATRIX_FMT matrix_fmt; + // these variables store information about the binary data that is read + // from the ildg-format header. if matrix_fmt==REDUCED + // then Grid will reconstruct the full matrix using the appropriate munger. + FloatingPointFormat fp_fmt; + MatrixFormat matrix_fmt; bool is_su; bool is_sp; - // Binary format std::string format; @@ -817,19 +809,17 @@ class IldgReader : public GridLimeReader { ////////////////////////////////////////////////////////////////// // If not BINARY_DATA read a string and parse ////////////////////////////////////////////////////////////////// - if ( strncmp(limeReaderType(LimeR), ILDG_BINARY_DATA,strlen(ILDG_BINARY_DATA) ) ) { + if ( strncmp(limeReaderType(LimeR), ILDG_BINARY_DATA,strlen(ILDG_BINARY_DATA) ) ) { // Copy out the string std::vector xmlc(nbytes+1,'\0'); limeReaderReadData((void *)&xmlc[0], &nbytes, LimeR); - //std::cout << GridLogMessage<< "Non binary record :" < munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && matrix_fmt==MATRIX_FMT::REDUCED ) { + } else if ( is_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else { GaugeSimpleMunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } - } else if ( fp_fmt == FP_FMT::IEEE32BIG) { + } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { format = std::string("IEEE32BIG"); - if ( is_su && matrix_fmt==MATRIX_FMT::REDUCED ) { + if ( is_su && matrix_fmt==MatrixFormat::REDUCED ) { Gauge3x2munger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && matrix_fmt==MATRIX_FMT::REDUCED ) { + } else if ( is_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else { diff --git a/documentation/manual.rst b/documentation/manual.rst index 0c463672c1..03b1e933e3 100644 --- a/documentation/manual.rst +++ b/documentation/manual.rst @@ -2104,21 +2104,26 @@ These are specialised to SciDAC writers, introducing facilities for generating t void readScidacFieldRecord(Lattice &field,userRecord &_userRecord); }; -They are also specialised to ILDG format writers, available and defined only for Gauge configurations:: +They are also specialised to ILDG format writers, available and defined only for Gauge configurations. Following the publication of Rev. 1.2 of the ILDG Binary File Format by the ILDG Metadata Working Group the ILDG format writer and reader have been updated to be able to write SU(2/3) and Sp(2N) lattice configurations in reduced formats and to read and reconstruct reduced format lattices.:: class IldgWriter : public ScidacWriter { - template - void writeConfiguration(Lattice > &Umu,int sequence,std::string LFN,std::string description); + template + void writeConfiguration(Lattice &Umu, int sequence, std::string LFN, std::string description) }; class IldgReader : public GridLimeReader { - template - void readConfiguration(Lattice > &Umu, FieldMetaData &FieldMetaData_) ; + + template + void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_); }; +The template parameters *group_name*, *MatrixFormat*, *FloatingPointFormat* specify the exact format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced or not reduced, and two options for the floating point format, 32 bit and 64 bit precision. If a lattice is not written out in reduced format then *group_name* has no practical effect on the writing process. + +Put a table here describing these constructs in more detail? + **Implementation detail** The lattice field routines internally use the above Binary routines to write the bulk data at an offset diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index f7a9a7c89a..c87cbf1183 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -9,6 +9,7 @@ Author: Azusa Yamaguchi Author: Peter Boyle Author: paboyle +Author: Gaurav Ray This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by @@ -31,7 +32,24 @@ Author: paboyle using namespace std; using namespace Grid; - ; + + +////////////////////////////////////////////// +// +// this template function generates a +// LatticeGaugeField of the chosen gauge +// group, passed as a template argument. +// +////////////////////////////////////////////// +template +inline void generateFieldHotConfiguration( LatticeGaugeField &Umu, GridParallelRNG &pRNG ) { + if constexpr( std::is_same_v == true ) { + SU::HotConfiguration(pRNG,Umu); + } else if constexpr ( std::is_same_v == true ) { + Sp::HotConfiguration(pRNG,Umu); + } else { static_assert(true, "Grid does not recognise the gauge group"); } + +} int main (int argc, char ** argv) @@ -43,8 +61,6 @@ int main (int argc, char ** argv) auto simd_layout = GridDefaultSimd(4,vComplex::Nsimd()); auto mpi_layout = GridDefaultMpi(); - //std::vector latt_size ({48,48,48,96}); - //std::vector latt_size ({32,32,32,32}); Coordinate latt_size ({8,8,8,16}); Coordinate clatt_size ({2,2,2,4}); int orthodir=3; @@ -54,9 +70,7 @@ int main (int argc, char ** argv) GridCartesian Coarse(clatt_size,simd_layout,mpi_layout); GridParallelRNG pRNGa(&Fine); - GridParallelRNG pRNGb(&Fine); GridSerialRNG sRNGa; - GridSerialRNG sRNGb; std::cout <({45,12,81,9})); @@ -67,18 +81,13 @@ int main (int argc, char ** argv) LatticeGaugeField Umu_diff(&Fine); LatticeGaugeField Umu_saved(&Fine); - //std::vector U(4,&Fine); - - SU::HotConfiguration(pRNGa,Umu); - using grpName = GroupName::SU; - - //Sp::HotConfiguration(pRNGa,Umu); - //using grpName = GroupName::Sp; // ideally want to automate this + using grpName = GroupName::Sp; + generateFieldHotConfiguration(Umu, pRNGa); - FP_FMT const fmt64 = FP_FMT::IEEE64BIG; - FP_FMT const fmt32 = FP_FMT::IEEE32BIG; - MATRIX_FMT const noGrpRdc = MATRIX_FMT::FULL; - MATRIX_FMT const GrpRdc = MATRIX_FMT::REDUCED; + FloatingPointFormat const fmt64 = FloatingPointFormat::IEEE64BIG; + FloatingPointFormat const fmt32 = FloatingPointFormat::IEEE32BIG; + MatrixFormat const noGrpRdc = MatrixFormat::FULL; + MatrixFormat const GrpRdc = MatrixFormat::REDUCED; FieldMetaData header; From a4357fa56e733c9ff06b0ad87fb8f8daae0d943a Mon Sep 17 00:00:00 2001 From: gray95 Date: Thu, 15 Jan 2026 19:28:07 +0000 Subject: [PATCH 08/67] one more draft pr to go --- Grid/parallelIO/IldgIO.h | 184 +++++++++++++++++----------- Grid/parallelIO/MetaData.h | 19 ++- tests/IO/Test_ildg_reducedfmt_io.cc | 52 ++++---- 3 files changed, 146 insertions(+), 109 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index f755870125..127e2d652c 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -659,41 +659,38 @@ class IldgWriter : public ScidacWriter { ildgFormat ildgfmt ; const std::string stNC = std::to_string( Nc ) ; - // find out which gauge group is going to be written - is_su su; - is_sp sp; - - if constexpr ( su.value ) { - std::cout << GridLogMessage << "is_su returned TRUE" << std::endl; - ildgfmt.field = std::string("su"+stNC+"gauge"); - } else if constexpr ( sp.value ) { - std::cout << GridLogMessage << "is_sp returned TRUE" << std::endl; - ildgfmt.field = std::string("sp"+stNC+"gauge"); - } else { - static_assert(su.value || sp.value, "Can't tell whether SU or Sp is being written!"); - } + // use the gauge group to populate the 'field' element of ildg header + if constexpr ( is_su::value ) { + std::cout << GridLogMessage << "writing SU(" << stNC << ") field" << std::endl; + ildgfmt.field = std::string("su"+stNC+"gauge"); + } else if constexpr ( is_sp::value ) { + std::cout << GridLogMessage << "writing Sp(" << stNC << ") field" << std::endl; + ildgfmt.field = std::string("sp"+stNC+"gauge"); + } else { + static_assert(1, "Unrecognised group; unable to determine field tag"); + } - // populate 'rows' element accordingly - if constexpr( matrix_fmt==MatrixFormat::REDUCED && su.value ) { - ildgfmt.rows = Nc-1 ; - } else if constexpr( matrix_fmt==MatrixFormat::REDUCED && sp.value ) { - ildgfmt.rows = Nc/2 ; + // populate 'rows' element of ildg header + if constexpr( matrix_fmt==MatrixFormat::REDUCED && is_su::value ) { + ildgfmt.rows = Nc-1 ; + } else if constexpr( matrix_fmt==MatrixFormat::REDUCED && is_sp::value ) { + ildgfmt.rows = Nc/2 ; } else if constexpr( matrix_fmt==MatrixFormat::FULL ) { - ildgfmt.rows = Nc; - } else { - static_assert(1, "Unknown MatrixFormat specified"); - } + ildgfmt.rows = Nc; + } else { + static_assert(1, "Unknown MatrixFormat specified"); + } // set fmt string for single/double precision - if constexpr( fp_fmt == FloatingPointFormat::IEEE32BIG ) { - header.floating_point = std::string("IEEE32BIG"); - ildgfmt.precision = 32; - } else if constexpr ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { - header.floating_point = std::string("IEEE64BIG"); - ildgfmt.precision = 64; - } else { - static_assert(1, "Unknown FloatingPointFormat specified"); - } + if constexpr( fp_fmt == FloatingPointFormat::IEEE32BIG ) { + header.floating_point = std::string("IEEE32BIG"); + ildgfmt.precision = 32; + } else if constexpr ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { + header.floating_point = std::string("IEEE64BIG"); + ildgfmt.precision = 64; + } else { + static_assert(1, "Unknown FloatingPointFormat specified"); + } ildgfmt.version = 1.2; ildgfmt.lx = header.dimension[0]; @@ -788,8 +785,8 @@ class IldgReader : public GridLimeReader { // then Grid will reconstruct the full matrix using the appropriate munger. FloatingPointFormat fp_fmt; MatrixFormat matrix_fmt; - bool is_su; - bool is_sp; + bool is_grp_su; + bool is_grp_sp; // Binary format std::string format; @@ -822,37 +819,47 @@ class IldgReader : public GridLimeReader { std::string xmlstring(&xmlc[0]); if ( !strncmp(limeReaderType(LimeR), ILDG_FORMAT,strlen(ILDG_FORMAT)) ) { - // if no rows element in ildg-format insert such an element with value of Nc - pugi::xml_document doc; - doc.load_string(xmlstring.c_str()); - - if(doc.child("ildgFormat").child("rows")) { - std::string rows = doc.child("ildgFormat").child("rows").child_value(); - std::cout << GridLogMessage << "rows element found in ildg-format and is " << rows << " so this lattice might be ildg 1.2 compliant." << std::endl; - } else { - std::cout << GridLogMessage << "rows element is not present in ildg-format - adding it after ." << std::endl; - pugi::xml_node ildgfmt = doc.child("ildgFormat"); - const std::string stNC = std::to_string(Nc); - ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); - - // write result back into xmlstring - std::ostringstream ss; - doc.save(ss); - xmlstring = ss.str(); - } + // if no rows element in ildg-format insert such an element with value = Nc + pugi::xml_document doc; + doc.load_string(xmlstring.c_str()); + + if(doc.child("ildgFormat").child("rows")) { + std::string rows = doc.child("ildgFormat").child("rows").child_value(); + std::cout << GridLogMessage << "rows element found in ildg-format and is " << rows << " so this lattice might be ildg 1.2 compliant." << std::endl; + } else { + std::cout << GridLogMessage << "rows element is not present in ildg-format - adding it after ." << std::endl; + pugi::xml_node ildgfmt = doc.child("ildgFormat"); + const std::string stNC = std::to_string(Nc); + ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); + + // write result back into xmlstring + std::ostringstream ss; + doc.save(ss); + xmlstring = ss.str(); + } XmlReader RD(xmlstring, true, ""); read(RD,"ildgFormat",ildgFormat_); - // what if data is in another format ieee64[?] if ( ildgFormat_.precision == 64 ) { fp_fmt = FloatingPointFormat::IEEE64BIG; } if ( ildgFormat_.precision == 32 ) { fp_fmt = FloatingPointFormat::IEEE32BIG; } - // do we need to reconstruct the full field with a munger? + // set vars to tell Grid if it needs to reconstruct the full field. std::cout << GridLogMessage << "ildgFormat_.rows is " << ildgFormat_.rows << std::endl; matrix_fmt = (ildgFormat_.rows < Nc) ? MatrixFormat::REDUCED : MatrixFormat::FULL; - is_su = (!strncmp(ildgFormat_.field.c_str(),"su",2)) ? true : false; - is_sp = (!strncmp(ildgFormat_.field.c_str(),"sp",2)) ? true : false; + is_grp_su = (!strncmp(ildgFormat_.field.c_str(),"su",2)) ? true : false; + is_grp_sp = (!strncmp(ildgFormat_.field.c_str(),"sp",2)) ? true : false; + + // check rows and Nc are related in the way we expect for each gauge group. + if (matrix_fmt==MatrixFormat::REDUCED && is_grp_su) { + assert( ildgFormat_.rows == Nc-1 ); + } + if (matrix_fmt==MatrixFormat::REDUCED && is_grp_sp) { + assert( ildgFormat_.rows == Nc/2 ); + } + if (matrix_fmt==MatrixFormat::FULL) { + assert( ildgFormat_.rows == Nc ); + } assert( ildgFormat_.lx == dims[0]); assert( ildgFormat_.ly == dims[1]); @@ -911,33 +918,64 @@ class IldgReader : public GridLimeReader { // std::cout << GridLogMessage << "ILDG Binary record found : " ILDG_BINARY_DATA << std::endl; uint64_t offset= ftello(File); - - if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { + + if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { format = std::string("IEEE64BIG"); - if ( is_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSpmunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else { - GaugeSimpleMunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { + Gauge3x2munger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { + GaugeSpmunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else { + GaugeSimpleMunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } + } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { + format = std::string("IEEE32BIG"); + if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { + Gauge3x2munger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { + GaugeSpmunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else { + GaugeSimpleMunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } + } else { assert("Unable to determine which readLatticeObject function template to instantiate."); } + + +/* if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { + format = std::string("IEEE64BIG"); + if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { + Gauge3x2munger munge; + using diskobj = dobjsuR; + } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { + GaugeSpmunger munge; + using diskobj = dobjspR; + } else { + GaugeSimpleMunger munge; + using diskobj = dobj; + } } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { - format = std::string("IEEE32BIG"); - if ( is_su && matrix_fmt==MatrixFormat::REDUCED ) { + format = std::string("IEEE32BIG"); + if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { Gauge3x2munger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_sp && matrix_fmt==MatrixFormat::REDUCED ) { + using diskobj = fobjsuR; + } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + using diskobj = fobjspR; } else { - GaugeSimpleMunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + GaugeSimpleMunger munge; + using diskobj = fobj; } } else { assert("Unable to determine which readLatticeObject function template to instantiate."); } + // read lattice into Umu using the appropriate munger. + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); +*/ + found_ildgBinary = 1; } diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 755cd8ddbf..9c2b1036ed 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -385,18 +385,17 @@ struct GaugeSpunmunger{ } }; -// add other options? eg. IEEE64LITTLE ? -enum struct FP_FMT { IEEE64BIG, IEEE32BIG }; -enum struct MATRIX_FMT { FULL, REDUCED }; +enum struct FloatingPointFormat { IEEE64BIG, IEEE32BIG }; +enum struct MatrixFormat { FULL, REDUCED }; // this struct is used to choose the appropriate // unmunger (for writing) at compile time. -template +template struct GaugeUnMunger; // no group reduction -template -struct GaugeUnMunger +template +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; @@ -409,8 +408,8 @@ struct GaugeUnMunger }; //group reduction for SU -template -struct GaugeUnMunger +template +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; using tmp_type = typename std::tuple_element_t(fp_fmt),std::tuple>; @@ -428,8 +427,8 @@ struct GaugeUnMunger }; //group reduction for Sp -template -struct GaugeUnMunger +template +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; using tmp_type = typename std::tuple_element_t(fp_fmt),std::tuple>; diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index c87cbf1183..47c7775bc0 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -33,25 +33,32 @@ Author: Gaurav Ray using namespace std; using namespace Grid; - ////////////////////////////////////////////// // // this template function generates a // LatticeGaugeField of the chosen gauge -// group, passed as a template argument. +// group passed as a template argument. // ////////////////////////////////////////////// template -inline void generateFieldHotConfiguration( LatticeGaugeField &Umu, GridParallelRNG &pRNG ) { +LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Fine, std::vector seed ) { + + GridParallelRNG pRNGa(&Fine); + std::cout < == true ) { - SU::HotConfiguration(pRNG,Umu); + SU::HotConfiguration(pRNGa,Umu); } else if constexpr ( std::is_same_v == true ) { - Sp::HotConfiguration(pRNG,Umu); - } else { static_assert(true, "Grid does not recognise the gauge group"); } + Sp::HotConfiguration(pRNGa,Umu); + } else { static_assert(true, "Grid does not recognise the gauge group"); } + return Umu; } - int main (int argc, char ** argv) { #ifdef HAVE_LIME @@ -62,32 +69,25 @@ int main (int argc, char ** argv) auto simd_layout = GridDefaultSimd(4,vComplex::Nsimd()); auto mpi_layout = GridDefaultMpi(); Coordinate latt_size ({8,8,8,16}); - Coordinate clatt_size ({2,2,2,4}); - int orthodir=3; - int orthosz =latt_size[orthodir]; - - GridCartesian Fine(latt_size,simd_layout,mpi_layout); - GridCartesian Coarse(clatt_size,simd_layout,mpi_layout); - GridParallelRNG pRNGa(&Fine); - GridSerialRNG sRNGa; + std::vector seed0 = {15,16,23,42}; + std::vector seed1 = {7,10,19,95}; + std::vector seed2 = {55,34,23,13}; - std::cout <({45,12,81,9})); - sRNGa.SeedFixedIntegers(std::vector({45,12,81,9})); - std::cout <(Umu, pRNGa); + LatticeGaugeField Umu = generateHotFieldConfiguration(Fine, seed0); + LatticeGaugeField Umu_diff = generateHotFieldConfiguration(Fine, seed1); + LatticeGaugeField Umu_saved = generateHotFieldConfiguration(Fine, seed2); + + // define enums for different field formats FloatingPointFormat const fmt64 = FloatingPointFormat::IEEE64BIG; FloatingPointFormat const fmt32 = FloatingPointFormat::IEEE32BIG; - MatrixFormat const noGrpRdc = MatrixFormat::FULL; - MatrixFormat const GrpRdc = MatrixFormat::REDUCED; + MatrixFormat const noGrpRdc = MatrixFormat::FULL; + MatrixFormat const GrpRdc = MatrixFormat::REDUCED; FieldMetaData header; From 71e362010f7ebd959fde31ec3dfbadcc5324d6c1 Mon Sep 17 00:00:00 2001 From: gray95 Date: Mon, 19 Jan 2026 17:04:26 +0000 Subject: [PATCH 09/67] added comments --- Grid/parallelIO/MetaData.h | 81 ++++++++++++++++++++++---------------- 1 file changed, 48 insertions(+), 33 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 9c2b1036ed..cad9f086df 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -217,6 +217,17 @@ inline void reconstruct3(LorentzColourMatrix & cm) } } +//////////////////////////////////////////////////////////////// +// this function takes a reduced format +// Sp(2N) field with N rows and 2N columns +// and reconstructs the full 2Nx2N matrix. +// the full matrix has block structure: +// | A B | +// |-B* A*| +// where A and B are NxN matrices. +// see appendix A.2 of +// https://www-zeuthen.desy.de/apewww/ILDG/specifications/ildg-file-format-1.2.pdf +//////////////////////////////////////////////////////////////// inline void reconstructSp(LorentzColourMatrix & cm) { assert( Nc%2==0 ); @@ -234,12 +245,12 @@ inline void reconstructSp(LorentzColourMatrix & cm) // Some data types for intermediate storage //////////////////////////////////////////////////////////////////////////////// template using iLorentzColour2x3 = iVector, Nc-1>, Nd >; -template using iLorentzColourNx2N = iVector, Nc/2>, Nd >; - typedef iLorentzColour2x3 LorentzColour2x3; typedef iLorentzColour2x3 LorentzColour2x3F; typedef iLorentzColour2x3 LorentzColour2x3D; +// intermediate types for Sp(2N) fields +template using iLorentzColourNx2N = iVector, Nc/2>, Nd >; typedef iLorentzColourNx2N LorentzColourNx2NF; typedef iLorentzColourNx2N LorentzColourNx2ND; @@ -363,86 +374,90 @@ struct GaugeSpmunger{ void operator() (fobj &in,sobj &out){ for(int mu=0;mu struct GaugeSpunmunger{ void operator() (sobj &in,fobj &out){ for(int mu=0;mu no group reduction - treat SU and Sp the same +// > group reduction for SU fields +// > group reduction for Sp fields +///////////////////////////////////////////////////////// template struct GaugeUnMunger; // no group reduction template -struct GaugeUnMunger +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; - using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; - BinarySimpleUnmunger unmunger; + using out_type = typename std::tuple_element_t(fp_fmt), std::tuple>; + BinarySimpleUnmunger unmunger; void operator() (in_type &in, out_type &out){ - unmunger(in,out); + unmunger(in,out); } - }; //group reduction for SU template -struct GaugeUnMunger +struct GaugeUnMunger { using in_type = typename vobj::scalar_object; - using tmp_type = typename std::tuple_element_t(fp_fmt),std::tuple>; - using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; + using tmp_type = typename std::tuple_element_t(fp_fmt), std::tuple>; + using out_type = typename std::tuple_element_t(fp_fmt), std::tuple>; - BinarySimpleUnmunger binary_unmunger; - Gauge3x2unmunger gauge_unmunger; + BinarySimpleUnmunger binary_unmunger; + Gauge3x2unmunger gauge_unmunger; void operator() (in_type &in, out_type &out){ - tmp_type tmp; + tmp_type tmp; binary_unmunger(in, tmp); - gauge_unmunger(tmp,out); + gauge_unmunger(tmp,out); } - }; //group reduction for Sp template -struct GaugeUnMunger +struct GaugeUnMunger { - using in_type = typename vobj::scalar_object; - using tmp_type = typename std::tuple_element_t(fp_fmt),std::tuple>; - using out_type = typename std::tuple_element_t(fp_fmt),std::tuple>; + using in_type = typename vobj::scalar_object; + using tmp_type = typename std::tuple_element_t(fp_fmt), std::tuple>; + using out_type = typename std::tuple_element_t(fp_fmt), std::tuple>; - BinarySimpleUnmunger binary_unmunger; - GaugeSpunmunger gauge_unmunger; + BinarySimpleUnmunger binary_unmunger; + GaugeSpunmunger gauge_unmunger; void operator() (in_type &in, out_type &out){ - tmp_type tmp; + tmp_type tmp; binary_unmunger(in, tmp); - gauge_unmunger(tmp, out); + gauge_unmunger(tmp, out); } - }; From 5adfc779dce0f0d669523e0dbab2dd035ae4746a Mon Sep 17 00:00:00 2001 From: gray95 Date: Mon, 19 Jan 2026 17:11:05 +0000 Subject: [PATCH 10/67] prepping for a git rebase --- tests/IO/Test_ildg_reducedfmt_io.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index 47c7775bc0..e545bae184 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -77,7 +77,7 @@ int main (int argc, char ** argv) GridCartesian Fine(latt_size,simd_layout,mpi_layout); // set the gauge group - using grpName = GroupName::Sp; + using grpName = GroupName::SU; LatticeGaugeField Umu = generateHotFieldConfiguration(Fine, seed0); LatticeGaugeField Umu_diff = generateHotFieldConfiguration(Fine, seed1); From b67aa6d3f863041361e406c16e6f20ef31efd3b6 Mon Sep 17 00:00:00 2001 From: gray95 Date: Mon, 19 Jan 2026 20:09:52 +0000 Subject: [PATCH 11/67] helper function to wrap BinaryIO::readLatticeObject --- Grid/parallelIO/IldgIO.h | 143 ++++++++++++++++++--------------------- 1 file changed, 67 insertions(+), 76 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 127e2d652c..dc05fcc830 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -440,7 +440,7 @@ class GridLimeWriter : public BinaryIO //////////////////////////////////////////// int err; uint32_t nersc_csum,scidac_csuma,scidac_csumb; - typedef typename GaugeUnMunger::out_type sobj; + typedef typename GaugeUnMunger::out_type sobj; uint64_t PayloadSize = sizeof(sobj) * grid->_gsites; if ( boss_node ) { createLimeRecordHeader(record_name, 0, 0, PayloadSize); @@ -465,7 +465,7 @@ class GridLimeWriter : public BinaryIO /////////////////////////////////////////// std::string format = getFormatString(); - GaugeUnMunger unmunger; + GaugeUnMunger unmunger; BinaryIO::writeLatticeObject(field, filename, unmunger, offset1, format, nersc_csum, scidac_csuma, scidac_csumb, control); @@ -734,6 +734,62 @@ class IldgWriter : public ScidacWriter { class IldgReader : public GridLimeReader { public: + ////////////////////////////////////////////////////////////////// + // this helper function wraps the logic for choosing + // the correct munger and intermediate lattice data type + // before then reading a lattice field from disk. + // This is a runtime choice as Grid won't know which munger/itype + // to use until it has read the header of a lattice cfg. Therefore + // templating this function on fp_fmt etc. does not work. + // At present we use a rather cumbersome if statement with 6 branches. + // The benefit of organising IldgReader this way is it makes + // readConfiguration clearer and more readable. + ////////////////////////////////////////////////////////////////// + template + void readLatticeBinaryObject(Lattice &Umu, std::string filename, FloatingPointFormat fp_fmt, MatrixFormat matrix_fmt, bool is_grp_su, bool is_grp_sp, uint64_t &offset, uint32_t &nersc_csum, uint32_t &scidac_csuma, uint32_t &scidac_csumb) + { + + typedef typename vobj::scalar_object sobj; + // need all the types we could possibly read from + // including the intermediate data types for + // reduced format lattices and single/double precision + typedef LorentzColourMatrixF fobj; + typedef LorentzColourMatrixD dobj; + typedef LorentzColour2x3F fobjsuR; + typedef LorentzColour2x3D dobjsuR; + typedef LorentzColourNx2NF fobjspR; + typedef LorentzColourNx2ND dobjspR; + + std::string format; + + if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { + format = std::string("IEEE64BIG"); + if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { + Gauge3x2munger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { + GaugeSpmunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else { + GaugeSimpleMunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } + } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { + format = std::string("IEEE32BIG"); + if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { + Gauge3x2munger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { + GaugeSpmunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } else { + GaugeSimpleMunger munge; + BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); + } + } else { assert("Unable to determine which readLatticeObject function template to instantiate."); } + + } + //////////////////////////////////////////////////////////////// // Read either Grid/SciDAC/ILDG configuration // Don't require scidac records EXCEPT checksum @@ -744,17 +800,6 @@ class IldgReader : public GridLimeReader { template void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_) { - typedef typename vobj::scalar_object sobj; - - // need all the types we could possibly read from, - // including the intermediate data types for reduced format lattices. - typedef LorentzColourMatrixF fobj; - typedef LorentzColourMatrixD dobj; - typedef LorentzColour2x3F fobjsuR; - typedef LorentzColour2x3D dobjsuR; - typedef LorentzColourNx2NF fobjspR; - typedef LorentzColourNx2ND dobjspR; - GridBase *grid = Umu.Grid(); Coordinate dims = Umu.Grid()->FullDimensions(); @@ -780,13 +825,13 @@ class IldgReader : public GridLimeReader { uint32_t scidac_csuma; uint32_t scidac_csumb; - // these variables store information about the binary data that is read - // from the ildg-format header. if matrix_fmt==REDUCED - // then Grid will reconstruct the full matrix using the appropriate munger. - FloatingPointFormat fp_fmt; - MatrixFormat matrix_fmt; - bool is_grp_su; - bool is_grp_sp; + // these variables store information about the binary data that is read + // from the ildg-format header. if matrix_fmt==REDUCED + // then Grid will reconstruct the full matrix using the appropriate munger. + FloatingPointFormat fp_fmt; + MatrixFormat matrix_fmt; + bool is_grp_su; + bool is_grp_sp; // Binary format std::string format; @@ -919,64 +964,10 @@ class IldgReader : public GridLimeReader { uint64_t offset= ftello(File); - if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { - format = std::string("IEEE64BIG"); - if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSpmunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else { - GaugeSimpleMunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } - } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { - format = std::string("IEEE32BIG"); - if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSpmunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } else { - GaugeSimpleMunger munge; - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); - } - } else { assert("Unable to determine which readLatticeObject function template to instantiate."); } - - -/* if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { - format = std::string("IEEE64BIG"); - if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; - using diskobj = dobjsuR; - } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSpmunger munge; - using diskobj = dobjspR; - } else { - GaugeSimpleMunger munge; - using diskobj = dobj; - } - } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { - format = std::string("IEEE32BIG"); - if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; - using diskobj = fobjsuR; - } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSpmunger munge; - using diskobj = fobjspR; - } else { - GaugeSimpleMunger munge; - using diskobj = fobj; - } - } else { assert("Unable to determine which readLatticeObject function template to instantiate."); } - - // read lattice into Umu using the appropriate munger. - BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); -*/ + readLatticeBinaryObject(Umu, filename, fp_fmt, matrix_fmt, is_grp_su, is_grp_sp, offset, nersc_csum, scidac_csuma, scidac_csumb); found_ildgBinary = 1; + } } From cd3a6ff7aeb3fa87cd74c72ea46f2f68b3fd9b21 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 20 Jan 2026 21:40:02 +0000 Subject: [PATCH 12/67] Test.cc refactored with the core code in its own function that is templated on the gauge group. --- tests/IO/Test_ildg_reducedfmt_io.cc | 148 +++++++++++++--------------- 1 file changed, 69 insertions(+), 79 deletions(-) diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index e545bae184..f873fc6c23 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -34,21 +34,19 @@ using namespace std; using namespace Grid; ////////////////////////////////////////////// -// -// this template function generates a +// this template function returns a // LatticeGaugeField of the chosen gauge // group passed as a template argument. -// ////////////////////////////////////////////// template -LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Fine, std::vector seed ) { +LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Grid, std::vector seed = {83,8,34,59} ) { - GridParallelRNG pRNGa(&Fine); + GridParallelRNG pRNGa(&Grid); std::cout < == true ) { SU::HotConfiguration(pRNGa,Umu); @@ -59,94 +57,86 @@ LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Fine, std::vecto return Umu; } -int main (int argc, char ** argv) -{ -#ifdef HAVE_LIME - Grid_init(&argc,&argv); - - std::cout < seed0 = {15,16,23,42}; - std::vector seed1 = {7,10,19,95}; - std::vector seed2 = {55,34,23,13}; - - GridCartesian Fine(latt_size,simd_layout,mpi_layout); - - // set the gauge group - using grpName = GroupName::SU; - - LatticeGaugeField Umu = generateHotFieldConfiguration(Fine, seed0); - LatticeGaugeField Umu_diff = generateHotFieldConfiguration(Fine, seed1); - LatticeGaugeField Umu_saved = generateHotFieldConfiguration(Fine, seed2); - - // define enums for different field formats - FloatingPointFormat const fmt64 = FloatingPointFormat::IEEE64BIG; - FloatingPointFormat const fmt32 = FloatingPointFormat::IEEE32BIG; - MatrixFormat const noGrpRdc = MatrixFormat::FULL; - MatrixFormat const GrpRdc = MatrixFormat::REDUCED; +/////////////////////////////////////////////////////////////// +// this template function generates writes a lattice +// field of a given gaugeGroup to disk. It can write in +// reduced format and single precision depending on +// the values of matrix_fmt and fp_fmt. +/////////////////////////////////////////////////////////////// +template +void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, std::string file) { + + // if the IldgReader class is not able to handle + // the binary format do nothing and exit the function + if constexpr( std::is_same_v && N%2==1) { + std::cout <= 4" << std::endl; + std::cout < && matrix_fmt==MatrixFormat::REDUCED && N>3) + { + std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); - Umu_saved = Umu; + + LatticeGaugeField Umu_saved = generateHotFieldConfiguration(Grid); + std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); - _IldgWriter1.close(); - Umu_saved = Umu; - std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); - _IldgWriter2.close(); - Umu_saved = Umu; - std::cout < seed0 = {15,16,23,42}; + std::vector seed1 = {95,63,71,66}; + + GridCartesian Grid(latt_size,simd_layout,mpi_layout); + + LatticeGaugeField Umu = generateHotFieldConfiguration(Grid, seed0); + LatticeGaugeField UmuSp = generateHotFieldConfiguration(Grid, seed1); + + // write/read SU fields + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_full64.4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_rdc64.4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_full32.4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_rdc32.4000"); + + // write/read Sp fields + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_full64.4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_rdc64.4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_full32.4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_rdc32.4000"); Grid_finalize(); #endif From 050f9236414c3501a646008a00996ca120a892e4 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 20 Jan 2026 22:05:29 +0000 Subject: [PATCH 13/67] added/cleaned up some comments --- Grid/parallelIO/IldgIO.h | 20 ++++++++++---------- Grid/parallelIO/MetaData.h | 10 +++++++--- 2 files changed, 17 insertions(+), 13 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index dc05fcc830..221398301d 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -737,22 +737,22 @@ class IldgReader : public GridLimeReader { ////////////////////////////////////////////////////////////////// // this helper function wraps the logic for choosing // the correct munger and intermediate lattice data type - // before then reading a lattice field from disk. - // This is a runtime choice as Grid won't know which munger/itype - // to use until it has read the header of a lattice cfg. Therefore - // templating this function on fp_fmt etc. does not work. - // At present we use a rather cumbersome if statement with 6 branches. - // The benefit of organising IldgReader this way is it makes + // before reading a lattice field from disk. + // This is a runtime choice as Grid has to first read the + // header of the lattice cfg. Therefore templating this function on + // FloatingPointFormat etc. does not work. It is rather + // cumbersome, it's if statement with 6 branches, + // but the benefit of organising IldgReader this way is it makes // readConfiguration clearer and more readable. ////////////////////////////////////////////////////////////////// template void readLatticeBinaryObject(Lattice &Umu, std::string filename, FloatingPointFormat fp_fmt, MatrixFormat matrix_fmt, bool is_grp_su, bool is_grp_sp, uint64_t &offset, uint32_t &nersc_csum, uint32_t &scidac_csuma, uint32_t &scidac_csumb) { - typedef typename vobj::scalar_object sobj; // need all the types we could possibly read from // including the intermediate data types for // reduced format lattices and single/double precision + typedef typename vobj::scalar_object sobj; typedef LorentzColourMatrixF fobj; typedef LorentzColourMatrixD dobj; typedef LorentzColour2x3F fobjsuR; @@ -825,15 +825,15 @@ class IldgReader : public GridLimeReader { uint32_t scidac_csuma; uint32_t scidac_csumb; - // these variables store information about the binary data that is read - // from the ildg-format header. if matrix_fmt==REDUCED + // these variables store information about the lattice that is read + // from its ildg-format header. if matrix_fmt==MatrixFormat::REDUCED // then Grid will reconstruct the full matrix using the appropriate munger. - FloatingPointFormat fp_fmt; MatrixFormat matrix_fmt; bool is_grp_su; bool is_grp_sp; // Binary format std::string format; + FloatingPointFormat fp_fmt; ////////////////////////////////////////////////////////////////////////// // Loop over all records diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index cad9f086df..5455dc0818 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -397,7 +397,9 @@ struct GaugeSpunmunger{ } }; -// use these as non-type template parameters in the GaugeUnMunger +// these are used as non-type template parameters when +// writing with GaugeUnMunger and as regular parameters +// when reading with IldgReader.readConfiguration enum struct FloatingPointFormat { IEEE64BIG, IEEE32BIG }; enum struct MatrixFormat { FULL, REDUCED }; ///////////////////////////////////////////////////////// @@ -407,6 +409,8 @@ enum struct MatrixFormat { FULL, REDUCED }; // > no group reduction - treat SU and Sp the same // > group reduction for SU fields // > group reduction for Sp fields +// It also exposes the intermediate out_type for use in +// IldgWriter.writeConfiguration ///////////////////////////////////////////////////////// template struct GaugeUnMunger; @@ -424,7 +428,7 @@ struct GaugeUnMunger } }; -//group reduction for SU +//template specialisation for SU template struct GaugeUnMunger { @@ -442,7 +446,7 @@ struct GaugeUnMunger } }; -//group reduction for Sp +//template specialisation for Sp template struct GaugeUnMunger { From cfedb495afdbafc1d37e6ce60d30d56030c50900 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 20 Jan 2026 22:17:34 +0000 Subject: [PATCH 14/67] small change to docs --- documentation/manual.rst | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/documentation/manual.rst b/documentation/manual.rst index 03b1e933e3..06251d752e 100644 --- a/documentation/manual.rst +++ b/documentation/manual.rst @@ -2104,7 +2104,7 @@ These are specialised to SciDAC writers, introducing facilities for generating t void readScidacFieldRecord(Lattice &field,userRecord &_userRecord); }; -They are also specialised to ILDG format writers, available and defined only for Gauge configurations. Following the publication of Rev. 1.2 of the ILDG Binary File Format by the ILDG Metadata Working Group the ILDG format writer and reader have been updated to be able to write SU(2/3) and Sp(2N) lattice configurations in reduced formats and to read and reconstruct reduced format lattices.:: +They are also specialised to ILDG format writers, available and defined only for Gauge configurations. Following the publication of Rev. 1.2 of the ILDG Binary File Format by the ILDG Metadata Working Group the ILDG format writer and reader have been updated to be able to read and write SU(2/3) and Sp(2N) lattice configurations in reduced formats.:: class IldgWriter : public ScidacWriter { @@ -2120,9 +2120,7 @@ They are also specialised to ILDG format writers, available and defined only for }; -The template parameters *group_name*, *MatrixFormat*, *FloatingPointFormat* specify the exact format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced or not reduced, and two options for the floating point format, 32 bit and 64 bit precision. If a lattice is not written out in reduced format then *group_name* has no practical effect on the writing process. - -Put a table here describing these constructs in more detail? +The template parameters *group_name*, *MatrixFormat*, *FloatingPointFormat* specify the format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced or not reduced, and two options for the floating point format, 32 bit and 64 bit precision. Reduced format Sp gauge fields take up half as much memory as the non-reduced fields while SU(3) fields take up 2/3 as much memory as their non-reduced counterparts. **Implementation detail** From f3221253b71eb715ac754ff0cff7fc0c3cab4590 Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 21 Jan 2026 18:17:10 +0000 Subject: [PATCH 15/67] added if-guards and another explanatory comment to the new ildg reduced format Test --- tests/IO/Test_ildg_reducedfmt_io.cc | 52 +++++++++++++++-------------- 1 file changed, 27 insertions(+), 25 deletions(-) diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index f873fc6c23..f6d556f947 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -33,13 +33,16 @@ Author: Gaurav Ray using namespace std; using namespace Grid; +// this test demonstrates and checks IldgWriter/Readers' ability to +// read/write ildg 1.2 compliant SU and Sp lattices, including in reduced format. + ////////////////////////////////////////////// // this template function returns a // LatticeGaugeField of the chosen gauge // group passed as a template argument. ////////////////////////////////////////////// template -LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Grid, std::vector seed = {83,8,34,59} ) { +LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Grid, std::vector seed ) { GridParallelRNG pRNGa(&Grid); std::cout < void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, std::string file) { - // if the IldgReader class is not able to handle - // the binary format do nothing and exit the function if constexpr( std::is_same_v && N%2==1) { std::cout <= 4" << std::endl; std::cout < && matrix_fmt==MatrixFormat::REDUCED && N>3) - { - std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); - LatticeGaugeField Umu_saved = generateHotFieldConfiguration(Grid); + if constexpr( std::is_same_v && matrix_fmt==MatrixFormat::REDUCED && N>3) + { + std::cout <(Grid, seed0); LatticeGaugeField UmuSp = generateHotFieldConfiguration(Grid, seed1); - // write/read SU fields - writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_full64.4000"); - writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_rdc64.4000"); - writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_full32.4000"); - writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su_rdc32.4000"); + // write and read SU fields + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su64_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su64_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su32_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su32_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); - // write/read Sp fields - writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_full64.4000"); - writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_rdc64.4000"); - writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_full32.4000"); - writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp_rdc32.4000"); + // write and read Sp fields + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp64_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp64_"+std::to_string(Nc/2)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp32_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp32_"+std::to_string(Nc/2)+"x"+std::to_string(Nc)+".4000"); Grid_finalize(); #endif From 245849fd2ef8abe1e8b8b3abe674ef6a9dff87d3 Mon Sep 17 00:00:00 2001 From: gray95 Date: Thu, 29 Jan 2026 18:44:33 +0000 Subject: [PATCH 16/67] started adding unit tests for Ildg IO classes. corrected pre-amble for end-to-end test. --- tests/IO/Test_IldgReaderWriter.cc | 76 +++++++++++++++++++++++++++++ tests/IO/Test_ildg_reducedfmt_io.cc | 3 +- 2 files changed, 77 insertions(+), 2 deletions(-) create mode 100644 tests/IO/Test_IldgReaderWriter.cc diff --git a/tests/IO/Test_IldgReaderWriter.cc b/tests/IO/Test_IldgReaderWriter.cc new file mode 100644 index 0000000000..fd44c1cdb8 --- /dev/null +++ b/tests/IO/Test_IldgReaderWriter.cc @@ -0,0 +1,76 @@ + /************************************************************************************* + + Grid physics library, www.github.com/paboyle/Grid + + Source file: ./Test_IldgWriter.cc + + Copyright (C) 2015 + +Author: Azusa Yamaguchi +Author: Peter Boyle +Author: Gaurav Ray + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License along + with this program; if not, write to the Free Software Foundation, Inc., + 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + + See the full license in the file "LICENSE" in the top level distribution directory + *************************************************************************************/ + /* END LEGAL */ +#include + +using namespace std; +using namespace Grid; + +// Unit tests to check the IldgWriter class. + +void checkWriteLimeIldgLFN(std::string &test_string); +void checkWriteConfiguration(LatticeGaugeField &Umu); + +// write a lime record and then read it back and check it matches. +void checkWriteLimeIldgLFN(std::string &test_string) { + auto simd_layout = GridDefaultSimd(4,vComplex::Nsimd()); + auto mpi_layout = GridDefaultMpi(); + Coordinate latt_size ({8,8,8,16}); + + GridCartesian Grid(latt_size,simd_layout,mpi_layout); + + std::string test_file = "test_file.lime"; + std::string check_string; + // open lime file + IldgWriter _ildgWriter(Grid.IsBoss()); + _ildgWriter.open(test_file); + _ildgWriter.writeLimeIldgLFN( test_string ); + _ildgWriter.close(); + // read back record + IldgReader _ildgReader; + _ildgReader.open(test_file); + _ildgReader.readLimeObject(check_string, ILDG_DATA_LFN); + std::cout < Author: Peter Boyle -Author: paboyle Author: Gaurav Ray This program is free software; you can redistribute it and/or modify From 3570f08cf2b31180d1794709683bfb5cf109c75f Mon Sep 17 00:00:00 2001 From: gray95 Date: Mon, 2 Feb 2026 20:38:54 +0000 Subject: [PATCH 17/67] added 1 unit test for Binary munger --- ...IldgReaderWriter.cc => Test_IldgWriter.cc} | 0 tests/IO/Test_io_mungers.cc | 95 +++++++++++++++++++ 2 files changed, 95 insertions(+) rename tests/IO/{Test_IldgReaderWriter.cc => Test_IldgWriter.cc} (100%) create mode 100644 tests/IO/Test_io_mungers.cc diff --git a/tests/IO/Test_IldgReaderWriter.cc b/tests/IO/Test_IldgWriter.cc similarity index 100% rename from tests/IO/Test_IldgReaderWriter.cc rename to tests/IO/Test_IldgWriter.cc diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc new file mode 100644 index 0000000000..2618c9c390 --- /dev/null +++ b/tests/IO/Test_io_mungers.cc @@ -0,0 +1,95 @@ + /************************************************************************************* + + Grid physics library, www.github.com/paboyle/Grid + + Source file: ./Test_IldgWriter.cc + + Copyright (C) 2015 + +Author: Azusa Yamaguchi +Author: Peter Boyle +Author: Gaurav Ray + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License along + with this program; if not, write to the Free Software Foundation, Inc., + 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + + See the full license in the file "LICENSE" in the top level distribution directory + *************************************************************************************/ + /* END LEGAL */ +#include + +using namespace std; +using namespace Grid; + +// Unit tests to check the various (un)munger classes in parallelIO/Metadata.h +void checkBinarySimpleMungers(); +void checkGaugeSimpleMungers(); +void checkGaugeDoubleStoredMungers(); +void checkGauge3x2mungers(); +void checkGaugeSpmungers(); + +// mungeing between the same type should yield the same result. +void checkBinarySimpleMungers() { + auto simd_layout = GridDefaultSimd(4,vComplex::Nsimd()); + auto mpi_layout = GridDefaultMpi(); + Coordinate latt_size = GridDefaultLatt(); + + GridCartesian Grid(latt_size,simd_layout,mpi_layout); + + LorentzColourMatrixF in_scalar_objectF; // single precision + LorentzColourMatrixF out_scalar_objectF; // single precision + LorentzColourMatrixD in_scalar_objectD; // double precision + LorentzColourMatrixD out_scalar_objectD; // double precision + + out_scalar_objectF = Zero(); + out_scalar_objectD = Zero(); + + in_scalar_objectF = 1.0; + in_scalar_objectD = 1.0; + + // BinarySimpleUnmunger + BinarySimpleUnmunger unmungerFF; + BinarySimpleUnmunger unmungerDF; + BinarySimpleUnmunger unmungerDD; + BinarySimpleUnmunger unmungerFD; + + unmungerFF(in_scalar_objectF, out_scalar_objectF); + assert(in_scalar_objectF==out_scalar_objectF); + assert(in_scalar_objectD==out_scalar_objectD); + + // BinarySimpleMunger + BinarySimpleMunger mungerFF; + BinarySimpleMunger mungerDF; + BinarySimpleMunger mungerDD; + BinarySimpleMunger mungerFD; + + +} + +int main (int argc, char ** argv) +{ +#ifdef HAVE_LIME + Grid_init(&argc,&argv); + + std::cout < Date: Tue, 3 Feb 2026 20:52:13 +0000 Subject: [PATCH 18/67] make comment on Spmunger clearer --- Grid/parallelIO/MetaData.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 5455dc0818..54f67d015b 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -383,7 +383,7 @@ struct GaugeSpmunger{ } }; -// write out Sp(2N) fields in reduced Nx2N format. +// transform Sp(2N) fields into reduced Nx2N format. template struct GaugeSpunmunger{ void operator() (sobj &in,fobj &out){ From 8c2cae5332d4eaa9fe14fdd9d86390b3d73f334e Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 3 Feb 2026 20:52:47 +0000 Subject: [PATCH 19/67] added unit test for Sp(un)munger --- tests/IO/Test_io_mungers.cc | 63 +++++++++++++++++++++++++++++++++++-- 1 file changed, 61 insertions(+), 2 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 2618c9c390..9c449de3b2 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -66,17 +66,76 @@ void checkBinarySimpleMungers() { unmungerFF(in_scalar_objectF, out_scalar_objectF); assert(in_scalar_objectF==out_scalar_objectF); - assert(in_scalar_objectD==out_scalar_objectD); + //assert(in_scalar_objectD==out_scalar_objectD); // BinarySimpleMunger BinarySimpleMunger mungerFF; BinarySimpleMunger mungerDF; BinarySimpleMunger mungerDD; BinarySimpleMunger mungerFD; + +} + +void checkGaugeSpmungers() { + ColourMatrix Sp_field; + + const Complex a(0.5, 0.5), abar(0.5, -0.5); + const Complex b(0.3, 0.9), bbar(0.3, -0.9); + + // fill top left + for(int i=0;i + GaugeSpunmunger unmunger; + LorentzColourNx2ND test_rdc = Zero(); + unmunger(test,test_rdc); + + //std::cout << test_rdc << std::endl; + + // reconstruct full field. GaugeSpmunger + GaugeSpmunger munger; + LorentzColourMatrixD test_recon; + munger(test_rdc,test_recon); + + std::cout << test_recon << std::endl; + + // round-trip test + assert(test==test_recon); } + int main (int argc, char ** argv) { #ifdef HAVE_LIME @@ -88,7 +147,7 @@ int main (int argc, char ** argv) //checkGaugeSimpleMungers(); //checkGaugeDoubleStoredMungers(); //checkGauge3x2mungers(); - //checkGaugeSpmungers(); + checkGaugeSpmungers(); Grid_finalize(); #endif From 623f4eef3ac9fc81762c6af46dc71eceb29591bc Mon Sep 17 00:00:00 2001 From: gray95 Date: Thu, 5 Feb 2026 19:48:39 +0000 Subject: [PATCH 20/67] added unit test for 3x2 mungers --- tests/IO/Test_io_mungers.cc | 57 +++++++++++++++++++++++++++++++++---- 1 file changed, 52 insertions(+), 5 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 9c449de3b2..d52119cd7d 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -32,7 +32,11 @@ Author: Gaurav Ray using namespace std; using namespace Grid; -// Unit tests to check the various (un)munger classes in parallelIO/Metadata.h +// Unit tests to check the various (un)munger classes +// and functions in parallelIO/Metadata.h +void check_reconstruct3(); +void check_reconstructSp(); +// void check reconstructSU(); void checkBinarySimpleMungers(); void checkGaugeSimpleMungers(); void checkGaugeDoubleStoredMungers(); @@ -76,6 +80,50 @@ void checkBinarySimpleMungers() { } +// generate SU(3) matrix -> reduce it -> reconstruct it -> all good? +void checkGauge3x2mungers() { + ColourMatrix Ta, Tb, Tc, arg; + SU3::generator(1, Ta); + SU3::generator(4, Tb); + SU3::generator(6, Tc); + + const RealD a(0.6), b(0.9), c(0.7); + + ColourMatrix SU3_field; + arg = timesI(a*Ta + b*Tb + c*Tc); + SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? + +// std::cout << SU3_field*adj(SU3_field); + LorentzColourMatrixD test=Zero(); + + for(int mu=0; mu + Gauge3x2unmunger unmunger; + LorentzColour2x3D test_rdc = Zero(); + unmunger(test,test_rdc); + + //std::cout << test_rdc; + + // reconstruct full field. Gauge3x2munger + Gauge3x2munger munger; + LorentzColourMatrixD test_recon; + munger(test_rdc,test_recon); + + std::cout << test_recon << std::endl; + + // round-trip test + std::cout << norm2(test_recon-test); + //assert(norm2(test_recon-test)<1e-05); + +} + +// round trip test - start with field in non-reduced format then reduce it and +// recover the non-reduced field before assert(final==start) void checkGaugeSpmungers() { ColourMatrix Sp_field; @@ -132,7 +180,6 @@ void checkGaugeSpmungers() { // round-trip test assert(test==test_recon); - } @@ -143,11 +190,11 @@ int main (int argc, char ** argv) std::cout < Date: Mon, 9 Feb 2026 19:48:26 +0000 Subject: [PATCH 21/67] adds unit test for reconstruct3(LorentzColourMatrix &cm) --- tests/IO/Test_io_mungers.cc | 186 +++++++++++++++++++++++++++++++++--- 1 file changed, 173 insertions(+), 13 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index d52119cd7d..0b0ab2f067 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -43,22 +43,111 @@ void checkGaugeDoubleStoredMungers(); void checkGauge3x2mungers(); void checkGaugeSpmungers(); +// check result is still in SU, determinant==1 +// need to mock up an SU matrix +void check_reconstruct3() { + + ColourMatrixD Ta, Tb, Tc, arg; + SU3::generator(1, Ta); + SU3::generator(4, Tb); + SU3::generator(6, Tc); + + const RealD a(0.6), b(0.9), c(0.7); + + ColourMatrixD SU3_field; + arg = timesI(a*Ta + b*Tb + c*Tc); + SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? + std::cout << GridLogMessage << SU3_field()()(Nc-1,1) << std::endl; + + //set last row equal to zero + for(int j=0;j(test, mu); + std::cout << GridLogMessage << new_cm()()(Nc-1,1) << std::endl; + auto det = Determinant(new_cm); + std::cout << GridLogMessage << "det: " << det << " norm2: " << norm2(det) << std::endl; + std::cout << GridLogMessage << "norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; + assert( abs(norm2(det)-1.0) < 1e-15 ); + } + + +} + +// is result an Sp matrix, Zero(), One()? +/*void check_reconstructSp() { + + ColourMatrix Sp_field = Zero(); + + LorentzColourMatrix cm = Zero(); + LorentzColourMatrix cm_zero = Zero(); + LorentzColourMatrix cm_one; + cm_one = 1.0; + + // should be zero matrix after reconstruction. + reconstructSp(cm); + assert(cm==cm_zero); + + + reconstructSp(cm); + + + const Complex a(0.5, 0.5), abar(0.5, -0.5); + const Complex b(0.3, 0.9), bbar(0.3, -0.9); + + // fill top left + for(int i=0;i mungerDF; BinarySimpleMunger mungerDD; BinarySimpleMunger mungerFD; + + mungerDD(in_scalar_objectD, out_scalar_objectD); + assert(in_scalar_objectD==out_scalar_objectD); + //assert(out_scalar_objectF==out_scalar_objectD); + + +} +// these are used in NerscIO.h +void checkGaugeSimpleMungers() { + + LorentzColourMatrixF in_scalar_objectF; // single precision + LorentzColourMatrixF out_scalar_objectF = Zero(); + LorentzColourMatrixD in_scalar_objectD; // double precision + LorentzColourMatrixD out_scalar_objectD = Zero(); + in_scalar_objectF = 1.0; + in_scalar_objectD = 1.0; + + // GaugeSimpleUnmunger + GaugeSimpleUnmunger unmungerFF; + GaugeSimpleUnmunger unmungerDF; + GaugeSimpleUnmunger unmungerDD; + GaugeSimpleUnmunger unmungerFD; + + unmungerFF(in_scalar_objectF, out_scalar_objectF); + assert(in_scalar_objectF==out_scalar_objectF); + + // GaugeSimpleMunger + GaugeSimpleMunger mungerFF; + GaugeSimpleMunger mungerDF; + GaugeSimpleMunger mungerDD; + GaugeSimpleMunger mungerFD; + + mungerDD(in_scalar_objectD, out_scalar_objectD); + assert(in_scalar_objectD==out_scalar_objectD); + +} + +void checkGaugeDoubleStoredMungers() { + + DoubleStoredColourMatrixF in_scalar_objectF; // single precision + DoubleStoredColourMatrixF out_scalar_objectF = Zero(); + DoubleStoredColourMatrixD in_scalar_objectD; // double precision + DoubleStoredColourMatrixD out_scalar_objectD = Zero(); + + in_scalar_objectF = 1.0; + in_scalar_objectD = 1.0; + + // GaugeDoubleStoredUnmunger + GaugeDoubleStoredUnmunger unmungerFF; + GaugeDoubleStoredUnmunger unmungerDF; + GaugeDoubleStoredUnmunger unmungerDD; + GaugeDoubleStoredUnmunger unmungerFD; + + unmungerFF(in_scalar_objectF, out_scalar_objectF); + assert(in_scalar_objectF==out_scalar_objectF); + + // GaugeDoubleStoredMunger + GaugeDoubleStoredMunger mungerFF; + GaugeDoubleStoredMunger mungerDF; + GaugeDoubleStoredMunger mungerDD; + GaugeDoubleStoredMunger mungerFD; + + mungerDD(in_scalar_objectD, out_scalar_objectD); + assert(in_scalar_objectD==out_scalar_objectD); + } + // generate SU(3) matrix -> reduce it -> reconstruct it -> all good? void checkGauge3x2mungers() { - ColourMatrix Ta, Tb, Tc, arg; + ColourMatrixD Ta, Tb, Tc, arg; SU3::generator(1, Ta); SU3::generator(4, Tb); SU3::generator(6, Tc); const RealD a(0.6), b(0.9), c(0.7); - ColourMatrix SU3_field; + ColourMatrixD SU3_field; arg = timesI(a*Ta + b*Tb + c*Tc); SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? @@ -190,10 +345,15 @@ int main (int argc, char ** argv) std::cout < Date: Wed, 11 Feb 2026 21:52:35 +0000 Subject: [PATCH 22/67] adds unit test for reconstructSp and adds if-constexpr guards for SU/Sp differentiation --- tests/IO/Test_io_mungers.cc | 110 +++++++++++++++++++----------------- 1 file changed, 57 insertions(+), 53 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 0b0ab2f067..043f0e53d7 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -34,9 +34,8 @@ using namespace Grid; // Unit tests to check the various (un)munger classes // and functions in parallelIO/Metadata.h -void check_reconstruct3(); +void check_reconstructSU(); void check_reconstructSp(); -// void check reconstructSU(); void checkBinarySimpleMungers(); void checkGaugeSimpleMungers(); void checkGaugeDoubleStoredMungers(); @@ -45,26 +44,26 @@ void checkGaugeSpmungers(); // check result is still in SU, determinant==1 // need to mock up an SU matrix -void check_reconstruct3() { +void check_reconstructSU() { ColourMatrixD Ta, Tb, Tc, arg; - SU3::generator(1, Ta); - SU3::generator(4, Tb); - SU3::generator(6, Tc); + SU::generator(1, Ta); + SU::generator(4, Tb); + SU::generator(6, Tc); const RealD a(0.6), b(0.9), c(0.7); ColourMatrixD SU3_field; arg = timesI(a*Ta + b*Tb + c*Tc); SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? - std::cout << GridLogMessage << SU3_field()()(Nc-1,1) << std::endl; + //std::cout << GridLogMessage << SU3_field()()(Nc-1,1) << std::endl; //set last row equal to zero for(int j=0;j(test, mu); - std::cout << GridLogMessage << new_cm()()(Nc-1,1) << std::endl; + //std::cout << GridLogMessage << new_cm()()(Nc-1,1) << std::endl; auto det = Determinant(new_cm); - std::cout << GridLogMessage << "det: " << det << " norm2: " << norm2(det) << std::endl; - std::cout << GridLogMessage << "norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; + //std::cout << GridLogMessage << "det: " << det << " norm2: " << norm2(det) << std::endl; + //std::cout << GridLogMessage << "norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; assert( abs(norm2(det)-1.0) < 1e-15 ); } - } -// is result an Sp matrix, Zero(), One()? -/*void check_reconstructSp() { - - ColourMatrix Sp_field = Zero(); +void check_reconstructSp() { LorentzColourMatrix cm = Zero(); LorentzColourMatrix cm_zero = Zero(); - LorentzColourMatrix cm_one; - cm_one = 1.0; // should be zero matrix after reconstruction. reconstructSp(cm); assert(cm==cm_zero); - - - reconstructSp(cm); + ColourMatrix Sp_field = Zero(); const Complex a(0.5, 0.5), abar(0.5, -0.5); const Complex b(0.3, 0.9), bbar(0.3, -0.9); @@ -121,24 +112,30 @@ void check_reconstruct3() { Sp_field()()(i,j) = b; } } - // fill bottom left - for(int i=2;i(test,mu); + for(int i=0;i reduce it -> reconstruct it -> all good? void checkGauge3x2mungers() { ColourMatrixD Ta, Tb, Tc, arg; - SU3::generator(1, Ta); - SU3::generator(4, Tb); - SU3::generator(6, Tc); + SU::generator(1, Ta); + SU::generator(4, Tb); + SU::generator(6, Tc); const RealD a(0.6), b(0.9), c(0.7); @@ -255,7 +252,7 @@ void checkGauge3x2mungers() { pokeLorentz(test,SU3_field,mu); } - std::cout << test << std::endl; + //std::cout << test << std::endl; // reduce field. Gauge3x2unmunger Gauge3x2unmunger unmunger; @@ -269,11 +266,10 @@ void checkGauge3x2mungers() { LorentzColourMatrixD test_recon; munger(test_rdc,test_recon); - std::cout << test_recon << std::endl; + //std::cout << test_recon << std::endl; // round-trip test - std::cout << norm2(test_recon-test); - //assert(norm2(test_recon-test)<1e-05); + assert(norm2(test_recon-test)<10e-7); } @@ -316,22 +312,18 @@ void checkGaugeSpmungers() { pokeLorentz(test,Sp_field,mu); } - std::cout << test << std::endl; + //std::cout << test << std::endl; // reduce field. GaugeSpunmunger GaugeSpunmunger unmunger; LorentzColourNx2ND test_rdc = Zero(); unmunger(test,test_rdc); - //std::cout << test_rdc << std::endl; - // reconstruct full field. GaugeSpmunger GaugeSpmunger munger; LorentzColourMatrixD test_recon; munger(test_rdc,test_recon); - std::cout << test_recon << std::endl; - // round-trip test assert(test==test_recon); @@ -348,13 +340,25 @@ int main (int argc, char ** argv) std::cout << GridLogMessage << "Nd is " << Nd << std::endl; std::cout << GridLogMessage << "Nc is " << Nc << std::endl; std::cout << GridLogMessage << "Nds is " << Nds << std::endl; - - check_reconstruct3(); - //checkBinarySimpleMungers(); - //checkGaugeSimpleMungers(); - //checkGaugeDoubleStoredMungers(); - //checkGauge3x2mungers(); - //checkGaugeSpmungers(); + + if constexpr(Nc>1 && Nc<4) { + std::cout << GridLogMessage << "Checking SU mungers " << std::endl; + check_reconstructSU(); + checkGauge3x2mungers(); + } + + if constexpr(Nc>2 && Nc%2==0) { + std::cout << GridLogMessage << "Checking Sp mungers " << std::endl; + check_reconstructSp(); + checkGaugeSpmungers(); + } + + std::cout << GridLogMessage << "Checking BinarySimple mungers " << std::endl; + checkBinarySimpleMungers(); + std::cout << GridLogMessage << "Checking GaugeSimple mungers " << std::endl; + checkGaugeSimpleMungers(); + std::cout << GridLogMessage << "Checking GaugeDoubleStored mungers " << std::endl; + checkGaugeDoubleStoredMungers(); Grid_finalize(); #endif From 6f112f9f9243eddd5b68e1b4e8a04df3295eb034 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 9 Dec 2025 19:20:07 +0000 Subject: [PATCH 23/67] resolve merge conflict from SU(N) branch --- Grid/parallelIO/MetaData.h | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 54f67d015b..782ca08f0d 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -217,6 +217,26 @@ inline void reconstruct3(LorentzColourMatrix & cm) } } +// insert multiline comment here +// template on vtype? +inline void reconstructSU(LorentzColourMatrix &cm) +{ + // what type does this expect? + using iLorentzColourNminus1xNminus1 = iVector, Nc-1>, Nd >; + using iLorentzColourNminus1xNminus1 = LorentzColourNmxNmD; + LorentzColourNmxNmD tmp(grid); + for(int mu=0;muk) ? j-1 : j; // for correct indexing of columns in tmp + tmp(mu)()(i,J) = cm(mu)()(i,j); + } + } + cm(mu)()(N,k) = adj(Determinant(tmp)); + } + } +} //////////////////////////////////////////////////////////////// // this function takes a reduced format // Sp(2N) field with N rows and 2N columns From facc2b38b4ecbc92fd2c56ffc387ff074fa2e295 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 17 Feb 2026 20:00:54 +0000 Subject: [PATCH 24/67] fixed reconstructSU and tested with N=3/4 --- Grid/parallelIO/MetaData.h | 33 +++++++++++++++++++++------------ 1 file changed, 21 insertions(+), 12 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 782ca08f0d..561d2ddadc 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -218,25 +218,34 @@ inline void reconstruct3(LorentzColourMatrix & cm) } // insert multiline comment here -// template on vtype? inline void reconstructSU(LorentzColourMatrix &cm) { - // what type does this expect? - using iLorentzColourNminus1xNminus1 = iVector, Nc-1>, Nd >; - using iLorentzColourNminus1xNminus1 = LorentzColourNmxNmD; - LorentzColourNmxNmD tmp(grid); - for(int mu=0;muk) ? j-1 : j; // for correct indexing of columns in tmp - tmp(mu)()(i,J) = cm(mu)()(i,j); + using ColourMatrixNm = iScalar > > ; + + ColourMatrixNm red; // for the Nc-1 by Nc-1 matrix + ColourMatrix old; // for the Nc-1 by N matrix read from disk + + for(int mu=0; mu(cm,mu); + for(int k=0; k(cm,old,mu); } } + //////////////////////////////////////////////////////////////// // this function takes a reduced format // Sp(2N) field with N rows and 2N columns From b8fdecea2130b44d280b983178c835bf2927f658 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 17 Feb 2026 21:49:06 +0000 Subject: [PATCH 25/67] implements GaugeSU(un)mungers --- Grid/parallelIO/MetaData.h | 37 ++++++++++++++++++++++++++++++++++++- 1 file changed, 36 insertions(+), 1 deletion(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 561d2ddadc..2ae2e9ba45 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -217,7 +217,14 @@ inline void reconstruct3(LorentzColourMatrix & cm) } } -// insert multiline comment here +// the elements of the final row of an SU(N) matrix +// can be obtained from the determinants of the N N-1 by N-1 matrices +// formed by deleting each column and the Nth row of the SU(N) matrix. +// 1) create a ColourMatrixNm object to store the N-1 dim matrix. +// 2) peek the colour matrix in a particular direction. +// 3) fill the ColourMatrixNm object from the peeked matrix. +// 4) take the Determinant and fill the corresponding element +// 5) repeat for each Lorentz index inline void reconstructSU(LorentzColourMatrix &cm) { using ColourMatrixNm = iScalar > > ; @@ -398,6 +405,33 @@ struct Gauge3x2unmunger{ } }; +template +struct GaugeSUmunger{ + void operator() (fobj &in,sobj &out){ + for(int mu=0;mu +struct GaugeSUunmunger{ + void operator() (sobj &in,fobj &out){ + for(int mu=0;mu struct GaugeSpmunger{ void operator() (fobj &in,sobj &out){ @@ -467,6 +501,7 @@ struct GaugeUnMunger BinarySimpleUnmunger binary_unmunger; Gauge3x2unmunger gauge_unmunger; + //GaugeSUunmunger gauge_unmunger; void operator() (in_type &in, out_type &out){ tmp_type tmp; From 48cc7bdb9f44e9fef2145c7d513cc8186fb6e7a3 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 17 Feb 2026 21:49:45 +0000 Subject: [PATCH 26/67] testing new reconstructSU function and GaugeSUmungers --- tests/IO/Test_io_mungers.cc | 115 +++++++++++++++++++++++++++++++----- 1 file changed, 101 insertions(+), 14 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 043f0e53d7..be5020ce20 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -32,26 +32,28 @@ Author: Gaurav Ray using namespace std; using namespace Grid; -// Unit tests to check the various (un)munger classes +// tests to check the various (un)munger classes // and functions in parallelIO/Metadata.h +void check_reconstruct3(); void check_reconstructSU(); void check_reconstructSp(); void checkBinarySimpleMungers(); void checkGaugeSimpleMungers(); void checkGaugeDoubleStoredMungers(); void checkGauge3x2mungers(); +void checkGaugeSUmungers(); void checkGaugeSpmungers(); // check result is still in SU, determinant==1 // need to mock up an SU matrix -void check_reconstructSU() { +void check_reconstruct3() { ColourMatrixD Ta, Tb, Tc, arg; SU::generator(1, Ta); SU::generator(4, Tb); SU::generator(6, Tc); - const RealD a(0.6), b(0.9), c(0.7); + const RealD a(0.1), b(0.2), c(0.3); ColourMatrixD SU3_field; arg = timesI(a*Ta + b*Tb + c*Tc); @@ -72,20 +74,56 @@ void check_reconstructSU() { reconstruct3(test); - ColourMatrixD new_cm = Zero(); - // check result is unitary and det==1 for(int mu=0;mu(test, mu); + auto new_cm = peekIndex(test, mu); //std::cout << GridLogMessage << new_cm()()(Nc-1,1) << std::endl; auto det = Determinant(new_cm); - //std::cout << GridLogMessage << "det: " << det << " norm2: " << norm2(det) << std::endl; - //std::cout << GridLogMessage << "norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; - assert( abs(norm2(det)-1.0) < 1e-15 ); + std::cout << GridLogMessage << "reconstruct3::norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; + //assert( abs(norm2(det)-1.0) < 1e-15 ); } } +void check_reconstructSU() { + + ColourMatrixD Ta, Tb, Tc, arg, SUN_field; + SU::generator(2, Ta); + SU::generator(4, Tb); + SU::generator(7, Tc); + + const RealD a(0.1), b(0.2), c(0.3); + + arg = timesI(a*Ta + b*Tb + c*Tc); + SUN_field = Exponentiate(arg, 2.0); // what is RealD alpha?? + //std::cout << GridLogMessage << "Original Field\n" << SUN_field << std::endl; + + //set last row equal to zero + for(int j=0;j(test, mu); + auto det = Determinant(new_cm); + //std::cout << GridLogMessage << "reconstructSU::norm2(det): " << norm2(det) << std::endl; + //std::cout << GridLogMessage << "reconstructSU::norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; + assert( abs( norm2(det)-1.0 ) < 1e-15 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); + } +} + + void check_reconstructSp() { LorentzColourMatrix cm = Zero(); @@ -269,7 +307,53 @@ void checkGauge3x2mungers() { //std::cout << test_recon << std::endl; // round-trip test - assert(norm2(test_recon-test)<10e-7); + std::cout << "SU3mungeing: " << norm2(test_recon-test) << std::endl; + //assert(norm2(test_recon-test)<10e-7); + +} + +void checkGaugeSUmungers() { + ColourMatrixD Ta, Tb, Tc, arg; + SU::generator(1, Ta); + SU::generator(4, Tb); + SU::generator(6, Tc); + + const RealD a(0.6), b(0.9), c(0.7); + + ColourMatrixD SUN_field; + arg = timesI(a*Ta + b*Tb + c*Tc); + SUN_field = Exponentiate(arg, 2.0); // what is RealD alpha?? + +// std::cout << SUN_field*adj(SUN_field); + LorentzColourMatrixD test=Zero(); + + for(int mu=0; mu + GaugeSUunmunger unmunger; + LorentzColour2x3D test_rdc = Zero(); + unmunger(test,test_rdc); + + //std::cout << test_rdc; + + // reconstruct full field. GaugeSUmunger + GaugeSUmunger munger; + LorentzColourMatrixD test_recon; + munger(test_rdc,test_recon); + + // round-trip test + std::cout << "SUmungeing:norm2 " << norm2(test_recon-test) << std::endl; + //assert(norm2(test_recon-test)<10e-7); + // check result is unitary and det==1 + for(int mu=0;mu(test_recon, mu); + auto det = Determinant(new_cm); + std::cout << GridLogMessage << "reconstructSU::norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; + //assert( abs(norm2(det)-1.0) < 1e-15 ); + } + } @@ -312,8 +396,6 @@ void checkGaugeSpmungers() { pokeLorentz(test,Sp_field,mu); } - //std::cout << test << std::endl; - // reduce field. GaugeSpunmunger GaugeSpunmunger unmunger; LorentzColourNx2ND test_rdc = Zero(); @@ -342,8 +424,8 @@ int main (int argc, char ** argv) std::cout << GridLogMessage << "Nds is " << Nds << std::endl; if constexpr(Nc>1 && Nc<4) { - std::cout << GridLogMessage << "Checking SU mungers " << std::endl; - check_reconstructSU(); + std::cout << GridLogMessage << "Checking SU(2/3) mungers " << std::endl; + check_reconstruct3(); checkGauge3x2mungers(); } @@ -360,6 +442,11 @@ int main (int argc, char ** argv) std::cout << GridLogMessage << "Checking GaugeDoubleStored mungers " << std::endl; checkGaugeDoubleStoredMungers(); + std::cout << GridLogMessage << "Checking reconstructSU" << std::endl; + check_reconstructSU(); + std::cout << GridLogMessage << "Checking GaugeSUmungers" << std::endl; + checkGaugeSUmungers(); + Grid_finalize(); #endif } From 070e03909da81c211140153312d6ed8a73dcf996 Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 18 Feb 2026 17:14:17 +0000 Subject: [PATCH 27/67] fixed reconstructSU - now pretty sure works with N=3/4 --- Grid/parallelIO/MetaData.h | 6 +----- 1 file changed, 1 insertion(+), 5 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 2ae2e9ba45..366ee2cbc2 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -241,13 +241,9 @@ inline void reconstructSU(LorentzColourMatrix &cm) for(int j=0; j(cm,old,mu); } From 5c6a57d88cd3eed5cec3166b60712d1aa638bfcb Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 18 Feb 2026 17:14:53 +0000 Subject: [PATCH 28/67] cleaned up unit tests and debug comments --- tests/IO/Test_io_mungers.cc | 53 +++++++++---------------------------- 1 file changed, 13 insertions(+), 40 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index be5020ce20..c7fd9d7eb4 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -2,7 +2,7 @@ Grid physics library, www.github.com/paboyle/Grid - Source file: ./Test_IldgWriter.cc + Source file: ./Test_io_mungers.cc Copyright (C) 2015 @@ -44,10 +44,9 @@ void checkGauge3x2mungers(); void checkGaugeSUmungers(); void checkGaugeSpmungers(); -// check result is still in SU, determinant==1 -// need to mock up an SU matrix void check_reconstruct3() { + // start by mocking up an SU matrix ColourMatrixD Ta, Tb, Tc, arg; SU::generator(1, Ta); SU::generator(4, Tb); @@ -57,15 +56,13 @@ void check_reconstruct3() { ColourMatrixD SU3_field; arg = timesI(a*Ta + b*Tb + c*Tc); - SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? - //std::cout << GridLogMessage << SU3_field()()(Nc-1,1) << std::endl; - + SU3_field = Exponentiate(arg, 2.0); + //set last row equal to zero for(int j=0;j(test, mu); - //std::cout << GridLogMessage << new_cm()()(Nc-1,1) << std::endl; auto det = Determinant(new_cm); - std::cout << GridLogMessage << "reconstruct3::norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; - //assert( abs(norm2(det)-1.0) < 1e-15 ); + assert( abs(norm2(det)-1.0) < 1e-15 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); } } @@ -95,15 +91,13 @@ void check_reconstructSU() { const RealD a(0.1), b(0.2), c(0.3); arg = timesI(a*Ta + b*Tb + c*Tc); - SUN_field = Exponentiate(arg, 2.0); // what is RealD alpha?? - //std::cout << GridLogMessage << "Original Field\n" << SUN_field << std::endl; + SUN_field = Exponentiate(arg, 2.0); //set last row equal to zero for(int j=0;j(test, mu); auto det = Determinant(new_cm); - //std::cout << GridLogMessage << "reconstructSU::norm2(det): " << norm2(det) << std::endl; - //std::cout << GridLogMessage << "reconstructSU::norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; assert( abs( norm2(det)-1.0 ) < 1e-15 ); assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); } @@ -194,7 +186,6 @@ void checkBinarySimpleMungers() { unmungerFF(in_scalar_objectF, out_scalar_objectF); assert(in_scalar_objectF==out_scalar_objectF); - //assert(in_scalar_objectD==out_scalar_objectD); // BinarySimpleMunger BinarySimpleMunger mungerFF; @@ -204,8 +195,6 @@ void checkBinarySimpleMungers() { mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); - //assert(out_scalar_objectF==out_scalar_objectD); - } // these are used in NerscIO.h @@ -269,8 +258,6 @@ void checkGaugeDoubleStoredMungers() { } - -// generate SU(3) matrix -> reduce it -> reconstruct it -> all good? void checkGauge3x2mungers() { ColourMatrixD Ta, Tb, Tc, arg; SU::generator(1, Ta); @@ -283,32 +270,24 @@ void checkGauge3x2mungers() { arg = timesI(a*Ta + b*Tb + c*Tc); SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? -// std::cout << SU3_field*adj(SU3_field); LorentzColourMatrixD test=Zero(); for(int mu=0; mu Gauge3x2unmunger unmunger; LorentzColour2x3D test_rdc = Zero(); unmunger(test,test_rdc); - //std::cout << test_rdc; - // reconstruct full field. Gauge3x2munger Gauge3x2munger munger; LorentzColourMatrixD test_recon; munger(test_rdc,test_recon); - //std::cout << test_recon << std::endl; - // round-trip test - std::cout << "SU3mungeing: " << norm2(test_recon-test) << std::endl; - //assert(norm2(test_recon-test)<10e-7); + assert(norm2(test_recon-test)<1e-10); } @@ -324,7 +303,6 @@ void checkGaugeSUmungers() { arg = timesI(a*Ta + b*Tb + c*Tc); SUN_field = Exponentiate(arg, 2.0); // what is RealD alpha?? -// std::cout << SUN_field*adj(SUN_field); LorentzColourMatrixD test=Zero(); for(int mu=0; mu GaugeSUmunger munger; LorentzColourMatrixD test_recon; munger(test_rdc,test_recon); - + // round-trip test - std::cout << "SUmungeing:norm2 " << norm2(test_recon-test) << std::endl; - //assert(norm2(test_recon-test)<10e-7); + assert(norm2(test_recon-test)<1e-15); + // check result is unitary and det==1 for(int mu=0;mu(test_recon, mu); auto det = Determinant(new_cm); - std::cout << GridLogMessage << "reconstructSU::norm2(UU^dag-1.0): " << norm2( new_cm*adj(new_cm)-1.0 ) << std::endl; - //assert( abs(norm2(det)-1.0) < 1e-15 ); + assert( abs(norm2(det)-1.0) < 1e-15 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); } - } -// round trip test - start with field in non-reduced format then reduce it and -// recover the non-reduced field before assert(final==start) void checkGaugeSpmungers() { ColourMatrix Sp_field; From 1810e7a12b1d7665eb7542dc690388ad7c11d815 Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 18 Feb 2026 17:53:37 +0000 Subject: [PATCH 29/67] swapped out Gauge3x2munger for shiny new GaugeSUmunger --- Grid/parallelIO/IldgIO.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 221398301d..3a13a60775 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -765,7 +765,7 @@ class IldgReader : public GridLimeReader { if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { format = std::string("IEEE64BIG"); if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; + GaugeSUmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; @@ -777,7 +777,7 @@ class IldgReader : public GridLimeReader { } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { format = std::string("IEEE32BIG"); if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - Gauge3x2munger munge; + GaugeSUmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; From c1eb79104bb5b175fcfb554f8c3714b29572a45c Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 18 Feb 2026 17:53:56 +0000 Subject: [PATCH 30/67] swapped out Gauge3x2munger for shiny new GaugeSUmunger --- Grid/parallelIO/MetaData.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 366ee2cbc2..f733ecb8e6 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -496,8 +496,8 @@ struct GaugeUnMunger using out_type = typename std::tuple_element_t(fp_fmt), std::tuple>; BinarySimpleUnmunger binary_unmunger; - Gauge3x2unmunger gauge_unmunger; - //GaugeSUunmunger gauge_unmunger; + //Gauge3x2unmunger gauge_unmunger; + GaugeSUunmunger gauge_unmunger; void operator() (in_type &in, out_type &out){ tmp_type tmp; From 9dfa530de5e0faf72454c3e8338211244e7915c5 Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 18 Feb 2026 17:55:14 +0000 Subject: [PATCH 31/67] removed if-guard so now can check SU(N>3) reduced fmt read/write --- tests/IO/Test_ildg_reducedfmt_io.cc | 9 --------- 1 file changed, 9 deletions(-) diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index 79ae8e315c..971e77bc83 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -84,15 +84,6 @@ void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, st _IldgWriter.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); - if constexpr( std::is_same_v && matrix_fmt==MatrixFormat::REDUCED && N>3) - { - std::cout < Date: Mon, 23 Feb 2026 02:43:54 +0000 Subject: [PATCH 32/67] fixed mock_SU_field - should now produce random SU fields --- tests/IO/Test_io_mungers.cc | 306 ++++++++++++++++++------------------ 1 file changed, 156 insertions(+), 150 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index c7fd9d7eb4..b11ec339c8 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -36,28 +36,52 @@ using namespace Grid; // and functions in parallelIO/Metadata.h void check_reconstruct3(); void check_reconstructSU(); +void checkGauge3x2mungers(); +void checkGaugeSUmungers(); void check_reconstructSp(); +void checkGaugeSpmungers(); void checkBinarySimpleMungers(); void checkGaugeSimpleMungers(); void checkGaugeDoubleStoredMungers(); -void checkGauge3x2mungers(); -void checkGaugeSUmungers(); -void checkGaugeSpmungers(); -void check_reconstruct3() { +void mock_SU_field(); + +////////////////////////////////////////////// +// helper function for generating SU matrices +////////////////////////////////////////////// +void mock_SU_field( ColourMatrixD &cm, std::vector seed ) { + GridSerialRNG sRNG; sRNG.SeedFixedIntegers(seed); - // start by mocking up an SU matrix - ColourMatrixD Ta, Tb, Tc, arg; - SU::generator(1, Ta); - SU::generator(4, Tb); - SU::generator(6, Tc); + ComplexD ca; + ColourMatrixD lie, la, ta; + ComplexD ci(0.0, 1.0); - const RealD a(0.1), b(0.2), c(0.3); + lie = Zero(); + + for (int a = 0; a < Nc*Nc-1; a++) { + random(sRNG, ca); + + ca = (ca + conjugate(ca)) * 0.5; + ca = ca - 0.5; + + SU::generator(a, ta); + + la = timesI(ca * ta); + + lie = lie + la; // e^{i la ta} + } + + cm = Exponentiate(lie, 2.0); + +} +//////////////////////////////////////////// + +void check_reconstruct3() { ColourMatrixD SU3_field; - arg = timesI(a*Ta + b*Tb + c*Tc); - SU3_field = Exponentiate(arg, 2.0); - + std::vector rng_seed = {1,2,3,4}; + mock_SU_field(SU3_field, rng_seed); + //set last row equal to zero for(int j=0;j::generator(2, Ta); - SU::generator(4, Tb); - SU::generator(7, Tc); - - const RealD a(0.1), b(0.2), c(0.3); - - arg = timesI(a*Ta + b*Tb + c*Tc); - SUN_field = Exponentiate(arg, 2.0); + ColourMatrixD SUN_field; + std::vector rnd_seed = {4,8,15,16}; + mock_SU_field(SUN_field, rnd_seed); //set last row equal to zero for(int j=0;j rnd_seed{20,6,15,8}; + mock_SU_field(SU3_field, rnd_seed); + + LorentzColourMatrixD test=Zero(); + for(int mu=0; mu + Gauge3x2unmunger unmunger; + LorentzColour2x3D test_rdc; + unmunger(test,test_rdc); + + // reconstruct full field. Gauge3x2munger + Gauge3x2munger munger; + LorentzColourMatrixD test_recon; + munger(test_rdc,test_recon); + + // round-trip test + //std::cout << "checkGauge3x2mungers: " << norm2(test_recon-test) << std::endl; + assert(norm2(test_recon-test)<1e-10); + +} + +void checkGaugeSUmungers() { + ColourMatrixD SUN_field; + std::vector rnd_seed{23,42,66,98}; + mock_SU_field(SUN_field, rnd_seed); + + LorentzColourMatrixD test=Zero(); + + for(int mu=0; mu + GaugeSUunmunger unmunger; + LorentzColour2x3D test_rdc = Zero(); + unmunger(test,test_rdc); + + // reconstruct full field. GaugeSUmunger + GaugeSUmunger munger; + LorentzColourMatrixD test_recon; + munger(test_rdc,test_recon); + + // round-trip test + //std::cout << "checkGaugeSUmungers: " << norm2(test_recon-test) << std::endl; + assert(norm2(test_recon-test)<1e-15); + + // check result is unitary and det==1 + for(int mu=0;mu(test_recon, mu); + auto det = Determinant(new_cm); + assert( abs(norm2(det)-1.0) < 1e-15 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); + } + +} + void check_reconstructSp() { @@ -166,6 +246,58 @@ void check_reconstructSp() { } +void checkGaugeSpmungers() { + ColourMatrix Sp_field; + + const Complex a(0.5, 0.5), abar(0.5, -0.5); + const Complex b(0.3, 0.9), bbar(0.3, -0.9); + + // fill top left + for(int i=0;i + GaugeSpunmunger unmunger; + LorentzColourNx2ND test_rdc = Zero(); + unmunger(test,test_rdc); + + // reconstruct full field. GaugeSpmunger + GaugeSpmunger munger; + LorentzColourMatrixD test_recon; + munger(test_rdc,test_recon); + + // round-trip test + assert(test==test_recon); + +} + // mungeing between the same type should yield the same result. void checkBinarySimpleMungers() { @@ -258,138 +390,12 @@ void checkGaugeDoubleStoredMungers() { } -void checkGauge3x2mungers() { - ColourMatrixD Ta, Tb, Tc, arg; - SU::generator(1, Ta); - SU::generator(4, Tb); - SU::generator(6, Tc); - - const RealD a(0.6), b(0.9), c(0.7); - - ColourMatrixD SU3_field; - arg = timesI(a*Ta + b*Tb + c*Tc); - SU3_field = Exponentiate(arg, 2.0); // what is RealD alpha?? - - LorentzColourMatrixD test=Zero(); - - for(int mu=0; mu - Gauge3x2unmunger unmunger; - LorentzColour2x3D test_rdc = Zero(); - unmunger(test,test_rdc); - - // reconstruct full field. Gauge3x2munger - Gauge3x2munger munger; - LorentzColourMatrixD test_recon; - munger(test_rdc,test_recon); - - // round-trip test - assert(norm2(test_recon-test)<1e-10); - -} - -void checkGaugeSUmungers() { - ColourMatrixD Ta, Tb, Tc, arg; - SU::generator(1, Ta); - SU::generator(4, Tb); - SU::generator(6, Tc); - - const RealD a(0.6), b(0.9), c(0.7); - - ColourMatrixD SUN_field; - arg = timesI(a*Ta + b*Tb + c*Tc); - SUN_field = Exponentiate(arg, 2.0); // what is RealD alpha?? - - LorentzColourMatrixD test=Zero(); - - for(int mu=0; mu - GaugeSUunmunger unmunger; - LorentzColour2x3D test_rdc = Zero(); - unmunger(test,test_rdc); - - // reconstruct full field. GaugeSUmunger - GaugeSUmunger munger; - LorentzColourMatrixD test_recon; - munger(test_rdc,test_recon); - - // round-trip test - assert(norm2(test_recon-test)<1e-15); - - // check result is unitary and det==1 - for(int mu=0;mu(test_recon, mu); - auto det = Determinant(new_cm); - assert( abs(norm2(det)-1.0) < 1e-15 ); - assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); - } - -} - -void checkGaugeSpmungers() { - ColourMatrix Sp_field; - - const Complex a(0.5, 0.5), abar(0.5, -0.5); - const Complex b(0.3, 0.9), bbar(0.3, -0.9); - - // fill top left - for(int i=0;i - GaugeSpunmunger unmunger; - LorentzColourNx2ND test_rdc = Zero(); - unmunger(test,test_rdc); - - // reconstruct full field. GaugeSpmunger - GaugeSpmunger munger; - LorentzColourMatrixD test_recon; - munger(test_rdc,test_recon); - - // round-trip test - assert(test==test_recon); - -} - int main (int argc, char ** argv) { #ifdef HAVE_LIME Grid_init(&argc,&argv); - + std::cout <1 && Nc<4) { - std::cout << GridLogMessage << "Checking SU(2/3) mungers " << std::endl; + std::cout << GridLogMessage << "Checking SU(" << Nc << ") mungers " << std::endl; check_reconstruct3(); checkGauge3x2mungers(); } From ee0d2ff17de4cc021462b4842a08c585850ebeba Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Mon, 23 Feb 2026 17:31:35 +0000 Subject: [PATCH 33/67] munger unit tests now return bool and output success msg --- tests/IO/Test_io_mungers.cc | 115 ++++++++++++++++++++++++------------ 1 file changed, 76 insertions(+), 39 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index b11ec339c8..24e81a66c7 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -34,15 +34,15 @@ using namespace Grid; // tests to check the various (un)munger classes // and functions in parallelIO/Metadata.h -void check_reconstruct3(); -void check_reconstructSU(); -void checkGauge3x2mungers(); -void checkGaugeSUmungers(); -void check_reconstructSp(); -void checkGaugeSpmungers(); -void checkBinarySimpleMungers(); -void checkGaugeSimpleMungers(); -void checkGaugeDoubleStoredMungers(); +bool check_reconstruct3(); +bool check_reconstructSU(); +bool checkGauge3x2mungers(); +bool checkGaugeSUmungers(); +bool check_reconstructSp(); +bool checkGaugeSpmungers(); +bool checkBinarySimpleMungers(); +bool checkGaugeSimpleMungers(); +bool checkGaugeDoubleStoredMungers(); void mock_SU_field(); @@ -76,7 +76,7 @@ void mock_SU_field( ColourMatrixD &cm, std::vector seed ) { } //////////////////////////////////////////// -void check_reconstruct3() { +bool check_reconstruct3() { ColourMatrixD SU3_field; std::vector rng_seed = {1,2,3,4}; @@ -103,9 +103,11 @@ void check_reconstruct3() { assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); } + return true; + } -void check_reconstructSU() { +bool check_reconstructSU() { ColourMatrixD SUN_field; std::vector rnd_seed = {4,8,15,16}; @@ -131,9 +133,12 @@ void check_reconstructSU() { assert( abs( norm2(det)-1.0 ) < 1e-15 ); assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); } + + return true; + } -void checkGauge3x2mungers() { +bool checkGauge3x2mungers() { ColourMatrixD SU3_field; std::vector rnd_seed{20,6,15,8}; @@ -158,9 +163,11 @@ void checkGauge3x2mungers() { //std::cout << "checkGauge3x2mungers: " << norm2(test_recon-test) << std::endl; assert(norm2(test_recon-test)<1e-10); + return true; + } -void checkGaugeSUmungers() { +bool checkGaugeSUmungers() { ColourMatrixD SUN_field; std::vector rnd_seed{23,42,66,98}; mock_SU_field(SUN_field, rnd_seed); @@ -193,10 +200,12 @@ void checkGaugeSUmungers() { assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); } + return true; + } -void check_reconstructSp() { +bool check_reconstructSp() { LorentzColourMatrix cm = Zero(); LorentzColourMatrix cm_zero = Zero(); @@ -244,9 +253,11 @@ void check_reconstructSp() { } } + return true; + } -void checkGaugeSpmungers() { +bool checkGaugeSpmungers() { ColourMatrix Sp_field; const Complex a(0.5, 0.5), abar(0.5, -0.5); @@ -296,11 +307,13 @@ void checkGaugeSpmungers() { // round-trip test assert(test==test_recon); + return true; + } // mungeing between the same type should yield the same result. -void checkBinarySimpleMungers() { +bool checkBinarySimpleMungers() { LorentzColourMatrixF in_scalar_objectF; // single precision LorentzColourMatrixF out_scalar_objectF = Zero(); @@ -328,9 +341,11 @@ void checkBinarySimpleMungers() { mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); + return true; + } // these are used in NerscIO.h -void checkGaugeSimpleMungers() { +bool checkGaugeSimpleMungers() { LorentzColourMatrixF in_scalar_objectF; // single precision LorentzColourMatrixF out_scalar_objectF = Zero(); @@ -358,9 +373,11 @@ void checkGaugeSimpleMungers() { mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); + return true; + } -void checkGaugeDoubleStoredMungers() { +bool checkGaugeDoubleStoredMungers() { DoubleStoredColourMatrixF in_scalar_objectF; // single precision DoubleStoredColourMatrixF out_scalar_objectF = Zero(); @@ -388,6 +405,8 @@ void checkGaugeDoubleStoredMungers() { mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); + return true; + } @@ -398,33 +417,51 @@ int main (int argc, char ** argv) std::cout <1 && Nc<4) { - std::cout << GridLogMessage << "Checking SU(" << Nc << ") mungers " << std::endl; - check_reconstruct3(); - checkGauge3x2mungers(); + if( check_reconstruct3() ) { + std::cout << GridLogMessage << "reconstruct3: PASS" << std::endl; + } + + if( checkGauge3x2mungers() ) { + std::cout << GridLogMessage << "Gauge3x2mungers: PASS" << std::endl; + } + } + + if( check_reconstructSU() ) { + std::cout << GridLogMessage << "reconstructSU: PASS" << std::endl; + } + if( checkGaugeSUmungers() ) { + std::cout << GridLogMessage << "GaugeSUmungers: PASS" << std::endl; } if constexpr(Nc>2 && Nc%2==0) { - std::cout << GridLogMessage << "Checking Sp mungers " << std::endl; - check_reconstructSp(); - checkGaugeSpmungers(); + std::cout < Date: Tue, 24 Feb 2026 01:54:33 +0000 Subject: [PATCH 34/67] fix 'Source file' in pre-amble --- Grid/parallelIO/IldgIOtypes.h | 2 +- Grid/parallelIO/MetaData.h | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/Grid/parallelIO/IldgIOtypes.h b/Grid/parallelIO/IldgIOtypes.h index 4fae9d9dac..4fa14f7b7d 100644 --- a/Grid/parallelIO/IldgIOtypes.h +++ b/Grid/parallelIO/IldgIOtypes.h @@ -2,7 +2,7 @@ Grid physics library, www.github.com/paboyle/Grid -Source file: ./lib/parallelIO/IldgIO.h +Source file: ./lib/parallelIO/IldgIOtypes.h Copyright (C) 2015 diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index f733ecb8e6..756a6eb8f0 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -2,7 +2,7 @@ Grid physics library, www.github.com/paboyle/Grid - Source file: ./lib/parallelIO/NerscIO.h + Source file: ./lib/parallelIO/Metadata.h Copyright (C) 2015 From ae2cc2dfcd3eebdbb82167a1349c3c98bb228463 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 24 Feb 2026 01:56:42 +0000 Subject: [PATCH 35/67] separated out and improved mocking of SU/Sp fields. fixed bug in reconstructSp for Nc!=4. added more tests for BinarySimple/GaugeSimple/GaugeDoubleStored mungers. --- tests/IO/Test_io_mungers.cc | 373 ++++++++++++++++++++++-------------- 1 file changed, 233 insertions(+), 140 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 24e81a66c7..f8b999c94c 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -45,16 +45,17 @@ bool checkGaugeSimpleMungers(); bool checkGaugeDoubleStoredMungers(); void mock_SU_field(); +void mock_Sp_field(); ////////////////////////////////////////////// -// helper function for generating SU matrices +// helper functions for generating (mock) SU/Sp matrices ////////////////////////////////////////////// void mock_SU_field( ColourMatrixD &cm, std::vector seed ) { + GridSerialRNG sRNG; sRNG.SeedFixedIntegers(seed); - ComplexD ca; ColourMatrixD lie, la, ta; - ComplexD ci(0.0, 1.0); + ComplexD ca; lie = Zero(); @@ -73,31 +74,57 @@ void mock_SU_field( ColourMatrixD &cm, std::vector seed ) { cm = Exponentiate(lie, 2.0); +} + +void mock_Sp_field( ColourMatrixD &cm, std::vector seed ) { + + GridSerialRNG sRNG; sRNG.SeedFixedIntegers(seed); + + std::vector elem(Nc*Nc/2); + for(auto&& e: elem) { + random(sRNG,e); + } + + // fill upper left and upper right blocks only + for(int i=0;i rng_seed = {1,2,3,4}; - mock_SU_field(SU3_field, rng_seed); + ColourMatrixD SU3_xfield, SU3_yfield, SU3_zfield, SU3_tfield; + + mock_SU_field(SU3_xfield, std::vector({30})); + mock_SU_field(SU3_yfield, std::vector({95})); + mock_SU_field(SU3_zfield, std::vector({7})); + mock_SU_field(SU3_tfield, std::vector({10})); //set last row equal to zero for(int j=0;j(test, mu); + auto new_cm = peekIndex(scalar, mu); auto det = Determinant(new_cm); assert( abs(norm2(det)-1.0) < 1e-15 ); assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); @@ -108,27 +135,34 @@ bool check_reconstruct3() { } bool check_reconstructSU() { - - ColourMatrixD SUN_field; - std::vector rnd_seed = {4,8,15,16}; - mock_SU_field(SUN_field, rnd_seed); + + ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; + + mock_SU_field(SUN_xfield, std::vector({4})); + mock_SU_field(SUN_yfield, std::vector({8})); + mock_SU_field(SUN_zfield, std::vector({15})); + mock_SU_field(SUN_tfield, std::vector({16})); //set last row equal to zero for(int j=0;j(test, mu); + auto new_cm = peekIndex(scalar, mu); auto det = Determinant(new_cm); assert( abs( norm2(det)-1.0 ) < 1e-15 ); assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); @@ -140,65 +174,65 @@ bool check_reconstructSU() { bool checkGauge3x2mungers() { - ColourMatrixD SU3_field; - std::vector rnd_seed{20,6,15,8}; - mock_SU_field(SU3_field, rnd_seed); + ColourMatrixD SU3_xfield, SU3_yfield, SU3_zfield, SU3_tfield; - LorentzColourMatrixD test=Zero(); - for(int mu=0; mu({66})); + mock_SU_field(SU3_yfield, std::vector({63})); + mock_SU_field(SU3_zfield, std::vector({71})); + mock_SU_field(SU3_tfield, std::vector({98})); + + LorentzColourMatrixD scalar = Zero(); + + pokeLorentz(scalar, SU3_xfield, 0); + pokeLorentz(scalar, SU3_yfield, 1); + pokeLorentz(scalar, SU3_zfield, 2); + pokeLorentz(scalar, SU3_tfield, 3); // reduce field. Gauge3x2unmunger Gauge3x2unmunger unmunger; - LorentzColour2x3D test_rdc; - unmunger(test,test_rdc); + LorentzColour2x3D scalar_rdc; + unmunger(scalar,scalar_rdc); // reconstruct full field. Gauge3x2munger Gauge3x2munger munger; - LorentzColourMatrixD test_recon; - munger(test_rdc,test_recon); + LorentzColourMatrixD scalar_recon; + munger(scalar_rdc,scalar_recon); // round-trip test - //std::cout << "checkGauge3x2mungers: " << norm2(test_recon-test) << std::endl; - assert(norm2(test_recon-test)<1e-10); + assert(norm2(scalar_recon-scalar)<1e-10); return true; } bool checkGaugeSUmungers() { - ColourMatrixD SUN_field; - std::vector rnd_seed{23,42,66,98}; - mock_SU_field(SUN_field, rnd_seed); - LorentzColourMatrixD test=Zero(); + ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; - for(int mu=0; mu({360})); + mock_SU_field(SUN_yfield, std::vector({804})); + mock_SU_field(SUN_zfield, std::vector({77})); + mock_SU_field(SUN_tfield, std::vector({65})); + + LorentzColourMatrixD scalar = Zero(); + + pokeLorentz(scalar, SUN_xfield, 0); + pokeLorentz(scalar, SUN_yfield, 1); + pokeLorentz(scalar, SUN_zfield, 2); + pokeLorentz(scalar, SUN_tfield, 3); // reduce field. GaugeSUunmunger GaugeSUunmunger unmunger; - LorentzColour2x3D test_rdc = Zero(); - unmunger(test,test_rdc); + LorentzColour2x3D scalar_rdc = Zero(); + unmunger(scalar,scalar_rdc); // reconstruct full field. GaugeSUmunger GaugeSUmunger munger; - LorentzColourMatrixD test_recon; - munger(test_rdc,test_recon); + LorentzColourMatrixD scalar_recon; + munger(scalar_rdc,scalar_recon); // round-trip test - //std::cout << "checkGaugeSUmungers: " << norm2(test_recon-test) << std::endl; - assert(norm2(test_recon-test)<1e-15); - - // check result is unitary and det==1 - for(int mu=0;mu(test_recon, mu); - auto det = Determinant(new_cm); - assert( abs(norm2(det)-1.0) < 1e-15 ); - assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); - } + assert(norm2(scalar_recon-scalar)<1e-15); return true; @@ -207,44 +241,28 @@ bool checkGaugeSUmungers() { bool check_reconstructSp() { - LorentzColourMatrix cm = Zero(); - LorentzColourMatrix cm_zero = Zero(); - - // should be zero matrix after reconstruction. - reconstructSp(cm); - assert(cm==cm_zero); - - ColourMatrix Sp_field = Zero(); - - const Complex a(0.5, 0.5), abar(0.5, -0.5); - const Complex b(0.3, 0.9), bbar(0.3, -0.9); + ColourMatrixD Sp_xfield, Sp_yfield, Sp_zfield, Sp_tfield; + Sp_xfield = Sp_yfield = Sp_zfield = Sp_tfield = Zero(); - // fill top left - for(int i=0;i({1066})); + mock_Sp_field(Sp_yfield, std::vector({1995})); + mock_Sp_field(Sp_zfield, std::vector({997})); + mock_Sp_field(Sp_tfield, std::vector({58})); - LorentzColourMatrixD test = Zero(); + LorentzColourMatrixD scalar = Zero(); - for(int mu=0; mu(test,mu); + auto Sp_cm = peekIndex(scalar,mu); for(int i=0;i({23})); + mock_Sp_field(Sp_yfield, std::vector({42})); + mock_Sp_field(Sp_zfield, std::vector({93})); + mock_Sp_field(Sp_tfield, std::vector({49})); - LorentzColourMatrixD test = Zero(); + LorentzColourMatrixD scalar = Zero(); - for(int mu=0; mu GaugeSpunmunger unmunger; - LorentzColourNx2ND test_rdc = Zero(); - unmunger(test,test_rdc); + LorentzColourNx2ND scalar_rdc = Zero(); + unmunger(scalar,scalar_rdc); // reconstruct full field. GaugeSpmunger GaugeSpmunger munger; - LorentzColourMatrixD test_recon; - munger(test_rdc,test_recon); + LorentzColourMatrixD scalar_recon; + munger(scalar_rdc,scalar_recon); // round-trip test - assert(test==test_recon); + assert(scalar==scalar_recon); return true; } -// mungeing between the same type should yield the same result. bool checkBinarySimpleMungers() { LorentzColourMatrixF in_scalar_objectF; // single precision LorentzColourMatrixF out_scalar_objectF = Zero(); LorentzColourMatrixD in_scalar_objectD; // double precision LorentzColourMatrixD out_scalar_objectD = Zero(); - - in_scalar_objectF = 1.0; - in_scalar_objectD = 1.0; - + + GridSerialRNG sRNG; sRNG.SeedFixedIntegers(std::vector({117})); + + random(sRNG,in_scalar_objectF); + random(sRNG,in_scalar_objectD); + // BinarySimpleUnmunger BinarySimpleUnmunger unmungerFF; - BinarySimpleUnmunger unmungerDF; + BinarySimpleUnmunger unmungerFD; BinarySimpleUnmunger unmungerDD; - BinarySimpleUnmunger unmungerFD; + BinarySimpleUnmunger unmungerDF; unmungerFF(in_scalar_objectF, out_scalar_objectF); + unmungerDD(in_scalar_objectD, out_scalar_objectD); + assert(in_scalar_objectF==out_scalar_objectF); + assert(in_scalar_objectD==out_scalar_objectD); + + random(sRNG,in_scalar_objectF); + unmungerFD(in_scalar_objectF, out_scalar_objectD); + unmungerDF(out_scalar_objectD, out_scalar_objectF); assert(in_scalar_objectF==out_scalar_objectF); + random(sRNG,in_scalar_objectD); + unmungerDF(in_scalar_objectD, out_scalar_objectF); + unmungerFD(out_scalar_objectF, out_scalar_objectD); + // lose exactness when going from double-->single precision + assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); + // BinarySimpleMunger BinarySimpleMunger mungerFF; BinarySimpleMunger mungerDF; BinarySimpleMunger mungerDD; BinarySimpleMunger mungerFD; + random(sRNG,in_scalar_objectF); + random(sRNG,in_scalar_objectD); + + mungerFF(in_scalar_objectF, out_scalar_objectF); mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); + assert(in_scalar_objectF==out_scalar_objectF); + random(sRNG,in_scalar_objectF); + mungerFD(in_scalar_objectF, out_scalar_objectD); + mungerDF(out_scalar_objectD, out_scalar_objectF); + assert(in_scalar_objectF==out_scalar_objectF); + + random(sRNG,in_scalar_objectD); + mungerDF(in_scalar_objectD, out_scalar_objectF); + mungerFD(out_scalar_objectF, out_scalar_objectD); + // lose exactness when going from double-->single precision + assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); + return true; } + // these are used in NerscIO.h bool checkGaugeSimpleMungers() { @@ -352,26 +383,57 @@ bool checkGaugeSimpleMungers() { LorentzColourMatrixD in_scalar_objectD; // double precision LorentzColourMatrixD out_scalar_objectD = Zero(); - in_scalar_objectF = 1.0; - in_scalar_objectD = 1.0; + GridSerialRNG sRNG; sRNG.SeedFixedIntegers(std::vector({57})); + + random(sRNG,in_scalar_objectF); + random(sRNG,in_scalar_objectD); // GaugeSimpleUnmunger GaugeSimpleUnmunger unmungerFF; - GaugeSimpleUnmunger unmungerDF; + GaugeSimpleUnmunger unmungerFD; GaugeSimpleUnmunger unmungerDD; - GaugeSimpleUnmunger unmungerFD; + GaugeSimpleUnmunger unmungerDF; unmungerFF(in_scalar_objectF, out_scalar_objectF); + unmungerDD(in_scalar_objectD, out_scalar_objectD); + assert(in_scalar_objectF==out_scalar_objectF); + assert(in_scalar_objectD==out_scalar_objectD); + + random(sRNG,in_scalar_objectF); + unmungerFD(in_scalar_objectF, out_scalar_objectD); + unmungerDF(out_scalar_objectD, out_scalar_objectF); assert(in_scalar_objectF==out_scalar_objectF); + random(sRNG,in_scalar_objectD); + unmungerDF(in_scalar_objectD, out_scalar_objectF); + unmungerFD(out_scalar_objectF, out_scalar_objectD); + // lose exactness when going from double-->single precision + assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); + // GaugeSimpleMunger GaugeSimpleMunger mungerFF; GaugeSimpleMunger mungerDF; GaugeSimpleMunger mungerDD; GaugeSimpleMunger mungerFD; - + + random(sRNG,in_scalar_objectF); + random(sRNG,in_scalar_objectD); + + mungerFF(in_scalar_objectF, out_scalar_objectF); mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); + assert(in_scalar_objectF==out_scalar_objectF); + + random(sRNG,in_scalar_objectF); + mungerFD(in_scalar_objectF, out_scalar_objectD); + mungerDF(out_scalar_objectD, out_scalar_objectF); + assert(in_scalar_objectF==out_scalar_objectF); + + random(sRNG,in_scalar_objectD); + mungerDF(in_scalar_objectD, out_scalar_objectF); + mungerFD(out_scalar_objectF, out_scalar_objectD); + // lose exactness when going from double-->single precision + assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); return true; @@ -383,27 +445,58 @@ bool checkGaugeDoubleStoredMungers() { DoubleStoredColourMatrixF out_scalar_objectF = Zero(); DoubleStoredColourMatrixD in_scalar_objectD; // double precision DoubleStoredColourMatrixD out_scalar_objectD = Zero(); - - in_scalar_objectF = 1.0; - in_scalar_objectD = 1.0; + + GridSerialRNG sRNG; sRNG.SeedFixedIntegers(std::vector({169})); + + random(sRNG,in_scalar_objectF); + random(sRNG,in_scalar_objectD); // GaugeDoubleStoredUnmunger GaugeDoubleStoredUnmunger unmungerFF; - GaugeDoubleStoredUnmunger unmungerDF; + GaugeDoubleStoredUnmunger unmungerFD; GaugeDoubleStoredUnmunger unmungerDD; - GaugeDoubleStoredUnmunger unmungerFD; + GaugeDoubleStoredUnmunger unmungerDF; unmungerFF(in_scalar_objectF, out_scalar_objectF); + unmungerDD(in_scalar_objectD, out_scalar_objectD); + assert(in_scalar_objectF==out_scalar_objectF); + assert(in_scalar_objectD==out_scalar_objectD); + + random(sRNG,in_scalar_objectF); + unmungerFD(in_scalar_objectF, out_scalar_objectD); + unmungerDF(out_scalar_objectD, out_scalar_objectF); assert(in_scalar_objectF==out_scalar_objectF); + random(sRNG,in_scalar_objectD); + unmungerDF(in_scalar_objectD, out_scalar_objectF); + unmungerFD(out_scalar_objectF, out_scalar_objectD); + // lose exactness when going from double-->single precision + assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); + // GaugeDoubleStoredMunger GaugeDoubleStoredMunger mungerFF; GaugeDoubleStoredMunger mungerDF; GaugeDoubleStoredMunger mungerDD; GaugeDoubleStoredMunger mungerFD; + random(sRNG,in_scalar_objectF); + random(sRNG,in_scalar_objectD); + + mungerFF(in_scalar_objectF, out_scalar_objectF); mungerDD(in_scalar_objectD, out_scalar_objectD); assert(in_scalar_objectD==out_scalar_objectD); + assert(in_scalar_objectF==out_scalar_objectF); + + random(sRNG,in_scalar_objectF); + mungerFD(in_scalar_objectF, out_scalar_objectD); + mungerDF(out_scalar_objectD, out_scalar_objectF); + assert(in_scalar_objectF==out_scalar_objectF); + + random(sRNG,in_scalar_objectD); + mungerDF(in_scalar_objectD, out_scalar_objectF); + mungerFD(out_scalar_objectF, out_scalar_objectD); + // lose exactness when going from double-->single precision + assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); return true; From e46ab610c8b97b00789c642598fc902d6fcb173a Mon Sep 17 00:00:00 2001 From: gray95 Date: Fri, 27 Feb 2026 20:29:59 +0000 Subject: [PATCH 36/67] loosens norm2 condition tolerance to 1e-12 --- tests/IO/Test_io_mungers.cc | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index f8b999c94c..fe6e284cf0 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -126,8 +126,8 @@ bool check_reconstruct3() { for(int mu=0;mu(scalar, mu); auto det = Determinant(new_cm); - assert( abs(norm2(det)-1.0) < 1e-15 ); - assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); + assert( abs(norm2(det)-1.0) < 1e-12 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-12 ); } return true; @@ -164,8 +164,8 @@ bool check_reconstructSU() { for(int mu=0;mu(scalar, mu); auto det = Determinant(new_cm); - assert( abs( norm2(det)-1.0 ) < 1e-15 ); - assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-15 ); + assert( abs( norm2(det)-1.0 ) < 1e-12 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-12 ); } return true; @@ -232,7 +232,7 @@ bool checkGaugeSUmungers() { munger(scalar_rdc,scalar_recon); // round-trip test - assert(norm2(scalar_recon-scalar)<1e-15); + assert(norm2(scalar_recon-scalar)<1e-12); return true; From 8d5d8f09feac02f6360ed82474c99aa219d2b1d5 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 18:53:06 +0000 Subject: [PATCH 37/67] add unique_su template parameter to readConfiguration and readLatticeBinaryObject --- Grid/parallelIO/IldgIO.h | 30 +++++++++++++++--------------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 3a13a60775..04420630fe 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -737,35 +737,35 @@ class IldgReader : public GridLimeReader { ////////////////////////////////////////////////////////////////// // this helper function wraps the logic for choosing // the correct munger and intermediate lattice data type - // before reading a lattice field from disk. + // when reading a lattice field from disk. // This is a runtime choice as Grid has to first read the // header of the lattice cfg. Therefore templating this function on - // FloatingPointFormat etc. does not work. It is rather - // cumbersome, it's if statement with 6 branches, + // FloatingPointFormat etc. does not work. It is a rather + // cumbersome if statement with 6 branches, // but the benefit of organising IldgReader this way is it makes // readConfiguration clearer and more readable. ////////////////////////////////////////////////////////////////// - template + template void readLatticeBinaryObject(Lattice &Umu, std::string filename, FloatingPointFormat fp_fmt, MatrixFormat matrix_fmt, bool is_grp_su, bool is_grp_sp, uint64_t &offset, uint32_t &nersc_csum, uint32_t &scidac_csuma, uint32_t &scidac_csumb) { - // need all the types we could possibly read from - // including the intermediate data types for - // reduced format lattices and single/double precision + // need all the types we could possibly read from + // including the intermediate data types for + // reduced format lattices and single/double precision typedef typename vobj::scalar_object sobj; typedef LorentzColourMatrixF fobj; typedef LorentzColourMatrixD dobj; - typedef LorentzColour2x3F fobjsuR; - typedef LorentzColour2x3D dobjsuR; - typedef LorentzColourNx2NF fobjspR; - typedef LorentzColourNx2ND dobjspR; + typedef LorentzColour2x3F fobjsuR; + typedef LorentzColour2x3D dobjsuR; + typedef LorentzColourNx2NF fobjspR; + typedef LorentzColourNx2ND dobjspR; std::string format; if ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { format = std::string("IEEE64BIG"); if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSUmunger munge; + GaugeSUmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; @@ -777,7 +777,7 @@ class IldgReader : public GridLimeReader { } else if ( fp_fmt == FloatingPointFormat::IEEE32BIG) { format = std::string("IEEE32BIG"); if ( is_grp_su && matrix_fmt==MatrixFormat::REDUCED ) { - GaugeSUmunger munge; + GaugeSUmunger munge; BinaryIO::readLatticeObject(Umu, filename, munge, offset, format,nersc_csum,scidac_csuma,scidac_csumb); } else if ( is_grp_sp && matrix_fmt==MatrixFormat::REDUCED ) { GaugeSpmunger munge; @@ -797,7 +797,7 @@ class IldgReader : public GridLimeReader { // Else use ILDG MetaData object if present. // Else use SciDAC MetaData object if present. //////////////////////////////////////////////////////////////// - template + template void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_) { GridBase *grid = Umu.Grid(); @@ -964,7 +964,7 @@ class IldgReader : public GridLimeReader { uint64_t offset= ftello(File); - readLatticeBinaryObject(Umu, filename, fp_fmt, matrix_fmt, is_grp_su, is_grp_sp, offset, nersc_csum, scidac_csuma, scidac_csumb); + readLatticeBinaryObject(Umu, filename, fp_fmt, matrix_fmt, is_grp_su, is_grp_sp, offset, nersc_csum, scidac_csuma, scidac_csumb); found_ildgBinary = 1; From 2a7a4c4910facd7e40b620ca65088e418421f3bd Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 18:55:47 +0000 Subject: [PATCH 38/67] adds unique_reconstructSU as alternative munger and adds corresponding template specialisation to GaugeSUmunger, specialising on unique_su --- Grid/parallelIO/MetaData.h | 135 ++++++++++++++++++++++++++++++------- 1 file changed, 110 insertions(+), 25 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 756a6eb8f0..f8e0d2c990 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -229,23 +229,88 @@ inline void reconstructSU(LorentzColourMatrix &cm) { using ColourMatrixNm = iScalar > > ; - ColourMatrixNm red; // for the Nc-1 by Nc-1 matrix - ColourMatrix old; // for the Nc-1 by N matrix read from disk + ColourMatrixNm Su; // for the Nc-1 by Nc-1 matrix + ColourMatrix SU; // for the Nc-1 by N matrix read from disk for(int mu=0; mu(cm,mu); + Su = Zero(); + SU = peekIndex(cm,mu); for(int k=0; k(cm,old,mu); + pokeIndex(cm,SU,mu); + } +} + + +bool is_perm_even(std::vector &v) { + + int n = v.size(); + std::vector a(n,0); + int c = 0; + + for(int j=0; j v(Nc); + std::iota(v.begin(), v.end(), 0); // {0,1,...,Nc-1} + + for(int mu=0; mu c = v; + c.erase(c.begin()+k); // {0,...,k-1,k+1,...Nc-1} + cm(mu)()(Nc-1,k) = 0; + // cycle through perms of c + do + { + std::vector P = c; + P.emplace_back(k); + bool E = is_perm_even(P); + ComplexD prod(1,0); + for(int i=0;i +struct GaugeSUmunger; +// two specialisations for two ways to reconstruct +// NcxNc matrix from Nc-1xNc matrix + template -struct GaugeSUmunger{ +struct GaugeSUmunger{ void operator() (fobj &in,sobj &out){ for(int mu=0;mu +struct GaugeSUmunger{ + void operator() (fobj &in,sobj &out){ + for(int mu=0;mu struct GaugeSUunmunger{ void operator() (sobj &in,fobj &out){ From 073511228d57574a0fd16bd243ebf52af3ebee44 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 18:57:35 +0000 Subject: [PATCH 39/67] adds unit tests for unique_reconstructSU and is_perm_even --- tests/IO/Test_io_mungers.cc | 82 ++++++++++++++++++++++++++++++++++--- 1 file changed, 77 insertions(+), 5 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index fe6e284cf0..3b5d300ebe 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -29,13 +29,14 @@ Author: Gaurav Ray /* END LEGAL */ #include -using namespace std; using namespace Grid; // tests to check the various (un)munger classes // and functions in parallelIO/Metadata.h bool check_reconstruct3(); bool check_reconstructSU(); +bool check_is_perm_even(); +bool check_unique_reconstructSU(); bool checkGauge3x2mungers(); bool checkGaugeSUmungers(); bool check_reconstructSp(); @@ -172,6 +173,65 @@ bool check_reconstructSU() { } +bool check_is_perm_even() { + + std::vector v = {0,1,2,3}; + int count = 0; + + // swap two elements and check if parity flipped. + for(auto& e: v) { + for(auto& d: v) { + if(e!=d) { + if(count%2==0) { assert( is_perm_even(v)==true ); } + else { assert( is_perm_even(v)==false ); } + std::swap(e, d); + count += 1; + } + } + } + + return true; + +} + +bool check_unique_reconstructSU() { + + ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; + + mock_SU_field(SUN_xfield, std::vector({2})); + mock_SU_field(SUN_yfield, std::vector({40})); + mock_SU_field(SUN_zfield, std::vector({28})); + mock_SU_field(SUN_tfield, std::vector({25})); + + //set last row equal to zero + for(int j=0;j(scalar, mu); + auto det = Determinant(new_cm); + assert( abs( norm2(det)-1.0 ) < 1e-12 ); + assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-12 ); + } + + return true; + +} + bool checkGauge3x2mungers() { ColourMatrixD SU3_xfield, SU3_yfield, SU3_zfield, SU3_tfield; @@ -205,6 +265,7 @@ bool checkGauge3x2mungers() { } +template bool checkGaugeSUmungers() { ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; @@ -227,9 +288,9 @@ bool checkGaugeSUmungers() { unmunger(scalar,scalar_rdc); // reconstruct full field. GaugeSUmunger - GaugeSUmunger munger; + GaugeSUmunger munger; LorentzColourMatrixD scalar_recon; - munger(scalar_rdc,scalar_recon); + munger(scalar_rdc, scalar_recon); // round-trip test assert(norm2(scalar_recon-scalar)<1e-12); @@ -510,6 +571,10 @@ int main (int argc, char ** argv) std::cout <1 && Nc<4) { if( check_reconstruct3() ) { @@ -524,8 +589,15 @@ int main (int argc, char ** argv) if( check_reconstructSU() ) { std::cout << GridLogMessage << "reconstructSU: PASS" << std::endl; } - if( checkGaugeSUmungers() ) { - std::cout << GridLogMessage << "GaugeSUmungers: PASS" << std::endl; + if( checkGaugeSUmungers() ) { + std::cout << GridLogMessage << "(unique_su=false) GaugeSUmungers: PASS" << std::endl; + } + + if( check_unique_reconstructSU() ) { + std::cout << GridLogMessage << "unique_reconstructSU: PASS" << std::endl; + } + if( checkGaugeSUmungers() ) { + std::cout << GridLogMessage << "(unique_su=true) GaugeSUmungers: PASS" << std::endl; } if constexpr(Nc>2 && Nc%2==0) { From 079fe8c4c1e116a3d31971fbde5b739a3b037614 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 18:59:41 +0000 Subject: [PATCH 40/67] adds examples showing the use of the unique_su template parameter in readConfiguration --- tests/IO/Test_ildg_reducedfmt_io.cc | 35 +++++++++++++++++------------ 1 file changed, 21 insertions(+), 14 deletions(-) diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index 971e77bc83..ac9ed3a8fd 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -29,7 +29,6 @@ Author: Gaurav Ray /* END LEGAL */ #include -using namespace std; using namespace Grid; // this test demonstrates and checks IldgWriter/Readers' ability to @@ -38,7 +37,7 @@ using namespace Grid; ////////////////////////////////////////////// // this template function returns a // LatticeGaugeField of the chosen gauge -// group passed as a template argument. +// group (passed as a template argument) ////////////////////////////////////////////// template LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Grid, std::vector seed ) { @@ -59,13 +58,14 @@ LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Grid, std::vecto return Umu; } -/////////////////////////////////////////////////////////////// -// this template function generates writes a lattice +///////////////////////////////////////////////////////////// +// this template function writes a lattice // field of a given gaugeGroup to disk. It can write in -// reduced format and single precision depending on -// the values of matrix_fmt and fp_fmt. -/////////////////////////////////////////////////////////////// -template +// reduced format and single/double precision depending on +// the values of matrix_fmt and fp_fmt. unique_su toggles +// between different functions when reading reduced fmt lattices. +///////////////////////////////////////////////////////////// +template void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, std::string file) { if constexpr( std::is_same_v && N%2==1) { @@ -93,7 +93,7 @@ void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, st std::cout <(Umu_saved, header); _IldgReader.close(); std::cout <(Umu, Grid, "./ckpoint_su64_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su64_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su64_unique_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su32_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); - writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su32_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su32_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(Umu, Grid, "./ckpoint_su32_unique_"+std::to_string(Nc-1)+"x"+std::to_string(Nc)+".4000"); - // write and read Sp fields - writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp64_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); + // write and read Sp lattices + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp64_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp64_"+std::to_string(Nc/2)+"x"+std::to_string(Nc)+".4000"); + writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp32_"+std::to_string(Nc)+"x"+std::to_string(Nc)+".4000"); writeReadIldgConfiguration(UmuSp, Grid, "./ckpoint_sp32_"+std::to_string(Nc/2)+"x"+std::to_string(Nc)+".4000"); From b7b572c1bebc31ee3a6fa925ad522f839da8de66 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 19:02:28 +0000 Subject: [PATCH 41/67] removes undefined function --- tests/IO/Test_IldgWriter.cc | 2 -- 1 file changed, 2 deletions(-) diff --git a/tests/IO/Test_IldgWriter.cc b/tests/IO/Test_IldgWriter.cc index fd44c1cdb8..3e0858e6bf 100644 --- a/tests/IO/Test_IldgWriter.cc +++ b/tests/IO/Test_IldgWriter.cc @@ -29,13 +29,11 @@ Author: Gaurav Ray /* END LEGAL */ #include -using namespace std; using namespace Grid; // Unit tests to check the IldgWriter class. void checkWriteLimeIldgLFN(std::string &test_string); -void checkWriteConfiguration(LatticeGaugeField &Umu); // write a lime record and then read it back and check it matches. void checkWriteLimeIldgLFN(std::string &test_string) { From 44157bdea18e443d6f81dba729fbcfa0440b3760 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 19:25:48 +0000 Subject: [PATCH 42/67] updates documentation with new ILDG I/O functionality --- documentation/manual.rst | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/documentation/manual.rst b/documentation/manual.rst index 06251d752e..41510fc765 100644 --- a/documentation/manual.rst +++ b/documentation/manual.rst @@ -2104,7 +2104,7 @@ These are specialised to SciDAC writers, introducing facilities for generating t void readScidacFieldRecord(Lattice &field,userRecord &_userRecord); }; -They are also specialised to ILDG format writers, available and defined only for Gauge configurations. Following the publication of Rev. 1.2 of the ILDG Binary File Format by the ILDG Metadata Working Group the ILDG format writer and reader have been updated to be able to read and write SU(2/3) and Sp(2N) lattice configurations in reduced formats.:: +They are also specialised to ILDG format writers, available and defined only for Gauge configurations. Following the publication of Rev. 1.2 of the ILDG Binary File Format by the ILDG Metadata Working Group the `IldgWriter` and `IldgReader` have been updated to be able to read and write SU(N) and Sp(2N) lattice configurations in reduced formats.:: class IldgWriter : public ScidacWriter { @@ -2115,12 +2115,12 @@ They are also specialised to ILDG format writers, available and defined only for class IldgReader : public GridLimeReader { - template + template void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_); }; -The template parameters *group_name*, *MatrixFormat*, *FloatingPointFormat* specify the format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced or not reduced, and two options for the floating point format, 32 bit and 64 bit precision. Reduced format Sp gauge fields take up half as much memory as the non-reduced fields while SU(3) fields take up 2/3 as much memory as their non-reduced counterparts. +`writeConfiguration`'s template parameters; *group_name*, *MatrixFormat*, and *FloatingPointFormat* specify the format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced and not reduced, and two options for the floating point format, 32 bit and 64 bit precision. Reduced format Sp gauge fields take up half as much memory as non-reduced fields while SU(N) fields take up N-1/N as much memory as their non-reduced counterparts. When reading reduced format lattices from disk there are multiple ways the full SU(N) field can be recovered. The `unique_su` template parameter of `readConfiguration` controls which of two alternative functions is used. Set `unique_su=true` if you want `Grid` to follow the prescription laid out in the ILDG Binary File Format specification when reconstructing reduced format lattices. **Implementation detail** From 681d5e21711a96720d63db71cc79812370a8d4ae Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 21:54:00 +0000 Subject: [PATCH 43/67] clean up GridLogMessage --- Grid/parallelIO/IldgIO.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 04420630fe..6c07d2c520 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -870,9 +870,9 @@ class IldgReader : public GridLimeReader { if(doc.child("ildgFormat").child("rows")) { std::string rows = doc.child("ildgFormat").child("rows").child_value(); - std::cout << GridLogMessage << "rows element found in ildg-format and is " << rows << " so this lattice might be ildg 1.2 compliant." << std::endl; + std::cout << GridLogMessage << " element present = " << rows << ". So this lattice might be ildg 1.2 compliant." << std::endl; } else { - std::cout << GridLogMessage << "rows element is not present in ildg-format - adding it after ." << std::endl; + std::cout << GridLogMessage << " element not present - adding it after ." << std::endl; pugi::xml_node ildgfmt = doc.child("ildgFormat"); const std::string stNC = std::to_string(Nc); ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); From fb32952ee8a7edc3ce83f53ab47d697bc4339512 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 3 Mar 2026 23:52:14 +0000 Subject: [PATCH 44/67] use assign multiply op so that product in unique_reconstruct is formed left-->right as in ildg spec --- Grid/parallelIO/MetaData.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index f8e0d2c990..b75e090cab 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -304,7 +304,7 @@ inline void unique_reconstructSU(LorentzColourMatrix &cm) bool E = is_perm_even(P); ComplexD prod(1,0); for(int i=0;i Date: Wed, 11 Mar 2026 16:58:22 +0000 Subject: [PATCH 45/67] address pr comments including adding a link to ildg binary spec 1.2 --- documentation/manual.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/documentation/manual.rst b/documentation/manual.rst index 41510fc765..b42bd42ce6 100644 --- a/documentation/manual.rst +++ b/documentation/manual.rst @@ -2104,7 +2104,7 @@ These are specialised to SciDAC writers, introducing facilities for generating t void readScidacFieldRecord(Lattice &field,userRecord &_userRecord); }; -They are also specialised to ILDG format writers, available and defined only for Gauge configurations. Following the publication of Rev. 1.2 of the ILDG Binary File Format by the ILDG Metadata Working Group the `IldgWriter` and `IldgReader` have been updated to be able to read and write SU(N) and Sp(2N) lattice configurations in reduced formats.:: +They are also specialised to ILDG format writers, available and defined only for Gauge configurations. These ``IldgWriter`` and ``IldgReader`` template classes can read and write ILDG-compliant SU(N) and Sp(2N) lattice configurations as full NxN and 2Nx2N matrix fields respectively. They also support reduced format read/write with (N-1)xN matrix fields for SU(N) and Nx2N matrix fields for Sp(2N), as specified by `Rev. 1.2 `_ of the ILDG Binary File Format.:: class IldgWriter : public ScidacWriter { @@ -2120,7 +2120,7 @@ They are also specialised to ILDG format writers, available and defined only for }; -`writeConfiguration`'s template parameters; *group_name*, *MatrixFormat*, and *FloatingPointFormat* specify the format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced and not reduced, and two options for the floating point format, 32 bit and 64 bit precision. Reduced format Sp gauge fields take up half as much memory as non-reduced fields while SU(N) fields take up N-1/N as much memory as their non-reduced counterparts. When reading reduced format lattices from disk there are multiple ways the full SU(N) field can be recovered. The `unique_su` template parameter of `readConfiguration` controls which of two alternative functions is used. Set `unique_su=true` if you want `Grid` to follow the prescription laid out in the ILDG Binary File Format specification when reconstructing reduced format lattices. +``writeConfiguration``'s template parameters; ``group_name``, ``MatrixFormat``, and ``FloatingPointFormat`` specify the format the lattices are to be written to disk with. At present there are two options for the group, SU and Sp, two options for the matrix format, reduced and not reduced, and two options for the floating point format, 32 bit and 64 bit precision. Reduced format Sp gauge fields take up half as much storage as non-reduced fields while SU(N) fields take up N-1/N as much storage as their non-reduced counterparts. When reading reduced format lattices from disk there are multiple ways the full SU(N) field can be recovered. The ``unique_su`` template parameter of ``readConfiguration`` controls which of two alternative functions is used. Set ``unique_su=true`` if you want ``Grid`` to follow the prescription laid out in the ILDG Binary File Format specification when reconstructing reduced format lattices. **Implementation detail** From a53f6f8dfd9a55cc7b44494ba73b68911d58981b Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 11 Mar 2026 17:15:06 +0000 Subject: [PATCH 46/67] minor formatting changes --- Grid/parallelIO/IldgIO.h | 54 ++++++++++++++++++++-------------------- 1 file changed, 27 insertions(+), 27 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 6c07d2c520..7b4b844b14 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -447,9 +447,9 @@ class GridLimeWriter : public BinaryIO fflush(File); } - // std::cout << "W sizeof(sobj)" <_gsites<_gsites<(); - GaugeUnMunger unmunger; - + GaugeUnMunger unmunger; + BinaryIO::writeLatticeObject(field, filename, unmunger, offset1, format, nersc_csum, scidac_csuma, scidac_csumb, control); /////////////////////////////////////////// @@ -659,38 +659,38 @@ class IldgWriter : public ScidacWriter { ildgFormat ildgfmt ; const std::string stNC = std::to_string( Nc ) ; - // use the gauge group to populate the 'field' element of ildg header + // use the gauge group to populate the 'field' element of ildg header if constexpr ( is_su::value ) { - std::cout << GridLogMessage << "writing SU(" << stNC << ") field" << std::endl; - ildgfmt.field = std::string("su"+stNC+"gauge"); - } else if constexpr ( is_sp::value ) { - std::cout << GridLogMessage << "writing Sp(" << stNC << ") field" << std::endl; - ildgfmt.field = std::string("sp"+stNC+"gauge"); - } else { - static_assert(1, "Unrecognised group; unable to determine field tag"); - } + std::cout << GridLogMessage << "writing SU(" << stNC << ") field" << std::endl; + ildgfmt.field = std::string("su"+stNC+"gauge"); + } else if constexpr ( is_sp::value ) { + std::cout << GridLogMessage << "writing Sp(" << stNC << ") field" << std::endl; + ildgfmt.field = std::string("sp"+stNC+"gauge"); + } else { + static_assert(1, "Unrecognised group; unable to determine field tag"); + } - // populate 'rows' element of ildg header + // populate 'rows' element of ildg header if constexpr( matrix_fmt==MatrixFormat::REDUCED && is_su::value ) { ildgfmt.rows = Nc-1 ; } else if constexpr( matrix_fmt==MatrixFormat::REDUCED && is_sp::value ) { ildgfmt.rows = Nc/2 ; } else if constexpr( matrix_fmt==MatrixFormat::FULL ) { - ildgfmt.rows = Nc; - } else { - static_assert(1, "Unknown MatrixFormat specified"); - } + ildgfmt.rows = Nc; + } else { + static_assert(1, "Unknown MatrixFormat specified"); + } - // set fmt string for single/double precision - if constexpr( fp_fmt == FloatingPointFormat::IEEE32BIG ) { - header.floating_point = std::string("IEEE32BIG"); + // set fmt string for single/double precision + if constexpr( fp_fmt == FloatingPointFormat::IEEE32BIG ) { + header.floating_point = std::string("IEEE32BIG"); ildgfmt.precision = 32; - } else if constexpr ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { - header.floating_point = std::string("IEEE64BIG"); + } else if constexpr ( fp_fmt == FloatingPointFormat::IEEE64BIG ) { + header.floating_point = std::string("IEEE64BIG"); ildgfmt.precision = 64; - } else { - static_assert(1, "Unknown FloatingPointFormat specified"); - } + } else { + static_assert(1, "Unknown FloatingPointFormat specified"); + } ildgfmt.version = 1.2; ildgfmt.lx = header.dimension[0]; From 59ba80c7f13aca62df6c1728a45387b6a5568ccd Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 11 Mar 2026 17:20:55 +0000 Subject: [PATCH 47/67] minor formatting changes - got rid of inconsistent tabs that were leading to indentation issues --- Grid/parallelIO/MetaData.h | 29 ++++++++++++++--------------- 1 file changed, 14 insertions(+), 15 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index b75e090cab..bac9e8508b 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -576,18 +576,17 @@ struct GaugeUnMunger template struct GaugeUnMunger { - using in_type = typename vobj::scalar_object; + using in_type = typename vobj::scalar_object; using tmp_type = typename std::tuple_element_t(fp_fmt), std::tuple>; using out_type = typename std::tuple_element_t(fp_fmt), std::tuple>; - BinarySimpleUnmunger binary_unmunger; - //Gauge3x2unmunger gauge_unmunger; - GaugeSUunmunger gauge_unmunger; + BinarySimpleUnmunger binary_unmunger; + GaugeSUunmunger gauge_unmunger; - void operator() (in_type &in, out_type &out){ - tmp_type tmp; - binary_unmunger(in, tmp); - gauge_unmunger(tmp,out); + void operator() (in_type &in, out_type &out){ + tmp_type tmp; + binary_unmunger(in, tmp); + gauge_unmunger(tmp,out); } }; @@ -599,14 +598,14 @@ struct GaugeUnMunger using tmp_type = typename std::tuple_element_t(fp_fmt), std::tuple>; using out_type = typename std::tuple_element_t(fp_fmt), std::tuple>; - BinarySimpleUnmunger binary_unmunger; - GaugeSpunmunger gauge_unmunger; + BinarySimpleUnmunger binary_unmunger; + GaugeSpunmunger gauge_unmunger; - void operator() (in_type &in, out_type &out){ - tmp_type tmp; - binary_unmunger(in, tmp); - gauge_unmunger(tmp, out); - } + void operator() (in_type &in, out_type &out){ + tmp_type tmp; + binary_unmunger(in, tmp); + gauge_unmunger(tmp, out); + } }; From ad959fd469ad681391f8bb72283380f0a64b43ed Mon Sep 17 00:00:00 2001 From: gray95 Date: Wed, 11 Mar 2026 17:57:33 +0000 Subject: [PATCH 48/67] removes redundant ' == true' after std::is_same_v<> --- tests/IO/Test_ildg_reducedfmt_io.cc | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index ac9ed3a8fd..9fbd0f1121 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -49,9 +49,9 @@ LatticeGaugeField generateHotFieldConfiguration( GridCartesian &Grid, std::vecto LatticeGaugeField Umu(&Grid); - if constexpr( std::is_same_v == true ) { + if constexpr( std::is_same_v ) { SU::HotConfiguration(pRNGa,Umu); - } else if constexpr ( std::is_same_v == true ) { + } else if constexpr ( std::is_same_v ) { Sp::HotConfiguration(pRNGa,Umu); } else { static_assert(true, "Grid does not recognise the gauge group"); } From 9f02b0da3f2167242aa2687c0216242f2b64da0b Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Thu, 12 Mar 2026 13:19:15 +0000 Subject: [PATCH 49/67] changes unit test return type to void --- tests/IO/Test_io_mungers.cc | 136 ++++++++++++++---------------------- 1 file changed, 52 insertions(+), 84 deletions(-) diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 3b5d300ebe..d0104f13bd 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -33,17 +33,17 @@ using namespace Grid; // tests to check the various (un)munger classes // and functions in parallelIO/Metadata.h -bool check_reconstruct3(); -bool check_reconstructSU(); -bool check_is_perm_even(); -bool check_unique_reconstructSU(); -bool checkGauge3x2mungers(); -bool checkGaugeSUmungers(); -bool check_reconstructSp(); -bool checkGaugeSpmungers(); -bool checkBinarySimpleMungers(); -bool checkGaugeSimpleMungers(); -bool checkGaugeDoubleStoredMungers(); +void check_reconstruct3(); +void check_reconstructSU(); +void check_is_perm_even(); +void check_unique_reconstructSU(); +void checkGauge3x2mungers(); +void checkGaugeSUmungers(); +void check_reconstructSp(); +void checkGaugeSpmungers(); +void checkBinarySimpleMungers(); +void checkGaugeSimpleMungers(); +void checkGaugeDoubleStoredMungers(); void mock_SU_field(); void mock_Sp_field(); @@ -97,7 +97,7 @@ void mock_Sp_field( ColourMatrixD &cm, std::vector seed ) { } //////////////////////////////////////////// -bool check_reconstruct3() { +void check_reconstruct3() { ColourMatrixD SU3_xfield, SU3_yfield, SU3_zfield, SU3_tfield; @@ -131,11 +131,9 @@ bool check_reconstruct3() { assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-12 ); } - return true; - } -bool check_reconstructSU() { +void check_reconstructSU() { ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; @@ -169,11 +167,9 @@ bool check_reconstructSU() { assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-12 ); } - return true; - } -bool check_is_perm_even() { +void check_is_perm_even() { std::vector v = {0,1,2,3}; int count = 0; @@ -189,12 +185,10 @@ bool check_is_perm_even() { } } } - - return true; } -bool check_unique_reconstructSU() { +void check_unique_reconstructSU() { ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; @@ -228,11 +222,9 @@ bool check_unique_reconstructSU() { assert( norm2( new_cm*adj(new_cm)-1.0 ) < 1e-12 ); } - return true; - } -bool checkGauge3x2mungers() { +void checkGauge3x2mungers() { ColourMatrixD SU3_xfield, SU3_yfield, SU3_zfield, SU3_tfield; @@ -261,12 +253,10 @@ bool checkGauge3x2mungers() { // round-trip test assert(norm2(scalar_recon-scalar)<1e-10); - return true; - } template -bool checkGaugeSUmungers() { +void checkGaugeSUmungers() { ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; @@ -295,12 +285,10 @@ bool checkGaugeSUmungers() { // round-trip test assert(norm2(scalar_recon-scalar)<1e-12); - return true; - } -bool check_reconstructSp() { +void check_reconstructSp() { ColourMatrixD Sp_xfield, Sp_yfield, Sp_zfield, Sp_tfield; Sp_xfield = Sp_yfield = Sp_zfield = Sp_tfield = Zero(); @@ -332,11 +320,9 @@ bool check_reconstructSp() { } } - return true; - } -bool checkGaugeSpmungers() { +void checkGaugeSpmungers() { ColourMatrixD Sp_xfield, Sp_yfield, Sp_zfield, Sp_tfield; Sp_xfield = Sp_yfield = Sp_zfield = Sp_tfield = Zero(); @@ -368,12 +354,10 @@ bool checkGaugeSpmungers() { // round-trip test assert(scalar==scalar_recon); - return true; - } -bool checkBinarySimpleMungers() { +void checkBinarySimpleMungers() { LorentzColourMatrixF in_scalar_objectF; // single precision LorentzColourMatrixF out_scalar_objectF = Zero(); @@ -432,12 +416,10 @@ bool checkBinarySimpleMungers() { // lose exactness when going from double-->single precision assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); - return true; - } // these are used in NerscIO.h -bool checkGaugeSimpleMungers() { +void checkGaugeSimpleMungers() { LorentzColourMatrixF in_scalar_objectF; // single precision LorentzColourMatrixF out_scalar_objectF = Zero(); @@ -496,11 +478,9 @@ bool checkGaugeSimpleMungers() { // lose exactness when going from double-->single precision assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); - return true; - } -bool checkGaugeDoubleStoredMungers() { +void checkGaugeDoubleStoredMungers() { DoubleStoredColourMatrixF in_scalar_objectF; // single precision DoubleStoredColourMatrixF out_scalar_objectF = Zero(); @@ -559,8 +539,6 @@ bool checkGaugeDoubleStoredMungers() { // lose exactness when going from double-->single precision assert(norm2( in_scalar_objectD - out_scalar_objectD ) < 1e-10); - return true; - } @@ -571,62 +549,52 @@ int main (int argc, char ** argv) std::cout <1 && Nc<4) { - if( check_reconstruct3() ) { - std::cout << GridLogMessage << "reconstruct3: PASS" << std::endl; - } + check_reconstruct3(); + std::cout << GridLogMessage << "reconstruct3: PASS" << std::endl; - if( checkGauge3x2mungers() ) { - std::cout << GridLogMessage << "Gauge3x2mungers: PASS" << std::endl; - } - } - - if( check_reconstructSU() ) { - std::cout << GridLogMessage << "reconstructSU: PASS" << std::endl; - } - if( checkGaugeSUmungers() ) { - std::cout << GridLogMessage << "(unique_su=false) GaugeSUmungers: PASS" << std::endl; + checkGauge3x2mungers(); + std::cout << GridLogMessage << "Gauge3x2mungers: PASS" << std::endl; } - if( check_unique_reconstructSU() ) { - std::cout << GridLogMessage << "unique_reconstructSU: PASS" << std::endl; - } - if( checkGaugeSUmungers() ) { - std::cout << GridLogMessage << "(unique_su=true) GaugeSUmungers: PASS" << std::endl; - } + check_reconstructSU(); + std::cout << GridLogMessage << "reconstructSU: PASS" << std::endl; + + checkGaugeSUmungers(); + std::cout << GridLogMessage << "(unique_su=false) GaugeSUmungers: PASS" << std::endl; + + check_unique_reconstructSU(); + std::cout << GridLogMessage << "unique_reconstructSU: PASS" << std::endl; + + checkGaugeSUmungers(); + std::cout << GridLogMessage << "(unique_su=true) GaugeSUmungers: PASS" << std::endl; if constexpr(Nc>2 && Nc%2==0) { std::cout < Date: Thu, 12 Mar 2026 15:06:01 +0000 Subject: [PATCH 50/67] minor formatting issues - tabs v spaces. added comment to more clearly explain the element logic when reading ildg lattices. --- Grid/parallelIO/IldgIO.h | 74 +++++++++++++++++++++------------------- 1 file changed, 38 insertions(+), 36 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 7b4b844b14..50aa700971 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -829,8 +829,8 @@ class IldgReader : public GridLimeReader { // from its ildg-format header. if matrix_fmt==MatrixFormat::REDUCED // then Grid will reconstruct the full matrix using the appropriate munger. MatrixFormat matrix_fmt; - bool is_grp_su; - bool is_grp_sp; + bool is_grp_su = false; + bool is_grp_sp = false; // Binary format std::string format; FloatingPointFormat fp_fmt; @@ -856,7 +856,7 @@ class IldgReader : public GridLimeReader { // Copy out the string std::vector xmlc(nbytes+1,'\0'); limeReaderReadData((void *)&xmlc[0], &nbytes, LimeR); - // std::cout << GridLogMessage<< "Non binary record :" < element present = " << rows << ". So this lattice might be ildg 1.2 compliant." << std::endl; - } else { - std::cout << GridLogMessage << " element not present - adding it after ." << std::endl; - pugi::xml_node ildgfmt = doc.child("ildgFormat"); - const std::string stNC = std::to_string(Nc); - ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); - - // write result back into xmlstring - std::ostringstream ss; - doc.save(ss); - xmlstring = ss.str(); - } + // 'older' format ildg lattices don't have a element in their header. + // Grid's XML parsing doesn't support optional parameters with default values. + // Rather than rewrite that, here we explicitly add the element with default + // value = Nc to the ildg header if it is missing, before passing it to XmlReader." + pugi::xml_document doc; + doc.load_string(xmlstring.c_str()); + + if(doc.child("ildgFormat").child("rows")) { + std::string rows = doc.child("ildgFormat").child("rows").child_value(); + std::cout << GridLogMessage << " element present = " << rows << ". So this lattice might be ildg 1.2 compliant." << std::endl; + } else { + std::cout << GridLogMessage << " element not present - adding it after ." << std::endl; + pugi::xml_node ildgfmt = doc.child("ildgFormat"); + const std::string stNC = std::to_string(Nc); + ildgfmt.insert_child_after("rows", ildgfmt.child("field")).text().set(stNC.c_str()); + + // write result back into xmlstring + std::ostringstream ss; + doc.save(ss); + xmlstring = ss.str(); + } XmlReader RD(xmlstring, true, ""); read(RD,"ildgFormat",ildgFormat_); @@ -892,19 +895,19 @@ class IldgReader : public GridLimeReader { // set vars to tell Grid if it needs to reconstruct the full field. std::cout << GridLogMessage << "ildgFormat_.rows is " << ildgFormat_.rows << std::endl; matrix_fmt = (ildgFormat_.rows < Nc) ? MatrixFormat::REDUCED : MatrixFormat::FULL; - is_grp_su = (!strncmp(ildgFormat_.field.c_str(),"su",2)) ? true : false; - is_grp_sp = (!strncmp(ildgFormat_.field.c_str(),"sp",2)) ? true : false; + if( !strncmp(ildgFormat_.field.c_str(),"su",2) ) { is_grp_su = true; } + if( !strncmp(ildgFormat_.field.c_str(),"sp",2) ) { is_grp_sp = true; } - // check rows and Nc are related in the way we expect for each gauge group. - if (matrix_fmt==MatrixFormat::REDUCED && is_grp_su) { - assert( ildgFormat_.rows == Nc-1 ); - } - if (matrix_fmt==MatrixFormat::REDUCED && is_grp_sp) { - assert( ildgFormat_.rows == Nc/2 ); - } - if (matrix_fmt==MatrixFormat::FULL) { - assert( ildgFormat_.rows == Nc ); - } + // check rows and Nc are related in the way we expect for each gauge group. + if (matrix_fmt==MatrixFormat::REDUCED && is_grp_su) { + assert( ildgFormat_.rows == Nc-1 ); + } + if (matrix_fmt==MatrixFormat::REDUCED && is_grp_sp) { + assert( ildgFormat_.rows == Nc/2 ); + } + if (matrix_fmt==MatrixFormat::FULL) { + assert( ildgFormat_.rows == Nc ); + } assert( ildgFormat_.lx == dims[0]); assert( ildgFormat_.ly == dims[1]); @@ -963,11 +966,10 @@ class IldgReader : public GridLimeReader { // std::cout << GridLogMessage << "ILDG Binary record found : " ILDG_BINARY_DATA << std::endl; uint64_t offset= ftello(File); - - readLatticeBinaryObject(Umu, filename, fp_fmt, matrix_fmt, is_grp_su, is_grp_sp, offset, nersc_csum, scidac_csuma, scidac_csumb); - found_ildgBinary = 1; + readLatticeBinaryObject(Umu, filename, fp_fmt, matrix_fmt, is_grp_su, is_grp_sp, offset, nersc_csum, scidac_csuma, scidac_csumb); + found_ildgBinary = 1; } } From 22501e27d90f1e38f1c9e68cd888befab91c2c55 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Thu, 12 Mar 2026 17:20:41 +0000 Subject: [PATCH 51/67] edits comments and var names for clarity --- Grid/parallelIO/MetaData.h | 49 ++++++++++++++++++++------------------ 1 file changed, 26 insertions(+), 23 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index bac9e8508b..cb27cfd679 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -220,16 +220,16 @@ inline void reconstruct3(LorentzColourMatrix & cm) // the elements of the final row of an SU(N) matrix // can be obtained from the determinants of the N N-1 by N-1 matrices // formed by deleting each column and the Nth row of the SU(N) matrix. -// 1) create a ColourMatrixNm object to store the N-1 dim matrix. +// 1) create a ColourSubMatrix object to store the N-1 dim matrix. // 2) peek the colour matrix in a particular direction. -// 3) fill the ColourMatrixNm object from the peeked matrix. +// 3) fill the ColourSubMatrix object from the peeked matrix. // 4) take the Determinant and fill the corresponding element // 5) repeat for each Lorentz index inline void reconstructSU(LorentzColourMatrix &cm) { - using ColourMatrixNm = iScalar > > ; + using ColourSubMatrix = iScalar > > ; - ColourMatrixNm Su; // for the Nc-1 by Nc-1 matrix + ColourSubMatrix Su; // for the Nc-1 by Nc-1 matrix ColourMatrix SU; // for the Nc-1 by N matrix read from disk for(int mu=0; mu {1,0,2} --> {0,1,2} implies {1,2,0} is even. bool is_perm_even(std::vector &v) { int n = v.size(); @@ -274,42 +277,42 @@ bool is_perm_even(std::vector &v) { } } - -/////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////// // // this function follows the prescription laid out by the -// ildg spec (linked in the comment to reconstrucSp), -// forming a sum of products in lexicographic order. +// ildg spec, forming a sum of products in lexicographic order. // // SU[N-1][k] = sum ( SU[0][i0].SU[1][i1]...SU[N-2][iN-2] )* // {i0...iN-2!=k} // -/////////////////////////////////////////////////////////// +// see appendix A.2 of +// https://www-zeuthen.desy.de/apewww/ILDG/specifications/ildg-file-format-1.2.pdf +///////////////////////////////////////////////////////////// inline void unique_reconstructSU(LorentzColourMatrix &cm) { - std::vector v(Nc); - std::iota(v.begin(), v.end(), 0); // {0,1,...,Nc-1} + std::vector indices(Nc); + std::iota(indices.begin(), indices.end(), 0); // {0,1,...,Nc-1} for(int mu=0; mu c = v; - c.erase(c.begin()+k); // {0,...,k-1,k+1,...Nc-1} + std::vector cols = indices; + cols.erase(cols.begin()+k); // {0,...,k-1,k+1,...Nc-1} cm(mu)()(Nc-1,k) = 0; - // cycle through perms of c + // construct sum of products + // as we cycle through permutations of cols do { - std::vector P = c; - P.emplace_back(k); - bool E = is_perm_even(P); + std::vector perm = cols; + perm.emplace_back(k); ComplexD prod(1,0); for(int i=0;i struct GaugeSUmunger; // two specialisations for two ways to reconstruct -// NcxNc matrix from Nc-1xNc matrix +// NcxNc matrix from (Nc-1)xNc matrix template struct GaugeSUmunger{ From f5d882ce47be9b1c2256ff58545d2ad0f432145a Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Thu, 12 Mar 2026 17:43:34 +0000 Subject: [PATCH 52/67] add assert statement to ensure either su or sp field read --- Grid/parallelIO/IldgIO.h | 2 ++ 1 file changed, 2 insertions(+) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index 50aa700971..cbac86df1b 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -897,6 +897,8 @@ class IldgReader : public GridLimeReader { matrix_fmt = (ildgFormat_.rows < Nc) ? MatrixFormat::REDUCED : MatrixFormat::FULL; if( !strncmp(ildgFormat_.field.c_str(),"su",2) ) { is_grp_su = true; } if( !strncmp(ildgFormat_.field.c_str(),"sp",2) ) { is_grp_sp = true; } + // check if field element corresponds to either su or sp + assert( is_grp_su || is_grp_sp ); // check rows and Nc are related in the way we expect for each gauge group. if (matrix_fmt==MatrixFormat::REDUCED && is_grp_su) { From 8a0cb02cf3d485245c7e7a3d7b340e3981f281f2 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 17 Mar 2026 00:34:31 +0000 Subject: [PATCH 53/67] augment CheckpointerParameters and ILDGCheckpointer to enable checkpointing with reduced format su/sp lattices --- Grid/qcd/hmc/checkpointers/BaseCheckpointer.h | 13 ++- Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h | 49 +++++++- .../sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc | 105 ++++++++++++++++++ 3 files changed, 157 insertions(+), 10 deletions(-) create mode 100644 tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc diff --git a/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h b/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h index 72f5e39bab..65b32afe4b 100644 --- a/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h @@ -39,15 +39,20 @@ class CheckpointerParameters : Serializable { std::string, rng_prefix, int, saveInterval, bool, saveSmeared, - std::string, format, ); + std::string, format, + std::string, group, + bool, reduced_matrix, + bool, unique_su, ); - CheckpointerParameters(std::string cf = "cfg", std::string sf="cfg_smr" , std::string rn = "rng", - int savemodulo = 1, const std::string &f = "IEEE64BIG") + CheckpointerParameters(std::string cf = "cfg", std::string sf="cfg_smr" , std::string rn = "rng", int savemodulo = 1, const std::string& f = "IEEE64BIG", std::string grp = "su", bool rdc_mat = true, bool unq_su = false) : config_prefix(cf), smeared_prefix(sf), rng_prefix(rn), saveInterval(savemodulo), - format(f){}; + format(f), + group(grp), + reduced_matrix(rdc_mat), + unique_su(unq_su) {}; template diff --git a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h index 50cca63d17..7b61cb4e79 100644 --- a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h @@ -65,7 +65,10 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { << "Allowed: IEEE32BIG | IEEE32 | IEEE64BIG | IEEE64" << std::endl; - exit(1); + if (!(Params.group=="su" || Params.group=="sp")) { + std::cout << GridLogError << "Unrecognized file format " << Params.group << std::endl; + } + exit(1); } } @@ -86,10 +89,41 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { << scidac_csumb << std::dec << std::endl; - + IldgWriter _IldgWriter(grid->IsBoss()); _IldgWriter.open(config); - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + + if(Params.format=="IEEE64BIG") { + if(Params.group=="su" && Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else if (Params.group=="su" && !Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else if (Params.group=="sp" && Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else if (Params.group=="sp" && !Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else { assert(0); } + } + else if (Params.format=="IEEE32BIG") { + if(Params.group=="su" && Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else if (Params.group=="su" && !Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else if (Params.group=="sp" && Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else if (Params.group=="sp" && !Params.reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + } + else { assert(0); } + } + _IldgWriter.close(); std::cout << GridLogMessage << "Written ILDG Configuration on " << config @@ -123,15 +157,18 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { this->check_filename(rng); this->check_filename(config); - - uint32_t nersc_csum,scidac_csuma,scidac_csumb; BinaryIO::readRNG(sRNG, pRNG, rng, 0,nersc_csum,scidac_csuma,scidac_csumb); FieldMetaData header; IldgReader _IldgReader; _IldgReader.open(config); - _IldgReader.readConfiguration(U,header); // format from the header + + if(Params.unique_su) { + _IldgReader.readConfiguration(U,header);// format from the header + } + else { _IldgReader.readConfiguration(U,header); } + _IldgReader.close(); std::cout << GridLogMessage << "Read ILDG Configuration from " << config diff --git a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc new file mode 100644 index 0000000000..475960d7b5 --- /dev/null +++ b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc @@ -0,0 +1,105 @@ +/************************************************************************************* + +Grid physics library, www.github.com/paboyle/Grid + +Source file: ./tests/Test_hmc_WilsonFermionGauge.cc + +Copyright (C) 2015 + +Author: Peter Boyle +Author: neo +Author: Guido Cossu + +This program is free software; you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation; either version 2 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License along +with this program; if not, write to the Free Software Foundation, Inc., +51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +See the full license in the file "LICENSE" in the top level distribution +directory +*************************************************************************************/ +/* END LEGAL */ +#include + +int main(int argc, char **argv) +{ + using namespace Grid; + + Grid_init(&argc, &argv); + GridLogLayout(); + + typedef GenericSpHMCRunner HMCWrapper; + HMCWrapper TheHMC; + + TheHMC.Resources.AddFourDimGrid("gauge"); + + // Checkpointer definition + CheckpointerParameters CPparams; + CPparams.config_prefix = "ckpoint/lat"; + CPparams.rng_prefix = "ckpoint/rng"; + CPparams.smeared_prefix = "ckpoint/smr"; + CPparams.saveInterval = 5; + CPparams.format = "IEEE64BIG"; + + // Ildg specific parameters + CPparams.group = "sp"; + CPparams.reduced_matrix = true; + CPparams.unique_su = false; + + TheHMC.Resources.LoadILDGCheckpointer(CPparams); + + RNGModuleParameters RNGpar; + RNGpar.serial_seeds = "12 22 32 42 52"; + RNGpar.parallel_seeds = "76 77 87 79 70"; + TheHMC.Resources.SetRNGSeeds(RNGpar); + + // Construct observables + // here there is too much indirection + typedef PlaquetteMod PlaqObs; + typedef TopologicalChargeMod QObs; + TheHMC.Resources.AddObservable(); + TopologyObsParameters TopParams; + TopParams.interval = 5; + TopParams.do_smearing = true; + TopParams.Smearing.init_step_size = 0.01; + TopParams.Smearing.tolerance = 1e-5; + //TopParams.Smearing.steps = 200; + //TopParams.Smearing.step_size = 0.01; + TopParams.Smearing.meas_interval = 50; + TopParams.Smearing.maxTau = 2.0; + TheHMC.Resources.AddObservable(TopParams); + ////////////////////////////////////////////// + + ///////////////////////////////////////////////////////////// + // Collect actions, here use more encapsulation + // need wrappers of the fermionic classes + // that have a complex construction + // standard + RealD beta = 8.0 ; + SpWilsonGaugeActionR Waction(beta); + + ActionLevel Level1(1); + Level1.push_back(&Waction); + //Level1.push_back(WGMod.getPtr()); + TheHMC.TheAction.push_back(Level1); + ///////////////////////////////////////////////////////////// + + // HMC parameters are serialisable + TheHMC.Parameters.MD.MDsteps = 10; + TheHMC.Parameters.MD.trajL = 1.0; + + TheHMC.ReadCommandLine(argc, argv); // these can be parameters from file + TheHMC.Run(); // no smearing + + Grid_finalize(); + +} // main From 348990c7d4f4ec514d21fa4a143b670987cfd8b1 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 17 Mar 2026 01:31:49 +0000 Subject: [PATCH 54/67] streamline test of reduced sp2n hmc runner --- tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc index 475960d7b5..10355b4ef0 100644 --- a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc +++ b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc @@ -48,18 +48,18 @@ int main(int argc, char **argv) CPparams.rng_prefix = "ckpoint/rng"; CPparams.smeared_prefix = "ckpoint/smr"; CPparams.saveInterval = 5; - CPparams.format = "IEEE64BIG"; + CPparams.format = "IEEE32BIG"; + //CPparams.saveSmeared = false; // Ildg specific parameters CPparams.group = "sp"; CPparams.reduced_matrix = true; - CPparams.unique_su = false; TheHMC.Resources.LoadILDGCheckpointer(CPparams); RNGModuleParameters RNGpar; - RNGpar.serial_seeds = "12 22 32 42 52"; - RNGpar.parallel_seeds = "76 77 87 79 70"; + RNGpar.serial_seeds = "55 1 37 3 97 8"; + RNGpar.parallel_seeds = "4 8 15 16 23 42"; TheHMC.Resources.SetRNGSeeds(RNGpar); // Construct observables @@ -72,8 +72,6 @@ int main(int argc, char **argv) TopParams.do_smearing = true; TopParams.Smearing.init_step_size = 0.01; TopParams.Smearing.tolerance = 1e-5; - //TopParams.Smearing.steps = 200; - //TopParams.Smearing.step_size = 0.01; TopParams.Smearing.meas_interval = 50; TopParams.Smearing.maxTau = 2.0; TheHMC.Resources.AddObservable(TopParams); @@ -89,7 +87,6 @@ int main(int argc, char **argv) ActionLevel Level1(1); Level1.push_back(&Waction); - //Level1.push_back(WGMod.getPtr()); TheHMC.TheAction.push_back(Level1); ///////////////////////////////////////////////////////////// @@ -98,7 +95,7 @@ int main(int argc, char **argv) TheHMC.Parameters.MD.trajL = 1.0; TheHMC.ReadCommandLine(argc, argv); // these can be parameters from file - TheHMC.Run(); // no smearing + TheHMC.Run(); Grid_finalize(); From 246a070ad5cf659e3152d11ac03f38b0cd665912 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 17 Mar 2026 01:32:53 +0000 Subject: [PATCH 55/67] adds default false value for CheckpointParameters.saveSmeared --- Grid/qcd/hmc/checkpointers/BaseCheckpointer.h | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h b/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h index 65b32afe4b..03972a8561 100644 --- a/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/BaseCheckpointer.h @@ -44,11 +44,12 @@ class CheckpointerParameters : Serializable { bool, reduced_matrix, bool, unique_su, ); - CheckpointerParameters(std::string cf = "cfg", std::string sf="cfg_smr" , std::string rn = "rng", int savemodulo = 1, const std::string& f = "IEEE64BIG", std::string grp = "su", bool rdc_mat = true, bool unq_su = false) + CheckpointerParameters(std::string cf = "cfg", std::string sf="cfg_smr" , std::string rn = "rng", int savemodulo = 1, bool save_smr = false, const std::string& f = "IEEE64BIG", std::string grp = "su", bool rdc_mat = false, bool unq_su = false) : config_prefix(cf), smeared_prefix(sf), rng_prefix(rn), saveInterval(savemodulo), + saveSmeared(save_smr), format(f), group(grp), reduced_matrix(rdc_mat), From 6d34f5f092c3c579f527671006620d0ac6102be6 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 17 Mar 2026 01:34:48 +0000 Subject: [PATCH 56/67] fix check on gauge group from CheckpointerParameters.group --- Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h index 7b61cb4e79..a0817daabc 100644 --- a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h @@ -64,11 +64,12 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { std::cout << GridLogError << "Allowed: IEEE32BIG | IEEE32 | IEEE64BIG | IEEE64" << std::endl; - - if (!(Params.group=="su" || Params.group=="sp")) { - std::cout << GridLogError << "Unrecognized file format " << Params.group << std::endl; } - exit(1); + + if ( !((Params.group == std::string("su")) || (Params.group == std::string("sp"))) ) { + std::cout << GridLogError << "Unrecognized gauge group " << Params.group << std::endl; + std::cout << GridLogError << "Allowed: su | sp" << std::endl; + exit(1); } } From ec440dc8031d8a831e8a33e5c6468d8ecf213a82 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 17 Mar 2026 23:53:18 +0000 Subject: [PATCH 57/67] change order of template parameters in read/writeConfiguration for backwards compatibility. this necessitates explicitly specifying an additional template param in the ildg reduced format test. --- Grid/parallelIO/IldgIO.h | 6 +++--- tests/IO/Test_ildg_reducedfmt_io.cc | 6 ++++-- 2 files changed, 7 insertions(+), 5 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index cbac86df1b..ef22f2aa66 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -632,7 +632,7 @@ class IldgWriter : public ScidacWriter { // Don't require scidac records EXCEPT checksum // Use Grid MetaData object if present. //////////////////////////////////////////////////////////////// - template + template void writeConfiguration(Lattice &Umu, int sequence, std::string LFN, std::string description) { GridBase * grid = Umu.Grid(); @@ -797,7 +797,7 @@ class IldgReader : public GridLimeReader { // Else use ILDG MetaData object if present. // Else use SciDAC MetaData object if present. //////////////////////////////////////////////////////////////// - template + template void readConfiguration(Lattice &Umu, FieldMetaData &FieldMetaData_) { GridBase *grid = Umu.Grid(); @@ -828,9 +828,9 @@ class IldgReader : public GridLimeReader { // these variables store information about the lattice that is read // from its ildg-format header. if matrix_fmt==MatrixFormat::REDUCED // then Grid will reconstruct the full matrix using the appropriate munger. - MatrixFormat matrix_fmt; bool is_grp_su = false; bool is_grp_sp = false; + MatrixFormat matrix_fmt; // Binary format std::string format; FloatingPointFormat fp_fmt; diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index 9fbd0f1121..50783055d6 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -76,12 +76,14 @@ void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, st return; } + using stats = PeriodicGaugeStatistics; + std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); LatticeGaugeField Umu_saved(&Grid); @@ -93,7 +95,7 @@ void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, st std::cout <(Umu_saved, header); + _IldgReader.readConfiguration(Umu_saved, header); _IldgReader.close(); std::cout <(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); + _IldgWriter.writeConfiguration(Umu,4000,std::string("dummy_ildg_LFN"),std::string("dummy_config")); _IldgWriter.close(); LatticeGaugeField Umu_saved(&Grid); @@ -93,7 +95,7 @@ void writeReadIldgConfiguration( LatticeGaugeField &Umu, GridCartesian &Grid, st std::cout <(Umu_saved, header); + _IldgReader.readConfiguration(Umu_saved, header); _IldgReader.close(); std::cout < 0){ + // If Integer is greater than the one it's pointing to and + if(elements[i] > elements[i-1] && elements[i] > elements[mobileIndex]) { + isMobile = true; // greater than the previously found mobile integer + mobileIndex = i; + } + } + else if(dirs[i] == ">" && i < Nc-1) { + if(elements[i] > elements[i+1] && elements[i] > elements[mobileIndex]) { + isMobile = true; + mobileIndex = i; + } } } - } - + // If there is a valid mobile integer swap with the one it's pointing to + if(isMobile) { + if(dirs[mobileIndex] == "<") { // Integers pointing left + std::iter_swap(elements.begin()+mobileIndex-1, elements.begin()+mobileIndex); + std::iter_swap(dirs.begin()+mobileIndex-1, dirs.begin()+mobileIndex); + // Change the direction of every integer greater than the mobile integer + for(int i = 0; i < Nc; i++) { + if(elements[i] > elements[mobileIndex - 1]){ + if(dirs[i] == "<") { dirs[i] = ">"; } + else { dirs[i] = "<"; } + } + } + } + else { // Integers pointing right + std::iter_swap(elements.begin()+mobileIndex, elements.begin()+mobileIndex+1); + std::iter_swap(dirs.begin()+mobileIndex, dirs.begin()+mobileIndex+1); + // Change the direction of every integer greater than the mobile integer + for(int i = 0; i < Nc; i++) { + if(elements[i] > elements[mobileIndex + 1]){ + if(dirs[i] == "<") { dirs[i] = ">"; } + else { dirs[i] = "<"; } + } + } + } + // test if permutation is even or odd + if(count%2==0) { assert(is_perm_even(elements) == true); } + else { assert(is_perm_even(elements) == false); } + count++; + } + } } + void check_unique_reconstructSU() { ColourMatrixD SUN_xfield, SUN_yfield, SUN_zfield, SUN_tfield; From 32b9cc302753ba7ddcc0cca9b2a89d51364996fc Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 24 Mar 2026 17:26:20 +0000 Subject: [PATCH 60/67] factor out logic for choosing IldgWriter template instantiation and remove little endian from initialisation check --- Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h | 107 ++++++++++-------- 1 file changed, 57 insertions(+), 50 deletions(-) diff --git a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h index a0817daabc..e2dd98e959 100644 --- a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h @@ -54,16 +54,15 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { // check here that the format is valid int ieee32big = (Params.format == std::string("IEEE32BIG")); - int ieee32 = (Params.format == std::string("IEEE32")); int ieee64big = (Params.format == std::string("IEEE64BIG")); - int ieee64 = (Params.format == std::string("IEEE64")); - if (!(ieee64big || ieee32 || ieee32big || ieee64)) { + if (!(ieee64big || ieee32big)) { std::cout << GridLogError << "Unrecognized file format " << Params.format << std::endl; std::cout << GridLogError - << "Allowed: IEEE32BIG | IEEE32 | IEEE64BIG | IEEE64" + << "Allowed: IEEE32BIG | IEEE64BIG" << std::endl; + exit(1); } if ( !((Params.group == std::string("su")) || (Params.group == std::string("sp"))) ) { @@ -73,79 +72,87 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { } } - void TrajectoryComplete(int traj, - ConfigurationBase &SmartConfig, - GridSerialRNG &sRNG, - GridParallelRNG &pRNG) { - if ((traj % Params.saveInterval) == 0) { - std::string config, rng, smr; - this->build_filenames(traj, Params, config, smr, rng); + // choose appropriate template instantiation here since the + // non-const checkpointer parameters cannot be used directly as + // template arguments to IldgWriter + void chooseIldgWriter( std::string format, std::string group, bool reduced_matrix, + std::string lat_obj, int traj, + ConfigurationBase &SmartConfig) { + GridBase *grid = SmartConfig.get_U(false).Grid(); - uint32_t nersc_csum,scidac_csuma,scidac_csumb; - BinaryIO::writeRNG(sRNG, pRNG, rng, 0,nersc_csum,scidac_csuma,scidac_csumb); - std::cout << GridLogMessage << "Written BINARY RNG " << rng - << " checksum " << std::hex - << nersc_csum<<"/" - << scidac_csuma<<"/" - << scidac_csumb - << std::dec << std::endl; - IldgWriter _IldgWriter(grid->IsBoss()); - _IldgWriter.open(config); + _IldgWriter.open(lat_obj); - if(Params.format=="IEEE64BIG") { - if(Params.group=="su" && Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + if(format=="IEEE64BIG") { + if(group=="su" && reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else if (Params.group=="su" && !Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (group=="su" && !reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else if (Params.group=="sp" && Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (group=="sp" && reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else if (Params.group=="sp" && !Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (group=="sp" && !reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else { assert(0); } } - else if (Params.format=="IEEE32BIG") { - if(Params.group=="su" && Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (format=="IEEE32BIG") { + if(group=="su" && reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else if (Params.group=="su" && !Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (group=="su" && !reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else if (Params.group=="sp" && Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (group=="sp" && reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else if (Params.group=="sp" && !Params.reduced_matrix) { - _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, config, config); + else if (group=="sp" && !reduced_matrix) { + _IldgWriter.writeConfiguration(SmartConfig.get_U(false), traj, lat_obj, lat_obj); } - else { assert(0); } } _IldgWriter.close(); + } + - std::cout << GridLogMessage << "Written ILDG Configuration on " << config + void TrajectoryComplete(int traj, + ConfigurationBase &SmartConfig, + GridSerialRNG &sRNG, + GridParallelRNG &pRNG) { + if ((traj % Params.saveInterval) == 0) { + std::string config, rng, smr; + this->build_filenames(traj, Params, config, smr, rng); + uint32_t nersc_csum,scidac_csuma,scidac_csumb; + BinaryIO::writeRNG(sRNG, pRNG, rng, 0,nersc_csum,scidac_csuma,scidac_csumb); + std::cout << GridLogMessage << "Written BINARY RNG " << rng << " checksum " << std::hex << nersc_csum<<"/" << scidac_csuma<<"/" << scidac_csumb << std::dec << std::endl; - if ( Params.saveSmeared ) { - IldgWriter _IldgWriter(grid->IsBoss()); - _IldgWriter.open(smr); - _IldgWriter.writeConfiguration(SmartConfig.get_U(true), traj, smr, smr); - _IldgWriter.close(); + chooseIldgWriter( Params.format, Params.group, Params.reduced_matrix, config, traj, + SmartConfig ); - std::cout << GridLogMessage << "Written ILDG Configuration on " << smr + std::cout << GridLogMessage << "Written ILDG Configuration on " << config << " checksum " << std::hex << nersc_csum<<"/" << scidac_csuma<<"/" << scidac_csumb << std::dec << std::endl; + + if ( Params.saveSmeared ) { + chooseIldgWriter( Params.format, Params.group, Params.reduced_matrix, smr, traj, + SmartConfig ); + + std::cout << GridLogMessage << "Written ILDG Configuration on " << smr + << " checksum " << std::hex + << nersc_csum<<"/" + << scidac_csuma<<"/" + << scidac_csumb + << std::dec << std::endl; } } @@ -166,9 +173,9 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { _IldgReader.open(config); if(Params.unique_su) { - _IldgReader.readConfiguration(U,header);// format from the header + _IldgReader.readConfiguration(U,header);// format from the header } - else { _IldgReader.readConfiguration(U,header); } + else { _IldgReader.readConfiguration(U,header); } _IldgReader.close(); From aa25bdf45bac159aff0d2d708f1a4c824d6c93f6 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 24 Mar 2026 17:28:54 +0000 Subject: [PATCH 61/67] re-order for clarity and show unique_su parameter is available for su checkpointing --- tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc index 10355b4ef0..c4a8230af5 100644 --- a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc +++ b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc @@ -49,11 +49,12 @@ int main(int argc, char **argv) CPparams.smeared_prefix = "ckpoint/smr"; CPparams.saveInterval = 5; CPparams.format = "IEEE32BIG"; - //CPparams.saveSmeared = false; + CPparams.saveSmeared = true; // Ildg specific parameters CPparams.group = "sp"; CPparams.reduced_matrix = true; + //CPparams.unique_su = true; // only applies to su fields TheHMC.Resources.LoadILDGCheckpointer(CPparams); @@ -65,8 +66,9 @@ int main(int argc, char **argv) // Construct observables // here there is too much indirection typedef PlaquetteMod PlaqObs; - typedef TopologicalChargeMod QObs; TheHMC.Resources.AddObservable(); + + typedef TopologicalChargeMod QObs; TopologyObsParameters TopParams; TopParams.interval = 5; TopParams.do_smearing = true; From 1107e6033eacef34e86b0f1d13338de1175af6a4 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Tue, 24 Mar 2026 21:13:09 +0000 Subject: [PATCH 62/67] add runtime check for even Nc when writing Sp fields --- Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h index e2dd98e959..29d7484af9 100644 --- a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h @@ -70,6 +70,13 @@ class ILDGHmcCheckpointer : public BaseHmcCheckpointer { std::cout << GridLogError << "Allowed: su | sp" << std::endl; exit(1); } + + if ( Params.group == std::string("sp") && Nc%2!=0 ) { + std::cout << GridLogError << "Nc=" << Nc; + std::cout << ", Sp fields require even Nc" << std::endl; + exit(1); + } + } // choose appropriate template instantiation here since the From a7b1e476b680693a8162ddec5c16fe23a9543fe6 Mon Sep 17 00:00:00 2001 From: gray95 Date: Tue, 31 Mar 2026 18:09:46 +0100 Subject: [PATCH 63/67] adds static_assert to catch users trying to write Sp fields with Nc odd --- Grid/parallelIO/IldgIO.h | 3 +++ 1 file changed, 3 insertions(+) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index ef22f2aa66..f150ff555f 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -636,6 +636,9 @@ class IldgWriter : public ScidacWriter { void writeConfiguration(Lattice &Umu, int sequence, std::string LFN, std::string description) { GridBase * grid = Umu.Grid(); + + // catch malformed (wrt group_name) template instantiations at compile-time + static_assert( std::is_same_v || (std::is_same_v && Nc%2==0), "IldgWriter supports SU(Nc) and Sp(Nc=2k) lattices. For Sp fields Nc must be even" ); //////////////////////////////////////// // fill the Grid header From 1a3e3a6aabb89263387d480f9912e9fd251fd24a Mon Sep 17 00:00:00 2001 From: Ed Bennett Date: Wed, 1 Apr 2026 16:04:43 +0100 Subject: [PATCH 64/67] revert change outside of target region --- Grid/parallelIO/MetaData.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index cb27cfd679..1f7d210bc3 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -363,7 +363,7 @@ struct BinarySimpleUnmunger { typedef typename getPrecision::real_scalar_type sobj_stype; void operator()(sobj &in, fobj &out) { - // take word by word and transform according to the status + // take word by word and transform accoding to the status fobj_stype *out_buffer = (fobj_stype *)&out; sobj_stype *in_buffer = (sobj_stype *)∈ size_t fobj_words = sizeof(out) / sizeof(fobj_stype); From e3540885adc48304f447466c71be830af90844e0 Mon Sep 17 00:00:00 2001 From: Ed Bennett Date: Wed, 1 Apr 2026 16:04:55 +0100 Subject: [PATCH 65/67] revert change outside target region --- Grid/parallelIO/MetaData.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 1f7d210bc3..f1c94b5aeb 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -401,8 +401,8 @@ struct GaugeSimpleMunger{ void operator()(fobj &in, sobj &out) { for (int mu = 0; mu < Nd; mu++) { for (int i = 0; i < Nc; i++) { -for (int j = 0; j < Nc; j++) { - out(mu)()(i, j) = in(mu)()(i, j); + for (int j = 0; j < Nc; j++) { + out(mu)()(i, j) = in(mu)()(i, j); }} } }; From a5905857b3d2cf70dba618b1cbb3e76eddf5d302 Mon Sep 17 00:00:00 2001 From: Ed Bennett Date: Wed, 1 Apr 2026 16:05:03 +0100 Subject: [PATCH 66/67] update copyright messages --- Grid/parallelIO/IldgIO.h | 6 +++++- Grid/parallelIO/MetaData.h | 3 ++- Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h | 3 ++- tests/IO/Test_IldgWriter.cc | 2 +- tests/IO/Test_ildg_reducedfmt_io.cc | 2 +- tests/IO/Test_io_mungers.cc | 2 +- tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc | 3 ++- 7 files changed, 14 insertions(+), 7 deletions(-) diff --git a/Grid/parallelIO/IldgIO.h b/Grid/parallelIO/IldgIO.h index f150ff555f..fc47fc4e52 100644 --- a/Grid/parallelIO/IldgIO.h +++ b/Grid/parallelIO/IldgIO.h @@ -4,7 +4,11 @@ Grid physics library, www.github.com/paboyle/Grid Source file: ./lib/parallelIO/IldgIO.h -Copyright (C) 2015 +Copyright (C) 2015, 2026 + +Author: Peter Boyle +Author: Guido Cossu +Author: Gaurav Ray This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index f1c94b5aeb..34cdfc3587 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -4,10 +4,11 @@ Source file: ./lib/parallelIO/Metadata.h - Copyright (C) 2015 + Copyright (C) 2015, 2026 Author: Peter Boyle Author: Jamie Hudspith + Author: Gaurav Ray This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by diff --git a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h index 29d7484af9..80e413ba11 100644 --- a/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h +++ b/Grid/qcd/hmc/checkpointers/ILDGCheckpointer.h @@ -4,9 +4,10 @@ Grid physics library, www.github.com/paboyle/Grid Source file: ./lib/qcd/hmc/ILDGCheckpointer.h -Copyright (C) 2016 +Copyright (C) 2016, 2026 Author: Guido Cossu +Author: Gaurav Ray This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by diff --git a/tests/IO/Test_IldgWriter.cc b/tests/IO/Test_IldgWriter.cc index 3e0858e6bf..003bb79717 100644 --- a/tests/IO/Test_IldgWriter.cc +++ b/tests/IO/Test_IldgWriter.cc @@ -4,7 +4,7 @@ Source file: ./Test_IldgWriter.cc - Copyright (C) 2015 + Copyright (C) 2026 Author: Azusa Yamaguchi Author: Peter Boyle diff --git a/tests/IO/Test_ildg_reducedfmt_io.cc b/tests/IO/Test_ildg_reducedfmt_io.cc index 50783055d6..138e6990d9 100644 --- a/tests/IO/Test_ildg_reducedfmt_io.cc +++ b/tests/IO/Test_ildg_reducedfmt_io.cc @@ -4,7 +4,7 @@ Source file: ./Test_ildg_reducedfmt_io.cc - Copyright (C) 2015 + Copyright (C) 2026 Author: Azusa Yamaguchi Author: Peter Boyle diff --git a/tests/IO/Test_io_mungers.cc b/tests/IO/Test_io_mungers.cc index 3d2873b991..0e13946ded 100644 --- a/tests/IO/Test_io_mungers.cc +++ b/tests/IO/Test_io_mungers.cc @@ -4,7 +4,7 @@ Source file: ./Test_io_mungers.cc - Copyright (C) 2015 + Copyright (C) 2026 Author: Azusa Yamaguchi Author: Peter Boyle diff --git a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc index c4a8230af5..cc1de6011c 100644 --- a/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc +++ b/tests/sp2n/Test_hmc_Sp_ildg_pureGaugeWilson.cc @@ -4,11 +4,12 @@ Grid physics library, www.github.com/paboyle/Grid Source file: ./tests/Test_hmc_WilsonFermionGauge.cc -Copyright (C) 2015 +Copyright (C) 2026 Author: Peter Boyle Author: neo Author: Guido Cossu +Author: Gaurav Ray This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by From e608ff5065425de31e50516ada2f4fd703b550b8 Mon Sep 17 00:00:00 2001 From: Gaurav Ray Date: Wed, 1 Apr 2026 21:15:45 +0100 Subject: [PATCH 67/67] mark is_perm_even as inline to avoid multiple re-definition errors --- Grid/parallelIO/MetaData.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Grid/parallelIO/MetaData.h b/Grid/parallelIO/MetaData.h index 34cdfc3587..4255f024ca 100644 --- a/Grid/parallelIO/MetaData.h +++ b/Grid/parallelIO/MetaData.h @@ -253,7 +253,7 @@ inline void reconstructSU(LorentzColourMatrix &cm) // is an even or odd parity permutation. even/odd if the number of // inversions needed to get back to the lexicographic 1st sequence is // even or odd. ex: {1,2,0} --> {1,0,2} --> {0,1,2} implies {1,2,0} is even. -bool is_perm_even(std::vector &v) { +inline bool is_perm_even(std::vector &v) { int n = v.size(); std::vector a(n,0);