From 41f5a022047049ce61d4761897810bf8620470f5 Mon Sep 17 00:00:00 2001 From: Peter Boyle Date: Tue, 11 Aug 2026 13:33:04 -0400 Subject: [PATCH] IO updates for AllToAllV aggregation of contiguous chunks. May finally clean up the poor MPI2 IO performance issue that has been persistent. --- CLAUDE.md | 76 +++- Grid/communicator/Communicator_base.h | 10 + Grid/communicator/Communicator_mpi3.cc | 20 ++ Grid/communicator/Communicator_none.cc | 11 + Grid/parallelIO/BinaryIO.cc | 4 + Grid/parallelIO/BinaryIO.h | 457 +++++++++++++++++++++++- tests/IO/Test_aggregate_io.cc | 458 +++++++++++++++++++++++++ 7 files changed, 1019 insertions(+), 17 deletions(-) create mode 100644 tests/IO/Test_aggregate_io.cc diff --git a/CLAUDE.md b/CLAUDE.md index 4b6acf2aa..11f221c6e 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -39,17 +39,39 @@ Platform recipes from `README.md`: Required external libs: GMP, MPFR, OpenSSL, zlib. -## Running Tests +### Use `systems/` for real machines + +`systems//` holds the known-good build for each production platform (`Frontier`, `Aurora`, `Perlmutter`, `Summit`, `Tursa`, `Lumi`, `Booster`, `Crusher`, `SDCC-*`, `mac-arm`, …). Each contains a `config-command` (the exact `../../configure` invocation) and a `sourceme.sh` (module loads and env). **Prefer copying/adapting these over hand-rolling configure flags** — they encode compiler workarounds, `LDFLAGS`, and shared-memory settings that are easy to get wrong. `systems/WorkArounds.txt` records known vendor bugs. + +Note the GPU builds use `--enable-simd=GPU --enable-gen-simd-width=64`, so `Nsimd` is *not* 1 on device (it is `64/sizeof(scalar)`). + +### Regenerating `Make.inc` — required after adding or deleting source files + +`Make.inc` files are generated, not tracked in git (`.gitignore`d). `scripts/filelist` walks `Grid/`, `tests/*`, `benchmarks/`, `examples/`, and `HMC/` and writes the file lists and per-test `bin_PROGRAMS` rules. Every new `.cc`/`.h` in `Grid/`, and every new `Test_*.cc` / `Benchmark_*.cc` / `Example_*.cc`, is invisible to the build until you run: + +```bash +./scripts/filelist # from the source root, then re-run configure/make +``` + +`bootstrap.sh` runs it for you on the first setup. + +## Running Tests and Benchmarks ```bash # From build directory -make check # root-level tests (Test_simd, Test_cshift, etc.) -make -C tests/ tests # build tests in a subdirectory +make check # runs only Test_simd, Test_cshift, Test_stencil, Test_dwf_mixedcg_prec +make -C tests/ tests # build (not run) the tests in a subdirectory +make bench # build benchmarks ./tests/core/Test_simd # run a single test binary directly +mpirun -n 4 ./tests/core/Test_cshift --grid 16.16.16.16 --mpi 1.1.1.4 ``` +`make check` is a thin smoke test — building a subdirectory with `make -C tests/ tests` and running the relevant binaries directly is the normal development loop. Test binaries take Grid's standard command-line arguments (`--grid`, `--mpi`, `--accelerator-threads`, `--threads`, `--debug-signals`, `--log`); see `Grid/util/Init.cc`. + Test subdirectories and their focus: `core` (SIMD, stencil, comms), `solver` (CG, GMRES, eigensolvers), `hmc` (MD integrators), `forces` (fermion forces), `lanczos`, `IO`, `smearing`, `sp2n`, `debug`. +Tests and benchmarks that need optional fermion representations are guarded by `disable_tests_without_instantiations.h` / `disable_benchmarks_without_instantiations.h`, so a `--disable-fermion-reps --disable-gparity` build silently compiles them to no-ops. + ## Architecture ### Layer stack (bottom to top) @@ -75,9 +97,35 @@ Test subdirectories and their focus: `core` (SIMD, stencil, comms), `solver` (CG - `smearing/` — APE, Stout, HEX, gradient flow - `observables/` — Polyakov loop, plaquette, topological charge -### GPU acceleration +### GPU acceleration and the view/memory-manager discipline -GPU support is injected via macros (`accelerator_for`, `accelerator_for2dNB`). The `Grid/simd/` SIMD types map to scalar on GPU device code; host code paths remain vectorised. Unified virtual memory is on by default (`--enable-unified=yes`); device-aware MPI (`--enable-accelerator-aware-mpi`) avoids device→host copies on transfers. +GPU support is injected via macros in `Grid/threads/Accelerator.h` — `accelerator_for(i, n, nsimd, {...})`, `accelerator_forNB` (non-blocking, must be followed by `accelerator_barrier()`), `accelerator_for2dNB`, and `accelerator_inline`. On a CPU build these degrade to `thread_for` (OpenMP). Unified virtual memory is on by default (`--enable-unified=yes`); device-aware MPI (`--enable-accelerator-aware-mpi`) avoids device→host copies on transfers. + +Lattice data is **not** directly addressable inside a kernel. You must open a view with the correct access mode so `Grid/allocator/MemoryManager.h` can move/mark the data: + +```cpp +autoView(out_v, out, AcceleratorWriteDiscard); // RAII; closes at end of scope +autoView(in_v, in, AcceleratorRead); +accelerator_for(ss, grid->oSites(), Nsimd, { + coalescedWrite(out_v[ss], coalescedRead(in_v[ss])); +}); +``` + +Modes are `AcceleratorRead/Write/WriteDiscard` and `CpuRead/Write/WriteDiscard`. Getting the mode wrong (e.g. `AcceleratorRead` on a field you write) produces stale-data bugs that only appear on GPU builds. Inside kernels use `coalescedRead`/`coalescedWrite` rather than raw `operator[]` — they map the SIMD lane onto `threadIdx.x` so accesses stay coalesced. + +### Repo-local debugging skills (`skills/`) + +`skills/` contains hard-won, Grid-specific playbooks written as invocable skill files. Consult them before debugging in these areas rather than reasoning from first principles: + +| File | Covers | +|---|---| +| `gpu-memory-performance.md` | `acceleratorThreads()`, LambdaApply thread mapping, `coalescedRead` idiom, fused vs staged HBM access | +| `gpu-runtime-correctness.md` | GPU runtime returning early from sync, silent wrong answers | +| `communication-overlap.md` | 7-phase halo pipeline, per-packet events, host-staging vs GPU-direct RDMA | +| `mpi-heterogeneous.md` | `MPI_Sendrecv` device-buffer aliasing, deterministic reductions | +| `compiler-validation.md` | Isolating GPU compiler codegen bugs, minimal reproducers | +| `correctness-verification.md` | Double-run fingerprinting, per-packet checksums, flight recorder | +| `hang-diagnosis.md` | Diagnosing MPI/accelerator hangs | ### Memory and I/O @@ -85,14 +133,24 @@ GPU support is injected via macros (`accelerator_for`, `accelerator_for2dNB`). T - `Grid/parallelIO/` — distributed parallel reader/writer for ILDG (via LIME), SciDAC, and native binary formats - `Grid/serialisation/` — text, binary, HDF5, XML/JSON serialisation of arbitrary Grid objects -### HMC applications +### Executables -`HMC/` contains production-ready HMC driver programmes (e.g. `Mobius2p1f.cc`, `DWF_plus_DSDR_nf2plus1_Shamir_Gparity.cc`). These are built separately from the library tests. +- `HMC/` — production HMC driver programmes (e.g. `Mobius2p1f.cc`, `DWF_plus_DSDR_nf2plus1_Shamir_Gparity.cc`) +- `benchmarks/` — `Benchmark_dwf`, `Benchmark_ITT`, `Benchmark_comms`, `Benchmark_memory_bandwidth`, … used to qualify a new machine +- `examples/` — small, readable programmes (`Example_plaquette.cc`, `Example_Mobius_spectrum.cc`) that are the best starting point for learning the API + +Each of these directories auto-builds every top-level `.cc` as its own binary via `scripts/filelist`. + +Every programme is wrapped in `Grid_init(&argc, &argv)` / `Grid_finalize()` (`Grid/util/Init.h`). ## Key Conventions - **C++17** is required throughout. -- Template structure: most classes are templated on `<_FImpl>` (fermion impl) or `` (gauge impl), which encode the representation and precision. Instantiation is controlled by `--enable-fermion-instantiations`. +- Template structure: most classes are templated on `<_FImpl>` (fermion impl) or `` (gauge impl), which encode the representation and precision. Instantiation is controlled by `--enable-fermion-instantiations`; the explicit instantiation `.cc` files live under `Grid/qcd/action/fermion/instantiation//` and are selected by `scripts/filelist`. - The `RealD`/`RealF`/`ComplexD`/`ComplexF` typedefs are used everywhere; avoid raw `double`/`float`. -- Logging uses `Grid_log`, `Grid_error` macros (from `Grid/log/`); performance-critical paths use the `GRID_TRACE` / timer macros from `Grid/perfmon/`. +- Use `GRID_ASSERT(cond)` (defined in `Grid/GridStd.h`), not bare `assert` — it prints a Grid-formatted message and aborts cleanly under MPI. +- Logging is stream-based, not macro-based: `std::cout << GridLogMessage << ... << std::endl;`. Channels declared in `Grid/log/Log.h` include `GridLogError`, `GridLogWarning`, `GridLogDebug`, `GridLogPerformance`, `GridLogIterative`, `GridLogSolver`, `GridLogHMC`, `GridLogComms`, `GridLogMemory`, `GridLogDslash`, `GridLogIRL`, `GridLogMG`. A subset is switched on at runtime with e.g. `--log Error,Warning,Message,Performance,Iterative,Integrator,Debug,Colours` (names given without the `GridLog` prefix). +- Performance-critical paths use `GRID_TRACE(name)` from `Grid/perfmon/Tracing.h` (compiled out unless `--enable-tracing` selects a backend) and the `GridStopWatch` timers in `Grid/perfmon/Timer.h`. - Reductions across MPI ranks go through `GridBase::GlobalSum` / `GlobalMax`; never reduce with bare MPI calls inside library code. +- Everything lives in `NAMESPACE_BEGIN(Grid)` / `NAMESPACE_END(Grid)` macros; follow the surrounding file rather than writing `namespace Grid { }`. +- British spelling is used in identifiers and comments (`colour`, `neighbour`, `serialisation`). diff --git a/Grid/communicator/Communicator_base.h b/Grid/communicator/Communicator_base.h index deb93faee..8c235c289 100644 --- a/Grid/communicator/Communicator_base.h +++ b/Grid/communicator/Communicator_base.h @@ -238,6 +238,16 @@ public: } void AllToAll(int dim ,void *in,void *out,uint64_t words,uint64_t bytes); void AllToAll(void *in,void *out,uint64_t words ,uint64_t bytes); + + //////////////////////////////////////////////////////////////////////////// + // Variable count all to all. Counts and displacements are in units of + // "bytes" sized words and are indexed by rank within this communicator. + // For exchanges that are a permutation but do not divide evenly between + // ranks; AllToAll above is the uniform count special case. + //////////////////////////////////////////////////////////////////////////// + void AllToAllV(void *in ,const std::vector &sendcounts,const std::vector &senddispls, + void *out,const std::vector &recvcounts,const std::vector &recvdispls, + uint64_t bytes); template void Broadcast(int root,obj &data) { diff --git a/Grid/communicator/Communicator_mpi3.cc b/Grid/communicator/Communicator_mpi3.cc index 3b8561d3c..d48c2aa9b 100644 --- a/Grid/communicator/Communicator_mpi3.cc +++ b/Grid/communicator/Communicator_mpi3.cc @@ -945,5 +945,25 @@ void CartesianCommunicator::AllToAll(void *in,void *out,uint64_t words,uint64_t MPI_Alltoall(in,iwords,object,out,iwords,object,communicator); MPI_Type_free(&object); } +void CartesianCommunicator::AllToAllV(void *in ,const std::vector &sendcounts,const std::vector &senddispls, + void *out,const std::vector &recvcounts,const std::vector &recvdispls, + uint64_t bytes) +{ + FlightRecorder::StepLog("AllToAllV"); + GRID_ASSERT(sendcounts.size()==(size_t)_Nprocessors); + GRID_ASSERT(senddispls.size()==(size_t)_Nprocessors); + GRID_ASSERT(recvcounts.size()==(size_t)_Nprocessors); + GRID_ASSERT(recvdispls.size()==(size_t)_Nprocessors); + // MPI counts are "int"; the caller sizes the word to keep them in range + int ibytes = bytes; + GRID_ASSERT(bytes == (uint64_t)ibytes); + MPI_Datatype object; + MPI_Type_contiguous(ibytes,MPI_BYTE,&object); + MPI_Type_commit(&object); + int ierr = MPI_Alltoallv(in ,(int *)&sendcounts[0],(int *)&senddispls[0],object, + out,(int *)&recvcounts[0],(int *)&recvdispls[0],object,communicator); + GRID_ASSERT(ierr==0); + MPI_Type_free(&object); +} NAMESPACE_END(Grid); diff --git a/Grid/communicator/Communicator_none.cc b/Grid/communicator/Communicator_none.cc index a7232fcdd..6e56bc211 100644 --- a/Grid/communicator/Communicator_none.cc +++ b/Grid/communicator/Communicator_none.cc @@ -113,6 +113,17 @@ void CartesianCommunicator::AllToAll(void *in,void *out,uint64_t words,uint64_t { bcopy(in,out,bytes*words); } +void CartesianCommunicator::AllToAllV(void *in ,const std::vector &sendcounts,const std::vector &senddispls, + void *out,const std::vector &recvcounts,const std::vector &recvdispls, + uint64_t bytes) +{ + // Single rank: the exchange degenerates to a copy of our own segment + GRID_ASSERT(sendcounts.size()==1); + GRID_ASSERT(recvcounts.size()==1); + GRID_ASSERT(sendcounts[0]==recvcounts[0]); + bcopy((char *)in +(uint64_t)senddispls[0]*bytes, + (char *)out+(uint64_t)recvdispls[0]*bytes,bytes*(uint64_t)sendcounts[0]); +} int CartesianCommunicator::RankWorld(void){return 0;} void CartesianCommunicator::Barrier(void){} diff --git a/Grid/parallelIO/BinaryIO.cc b/Grid/parallelIO/BinaryIO.cc index ef1b6683a..c2fda3a29 100644 --- a/Grid/parallelIO/BinaryIO.cc +++ b/Grid/parallelIO/BinaryIO.cc @@ -2,3 +2,7 @@ int Grid::BinaryIO::latticeWriteMaxRetry = -1; Grid::BinaryIO::IoPerf Grid::BinaryIO::lastPerf; + +// Target size of a single contiguous file extent under BINARYIO_AGGREGATE. +// 4MB is around the knee for Lustre; exposed so it can be swept at runtime. +uint64_t Grid::BinaryIO::aggregateTargetBytes = 4*1024*1024; diff --git a/Grid/parallelIO/BinaryIO.h b/Grid/parallelIO/BinaryIO.h index 4df9fdf9e..af28fe3ac 100644 --- a/Grid/parallelIO/BinaryIO.h +++ b/Grid/parallelIO/BinaryIO.h @@ -39,6 +39,7 @@ #endif #include +#include #include NAMESPACE_BEGIN(Grid); @@ -87,6 +88,7 @@ class BinaryIO { static IoPerf lastPerf; static int latticeWriteMaxRetry; + static uint64_t aggregateTargetBytes; ///////////////////////////////////////////////////////////////////////////// // more byte manipulation helpers @@ -253,12 +255,392 @@ class BinaryIO { // Read or Write distributed lexico array of ANY object to a specific location in file ////////////////////////////////////////////////////////////////////////////////////// + static const int BINARYIO_AGGREGATE = 0x20; static const int BINARYIO_MASTER_APPEND = 0x10; static const int BINARYIO_UNORDERED = 0x08; static const int BINARYIO_LEXICOGRAPHIC = 0x04; static const int BINARYIO_READ = 0x02; static const int BINARYIO_WRITE = 0x01; +#ifdef USE_MPI_IO + ///////////////////////////////////////////////////////////////////////////// + // Aggregation: self controlled transposition onto an I/O friendly layout. + // + // Under BINARYIO_LEXICOGRAPHIC the subarray file view handed to MPI-IO has + // contiguous runs of only lLattice[0]*sizeof(fobj) bytes -- a few KB for + // typical local volumes. Rather than rely on collective buffering to repair + // that, redistribute the payload ourselves so every rank owns a contiguous + // range of the global lexicographic site ordering, then issue large plain + // contiguous writes. + // + // "Un-splitting" the nunsplit fastest dimensions means the row of ranks + // sharing the remaining process coordinates collectively owns whole global + // hyperplanes. All data movement is then confined to that row communicator. + // Every rank still owns exactly lSites() sites afterwards, so the exchange is + // a pure permutation and needs no divisibility condition on the process grid. + ///////////////////////////////////////////////////////////////////////////// + struct AggregationPlan { + int nunsplit{0}; // number of fastest dimensions un-split + int rowsize{0}; // ranks in the aggregation (row) communicator + int rowrank{0}; // our logical (lexicographic) index within the row + uint64_t lsites{0}; // sites per rank -- invariant under the permutation + uint64_t chunk{0}; // sites in one globally contiguous run owned by the row + std::unique_ptr rowcomm; + // counts and displacements are indexed by rank within rowcomm + std::vector sendcounts, senddispls, recvcounts, recvdispls; + std::vector scatter; // recv slot -> slot in the aggregated buffer + std::vector extentGsite; // global lex site index of extent start + std::vector extentLocal; // offset of extent within aggregated buffer + std::vector extentSites; // sites in this extent + }; + + static inline void BuildAggregationPlan(GridBase *grid,uint64_t fobjSize,AggregationPlan &p) + { + int ndim = grid->Dimensions(); + Coordinate psizes = grid->ProcessorGrid(); + Coordinate pcoor = grid->ThisProcessorCoor(); + Coordinate gLattice= grid->GlobalDimensions(); + Coordinate lLattice= grid->LocalDimensions(); + Coordinate lstart = grid->LocalStarts(); + + uint64_t lsites = grid->lSites(); + p.lsites = lsites; + + ////////////////////////////////////////////////////////////////////////// + // Un-splitting dims 0..k-1 gives the row a contiguous run of + // chunk(k) = prod_{d= aggregateTargetBytes ) { k = trial; break; } + } + p.nunsplit = k; + + ////////////////////////////////////////////////////////////////////////// + // The box the row collectively owns, expressed in global coordinates. + // Restricting the global lexicographic order to this box preserves the + // ordering, so the row index below is monotone in the global index. + ////////////////////////////////////////////////////////////////////////// + Coordinate B(ndim), S(ndim); + for(int d=0; dProcessorCount(); + + ////////////////////////////////////////////////////////////////////////// + // Our logical index in the row is the forward lexicographic index of the + // un-split process coordinates, so that increasing logical index means + // increasing global lexicographic position in the file. The communicator + // numbers its own ranks by the reversed (MPI) convention, so build the map + // between the two rather than assuming either. + ////////////////////////////////////////////////////////////////////////// + int64_t logical=0, lstride=1; + for(int d=0; d commOf(p.rowsize,0); + commOf[p.rowrank] = (uint64_t)p.rowcomm->ThisRank(); + p.rowcomm->GlobalSumVector(&commOf[0],p.rowsize); + + uint64_t mystart = (uint64_t)p.rowrank * lsites; + uint64_t myend = mystart + lsites; + + Coordinate lcoor(ndim), bcoor(ndim), gcoor(ndim); + + ////////////////////////////////////////////////////////////////////////// + // Send side. Walking our local sites in local lexicographic order walks + // the row index monotonically, so the send buffer is iodata untouched and + // we need only the per destination counts. + ////////////////////////////////////////////////////////////////////////// + std::vector sendLogical(p.rowsize,0); + for(uint64_t L=0; L recvLogical(p.rowsize,0), recvDisplLogical(p.rowsize,0); + std::vector source(lsites); + for(uint64_t pos=0; pos fill(p.rowsize,0); + for(uint64_t pos=0; pos sendc(p.rowsize),recvc(p.rowsize); + for(int c=0;cAllToAll(&sendc[0],&recvc[0],1,sizeof(uint64_t)); + for(int c=0;cIsBoss() ) return; + std::cout << GridLogMessage << "IOobject: aggregate " << what + << " un-splitting " << p.nunsplit << " fastest dimensions, row of " + << p.rowsize << " ranks" << std::endl; + std::cout << GridLogMessage << "IOobject: aggregate " << p.extentSites.size() + << " extent(s)/rank, first " << p.extentSites[0]*fobjSize/1024./1024. << " MB" + << " (target " << aggregateTargetBytes/1024./1024. << " MB)" << std::endl; + std::cout << GridLogMessage << "IOobject: aggregate buffer overhead " + << p.lsites*fobjSize/1024./1024. << " MB/rank" << std::endl; + } + + //////////////////////////////////////////////////////////////////////////// + // Stage timings. The interesting quantity is the slowest rank, since every + // stage is followed sooner or later by a synchronisation, so reduce with + // GlobalMax rather than reporting whatever the boss happened to see. + //////////////////////////////////////////////////////////////////////////// + static inline void ReportStages(GridBase *grid,const char *what, + const std::vector &names, + std::vector &useconds) + { + GRID_ASSERT(names.size()==useconds.size()); + for(uint64_t i=0;iGlobalMax(useconds[i]); + if ( grid->IsBoss() ) { + std::cout << GridLogMessage << "IOobject: aggregate " << what << " stages (max over ranks, s):"; + for(uint64_t i=0;i + static inline void AggregateExchange(GridBase *grid,AggregationPlan &p,std::vector &iodata, + std::vector &aggregated,int forward) + { + uint64_t lsites = p.lsites; + GridStopWatch talloc,tperm,tcomm; + + talloc.Start(); + std::vector tmp(lsites); + talloc.Stop(); + + if ( forward ) { // iodata (local order) -> aggregated (lexicographic order) + tcomm.Start(); + p.rowcomm->AllToAllV(&iodata[0],p.sendcounts,p.senddispls, + &tmp[0], p.recvcounts,p.recvdispls,sizeof(fobj)); + tcomm.Stop(); + tperm.Start(); + thread_for(s,lsites,{ aggregated[p.scatter[s]] = tmp[s]; }); + tperm.Stop(); + } else { // aggregated -> iodata, the exact mirror + tperm.Start(); + thread_for(s,lsites,{ tmp[s] = aggregated[p.scatter[s]]; }); + tperm.Stop(); + tcomm.Start(); + p.rowcomm->AllToAllV(&tmp[0], p.recvcounts,p.recvdispls, + &iodata[0],p.sendcounts,p.senddispls,sizeof(fobj)); + tcomm.Stop(); + } + + std::vector us = { (RealD)talloc.useconds(), (RealD)tperm.useconds(), (RealD)tcomm.useconds() }; + ReportStages(grid,forward?"exchange (write)":"exchange (read)", + {"alloc","permute","alltoallv"},us); + } + + template + static inline void AggregateWrite(GridBase *grid,AggregationPlan &p,std::vector &aggregated, + std::string file,uint64_t offset) + { + ////////////////////////////////////////////////////////////////////////// + // All ranks write concurrently into a shared file, so the file must exist + // before any of them open it for update, but it does NOT have to be the + // right length first: the extents tile the record exactly, so writing them + // extends a short file to precisely offset+payload. + // + // Records are created in sequence, so this payload ends the file: the + // length must end up precisely offset+payload. Anything beyond is left + // over from whatever the file previously held and must not survive -- a + // shorter new record written over a longer old one would otherwise leave + // a trailing fragment of the previous contents masquerading as data. + // That is the only case needing a truncate, so stat first and truncate + // afterwards only when the size actually came out wrong. Measured on + // Frontier, an unconditional truncate up front cost 0.22 to 5.4 s per + // record -- 15 to 25% of a 19 GB write and 100% of a small one -- while + // create, open and close together cost a few milliseconds. It is per + // record, so multi record files do not amortise it away. + // + // ::truncate is used because the C++ standard library cannot express this. + // std::filebuf has no length operation at all; ios::trunc only truncates to + // zero; seeking past the end and writing a byte can grow a file but never + // shrink one; and there is no portable way to recover a descriptor from a + // stream in order to call ftruncate. C++17 does finally offer + // std::filesystem::resize_file, but that would be Grid's first + // dependency and needs -lstdc++fs on the older toolchains still in use. + ////////////////////////////////////////////////////////////////////////// + GridStopWatch tcreate,ttrunc,tbar,topen,twrite,tclose,tskew; + uint64_t need = offset + (uint64_t)grid->_gsites*sizeof(fobj); + + tcreate.Start(); + if ( grid->IsBoss() ) { + // opening for update needs the file to exist; create one only if not + std::fstream probe(file,std::ios::binary|std::ios::out|std::ios::in); + if ( !probe.is_open() ) { + std::ofstream create(file,std::ios::binary|std::ios::out); + create.close(); + } + } + tcreate.Stop(); + + tbar.Start(); + grid->Barrier(); + tbar.Stop(); + + std::ofstream fout; + fout.exceptions( std::fstream::failbit | std::fstream::badbit ); + try { + topen.Start(); + fout.open(file,std::ios::binary|std::ios::out|std::ios::in); + topen.Stop(); + twrite.Start(); + for(uint64_t e=0;eBarrier(); // every extent must be on its way first + tskew.Stop(); + + ttrunc.Start(); + if ( grid->IsBoss() ) { + struct stat sb; + int ierr = ::stat(file.c_str(),&sb); + GRID_ASSERT(ierr==0); + if ( (uint64_t)sb.st_size != need ) { // only when a longer record preceded us + ierr = ::truncate(file.c_str(),(off_t)need); + GRID_ASSERT(ierr==0); + } + } + grid->Barrier(); + ttrunc.Stop(); + + std::vector us = { (RealD)tcreate.useconds(), (RealD)tbar.useconds(), + (RealD)topen.useconds(), (RealD)twrite.useconds(), + (RealD)tclose.useconds(), (RealD)tskew.useconds(), + (RealD)ttrunc.useconds() }; + ReportStages(grid,"write",{"create","barrier","open","seek+write","close","skew","stat+truncate"},us); + } + + template + static inline void AggregateRead(GridBase *grid,AggregationPlan &p,std::vector &aggregated, + std::string file,uint64_t offset) + { + GridStopWatch topen,tread,tclose; + std::ifstream fin; + topen.Start(); + fin.open(file,std::ios::binary|std::ios::in); + topen.Stop(); + tread.Start(); + for(uint64_t e=0;e us = { (RealD)topen.useconds(), (RealD)tread.useconds(), (RealD)tclose.useconds() }; + ReportStages(grid,"read",{"open","seek+read","close"},us); + } +#endif + template static inline void IOobject(word w, GridBase *grid, @@ -302,6 +684,18 @@ class BinaryIO { lStart[d] = 0; } + ////////////////////////////////////////////////////////////////////////////// + // Aggregate the lexicographic layout onto contiguous per rank extents + // ourselves rather than leaving it to MPI-IO collective buffering + ////////////////////////////////////////////////////////////////////////////// + int aggregate = (control & BINARYIO_AGGREGATE) + && (control & BINARYIO_LEXICOGRAPHIC) + && !(control & BINARYIO_MASTER_APPEND) + && (nrank > 1); +#ifndef USE_MPI_IO + GRID_ASSERT(aggregate==0); // BINARYIO_AGGREGATE requires MPI +#endif + #ifdef USE_MPI_IO std::vector distribs(ndim,MPI_DISTRIBUTE_BLOCK); std::vector dargs (ndim,MPI_DISTRIBUTE_DFLT_DARG); @@ -329,6 +723,8 @@ class BinaryIO { ierr = MPI_Type_contiguous(numword,mpiword,&mpiObject); GRID_ASSERT(ierr==0); ierr = MPI_Type_commit(&mpiObject); + // The subarray view is what aggregation exists to avoid; do not build it + if ( !aggregate ) { ////////////////////////////////////////////////////////////////////////////// // File global array data type ////////////////////////////////////////////////////////////////////////////// @@ -340,6 +736,7 @@ class BinaryIO { ////////////////////////////////////////////////////////////////////////////// ierr=MPI_Type_create_subarray(ndim,&lLattice[0],&lLattice[0],&lStart[0],MPI_ORDER_FORTRAN, mpiObject,&localArray); GRID_ASSERT(ierr==0); ierr=MPI_Type_commit(&localArray); GRID_ASSERT(ierr==0); + } #endif ////////////////////////////////////////////////////////////////////////////// @@ -358,7 +755,19 @@ class BinaryIO { timer.Start(); - if ( (control & BINARYIO_LEXICOGRAPHIC) && (nrank > 1) ) { + if ( aggregate ) { +#ifdef USE_MPI_IO + std::cout<< GridLogMessage<<"IOobject: aggregate read I/O "<< file<< std::endl; + AggregationPlan plan; + BuildAggregationPlan(grid,sizeof(fobj),plan); + ReportAggregationPlan(grid,plan,sizeof(fobj),"read"); + std::vector aggregated(lsites); + AggregateRead(grid,plan,aggregated,file,offset); + AggregateExchange(grid,plan,iodata,aggregated,0); +#else + GRID_ASSERT(0); +#endif + } else if ( (control & BINARYIO_LEXICOGRAPHIC) && (nrank > 1) ) { #ifdef USE_MPI_IO std::cout<< GridLogMessage<<"IOobject: MPI read I/O "<< file<< std::endl; ierr=MPI_File_open(grid->communicator,(char *) file.c_str(), MPI_MODE_RDONLY, MPI_INFO_NULL, &fh); GRID_ASSERT(ierr==0); @@ -416,7 +825,25 @@ class BinaryIO { grid->Barrier(); timer.Start(); - if ( (control & BINARYIO_LEXICOGRAPHIC) && (nrank > 1) ) { + if ( aggregate ) { +#ifdef USE_MPI_IO + std::cout << GridLogMessage <<"IOobject: aggregate write I/O " << file << std::endl; + AggregationPlan plan; + BuildAggregationPlan(grid,sizeof(fobj),plan); + ReportAggregationPlan(grid,plan,sizeof(fobj),"write"); + std::vector aggregated(lsites); + AggregateExchange(grid,plan,iodata,aggregated,1); + AggregateWrite(grid,plan,aggregated,file,offset); + //////////////////////////////////////////////////////////////////////// + // Not every rank ends at the end of the payload, so the position can + // not be recovered from a file handle. Callers (Lime record chaining) + // rely on this being the first byte past the record. + //////////////////////////////////////////////////////////////////////// + offset = offset + (uint64_t)grid->_gsites*sizeof(fobj); +#else + GRID_ASSERT(0); +#endif + } else if ( (control & BINARYIO_LEXICOGRAPHIC) && (nrank > 1) ) { #ifdef USE_MPI_IO std::cout << GridLogMessage <<"IOobject: MPI write I/O " << file << std::endl; ierr = MPI_File_open(grid->communicator, (char *)file.c_str(), MPI_MODE_RDWR | MPI_MODE_CREATE, MPI_INFO_NULL, &fh); @@ -461,12 +888,26 @@ class BinaryIO { std::ofstream fout; fout.exceptions ( std::fstream::failbit | std::fstream::badbit ); + + //////////////////////////////////////////////////////////////////// + // Grid's model is that the boss rank performs the metadata + // operations and every other rank only seeks and writes into a file + // that already exists. Opening with ios::out on all ranks broke that: + // it is O_TRUNC, so a rank opening late truncated the file back to + // zero after an earlier rank had written its segment, leaving a hole + // in its place. The barriers around this block are outside it and do + // not order the opens against the writes. Let the boss create and + // empty the file, then everyone opens for update only. Same resulting + // length, one metadata operation instead of one per rank, no race. + //////////////////////////////////////////////////////////////////// + if ( !offset && grid->IsBoss() ) { // offset zero: this record starts the file + std::ofstream create(file,std::ios::binary|std::ios::out); + create.close(); + } + grid->Barrier(); + try { - if (offset) { // Must already exist and contain data - fout.open(file,std::ios::binary|std::ios::out|std::ios::in); - } else { // Allow create - fout.open(file,std::ios::binary|std::ios::out); - } + fout.open(file,std::ios::binary|std::ios::out|std::ios::in); } catch (const std::fstream::failure& exc) { std::cout << GridLogError << "Error in opening the file " << file << " for output" < + + 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 */ + +// Correctness and performance test for BINARYIO_AGGREGATE. +// +// Correctness, per aggregateTargetBytes: +// 1. write via the MPI-IO lexicographic path -> ref.bin +// 2. write via the aggregate path -> agg.bin +// 3. the two files must be byte identical <- proves the layout matches +// 4. read agg.bin back through the aggregate path <- proves the mirror inverts +// 5. write and read the non-lexicographic path -> raw.bin. Its layout is +// different by construction (each rank owns one contiguous segment in +// rank order) so it cannot be compared byte for byte, but the NERSC and +// SciDAC checksums are computed from the global site index and are +// therefore layout independent: they must match the other two paths. +// 6. a record written over a longer pre-existing file must leave the file at +// exactly offset+payload, with no trailing fragment of the old contents +// +// Performance: three paths, both directions, timed with the client page cache +// dropped before every read so that a read back reports filesystem bandwidth +// rather than memory bandwidth. The non-lexicographic path is the zero +// overhead reference: no transposition, no layout independence, one disjoint +// contiguous segment per rank, which is the arrangement that reaches full +// filesystem bandwidth on a leadership machine. It is the upper bound the +// other two are trying to approach. +// +// Options: +// --aggregate-target sweep this one target only (default: sweep +// 1, 1024, 64K, 4M) +// --io-reps repetitions in the performance section +// (default 3; 0 disables it) +// --io-no-correctness skip the correctness section, which reads the +// whole file on one rank and is not affordable +// at very large volume +// --io-read-only time reads only, of files left in place by an +// earlier job. Reading back what this job just +// wrote measures the client page cache; a fresh +// allocation pointed at the same directory is +// the only way to get a cold read without root. +// +// The exchange is only meaningfully exercised when the fast dimensions are +// split across ranks; --mpi 1.1.X.Y leaves the rows of size one and the test +// then passes vacuously. Non-uniform AllToAllV counts additionally need an +// odd process factor in a fast dimension and a small local volume. + +#include +#include + +using namespace Grid; + +///////////////////////////////////////////////////////////////////////////// +// Compare in chunks. Slurping both files into memory is fine for a few MB +// and fatal for the multi-GB records this test is meant to reach. +///////////////////////////////////////////////////////////////////////////// +static bool FilesIdentical(std::string a,std::string b) +{ + std::ifstream fa(a,std::ios::binary), fb(b,std::ios::binary); + if ( !fa.good() || !fb.good() ) { + std::cout< va(chunk), vb(chunk); + uint64_t done=0; + while ( done < sa ) { + uint64_t n = std::min(chunk,sa-done); + fa.read(&va[0],n); + fb.read(&vb[0],n); + for(uint64_t i=0;iBarrier(); + int fd = ::open(file.c_str(),O_RDONLY); + if ( fd >= 0 ) { +#ifdef POSIX_FADV_DONTNEED + ::posix_fadvise(fd,0,0,POSIX_FADV_DONTNEED); +#endif + ::close(fd); + } + grid->Barrier(); +} + +static uint64_t OptionU64(int argc,char **argv,const char *opt,uint64_t def) +{ + if ( GridCmdOptionExists(argv,argv+argc,opt) ) { + std::string arg = GridCmdOptionPayload(argv,argv+argc,opt); + return (uint64_t)std::stoull(arg); + } + return def; +} + +int main(int argc,char **argv) +{ + Grid_init(&argc,&argv); + + Coordinate latt = GridDefaultLatt(); + Coordinate simd = GridDefaultSimd(Nd,vComplexD::Nsimd()); + Coordinate mpi = GridDefaultMpi(); + GridCartesian grid(latt,simd,mpi); + + typedef vLorentzColourMatrixD vobj; + typedef LorentzColourMatrixD sobj; + + GridParallelRNG pRNG(&grid); + pRNG.SeedFixedIntegers(std::vector({1,2,3,4})); + LatticeGaugeFieldD Umu(&grid); + random(pRNG,Umu); + + BinarySimpleMunger munge; + const std::string format("IEEE64BIG"); + + const int lex = BinaryIO::BINARYIO_LEXICOGRAPHIC; + const int agg = BinaryIO::BINARYIO_LEXICOGRAPHIC|BinaryIO::BINARYIO_AGGREGATE; + const int raw = 0; // no BINARYIO_LEXICOGRAPHIC: contiguous segment per rank + + uint64_t payload = (uint64_t)grid._gsites*sizeof(sobj); + + std::vector targets = {1, 1024, 64*1024, 4*1024*1024}; + if ( GridCmdOptionExists(argv,argv+argc,"--aggregate-target") ) { + targets.clear(); + targets.push_back(OptionU64(argc,argv,"--aggregate-target",4*1024*1024)); + } + uint64_t reps = OptionU64(argc,argv,"--io-reps",3); + bool correctness = !GridCmdOptionExists(argv,argv+argc,"--io-no-correctness"); + // Read only: time reads of files left by an earlier job. The only way to + // get a cold client cache without root is to read on an allocation that did + // not write the data, so run one job to write and a second, pointed at the + // same directory, with this flag. + bool readonly = GridCmdOptionExists(argv,argv+argc,"--io-read-only"); + if ( readonly ) correctness = false; + + std::cout<(Umu,"ref.bin",munge,off,format,n1,a1,b1,lex); + + BinaryIO::aggregateTargetBytes = target; + off=0; + BinaryIO::writeLatticeObject(Umu,"agg.bin",munge,off,format,n2,a2,b2,agg); + + grid.Barrier(); + + if ( grid.IsBoss() ) { + if ( !FilesIdentical("ref.bin","agg.bin") ) { + std::cout<(Uchk,"agg.bin",munge,off,format,n3,a3,b3,agg); + if ( (n3!=n1)||(a3!=a1)||(b3!=b1) ) { + std::cout<(Umu,"raw.bin",munge,off,format,n4,a4,b4,raw); + grid.Barrier(); + if ( (n4!=n1)||(a4!=a1)||(b4!=b1) ) { + std::cout<(Uchk,"raw.bin",munge,off,format,n5,a5,b5,raw); + if ( (n5!=n1)||(a5!=a1)||(b5!=b1) ) { + std::cout<(Umu,"trunc.bin",munge,off,format,n,a,b,agg); + grid.Barrier(); + + if ( grid.IsBoss() ) { + std::ifstream f("trunc.bin",std::ios::binary|std::ios::ate); + uint64_t sz = (uint64_t)f.tellg(); + f.close(); + if ( sz != expect ) { + std::cout< wref,wagg,wraw,rref,ragg,rraw; + + for(uint64_t rep=0;rep(Umu,"ref.bin",munge,off,format,n,a,b,lex); + wref.push_back(BinaryIO::lastPerf.mbytesPerSecond); + cn[0]=n; ca[0]=a; cb[0]=b; + + off=0; + BinaryIO::writeLatticeObject(Umu,"agg.bin",munge,off,format,n,a,b,agg); + wagg.push_back(BinaryIO::lastPerf.mbytesPerSecond); + cn[1]=n; ca[1]=a; cb[1]=b; + + off=0; + BinaryIO::writeLatticeObject(Umu,"raw.bin",munge,off,format,n,a,b,raw); + wraw.push_back(BinaryIO::lastPerf.mbytesPerSecond); + cn[2]=n; ca[2]=a; cb[2]=b; + if ( !agreeing("write",0,2) ) failures++; + + } // !readonly + + DropCache(&grid,"ref.bin"); + off=0; + BinaryIO::readLatticeObject(Uio,"ref.bin",munge,off,format,n,a,b,lex); + rref.push_back(BinaryIO::lastPerf.mbytesPerSecond); + cn[3]=n; ca[3]=a; cb[3]=b; + + DropCache(&grid,"agg.bin"); + off=0; + BinaryIO::readLatticeObject(Uio,"agg.bin",munge,off,format,n,a,b,agg); + ragg.push_back(BinaryIO::lastPerf.mbytesPerSecond); + cn[4]=n; ca[4]=a; cb[4]=b; + + DropCache(&grid,"raw.bin"); + off=0; + BinaryIO::readLatticeObject(Uio,"raw.bin",munge,off,format,n,a,b,raw); + rraw.push_back(BinaryIO::lastPerf.mbytesPerSecond); + cn[5]=n; ca[5]=a; cb[5]=b; + if ( !agreeing("read back",readonly?3:0,5) ) failures++; + } + + if ( grid.IsBoss() ) { + auto report = [&](const char *name,std::vector &v){ + if ( v.empty() ) return; + RealD best=0, sum=0; + for(auto x : v){ if(x>best) best=x; sum+=x; } + // First sample includes file creation, later ones do not; quote both + // rather than a mean that mixes the two. + std::cout< &v){ RealD m=0; for(auto x:v) if(x>m) m=x; return m; }; + if ( !wraw.empty() && best(wraw) > 0 ) { + std::cout< 0 ) { + std::cout<