diff --git a/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h b/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h index 14344d708..24451cc6e 100644 --- a/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h +++ b/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h @@ -88,16 +88,14 @@ public: double telLeafMaxInv; uint64_t nLeaf; // BIG LEAVES. Below span s blocks a sub-block lives on <= s of the Pr - // process rows / s of the Pc columns; when s is small RELATIVE TO THE - // GRID the SUMMA rings run on a few ranks while the rest block in their - // next SendToRecvFrom (histogram 2026-08-27: 93% of ring time in the - // 3.7 MB single-block panels of exactly these levels). Instead: gather - // the (s*nb)^2 sub-block to one rank, invert locally (inverseLU), scatter - // back. Fires only while span < min(Pr,Pc) -- when the whole grid - // participates the rings are not degenerate and gathering would only - // concentrate memory (and at the top of the tree, gather the whole - // matrix). Default 9: banked on Frontier 2026-08-27 (invert 132 -> 27.6 s - // at N=138240 on a 16x18 grid; s=18 gave 30.6, s=4 37.0). + // process rows / s of the Pc columns; when s is small relative to the grid + // the SUMMA rings run on a few ranks while the rest block in their next + // SendToRecvFrom. Instead: gather the (s*nb)^2 sub-block to one rank, + // invert locally (inverseLU), scatter back. Fires only while + // span < min(Pr,Pc) -- when the whole grid participates the rings are not + // degenerate and gathering would only concentrate memory (and at the top + // of the tree would gather the whole matrix). Larger s concentrates more + // memory on the root; smaller s loses ring parallelism. int leafSpan = 9; uint64_t nBigLeaf = 0; int64_t maxBigW = 0; double tBigGather = 0, tBigInv = 0, tBigScatter = 0; diff --git a/Grid/algorithms/multigrid/DenseCoarseMatrix.h b/Grid/algorithms/multigrid/DenseCoarseMatrix.h index 51bbbecc3..67cd47fc2 100644 --- a/Grid/algorithms/multigrid/DenseCoarseMatrix.h +++ b/Grid/algorithms/multigrid/DenseCoarseMatrix.h @@ -51,9 +51,7 @@ NAMESPACE_BEGIN(Grid); // complement (BlockCyclicSchurInverse): fp64 rank-major import -> // RowsToCyclic -> in-place recursion (pure point-to-point SUMMA rings and // local leaves; bitwise reproducible) -> CyclicToRows -> ONE terminal -// rounding into the fp32 apply slab. Distributed at every N and P; banked -// 3.87x against its retired 1D predecessor and >10x against SLATE at -// N=138240 on 288 GCDs. +// rounding into the fp32 apply slab. Distributed at every N and P. // // - Split-K apply through GridBLAS.gemmBatched with EXPLICIT leading dimensions // (arXiv:2409.03904 fig 11): the tiny-output/huge-K GEMM Y = slab^T X becomes @@ -213,8 +211,7 @@ public: GRID_ASSERT( l2r[myLex] == grid->ThisRank() ); // x and the slab columns are in GLOBAL-SITE order (hX[myGsite*nbasis+b]); // the gathered blocks are in RANK-MAJOR order (rank*nrows + ss*nbasis + b). - // On one rank the two coincide, which is how the laptop passed while - // Frontier VERIFYed 0.9965 (2026-08-27). Scatter through the inverse map. + // The two coincide only on one rank, so scatter through the inverse map. std::vector g2rm; BuildRankMajorMap(g2rm); std::vector rm2g(N); for(int64_t g=0; gThisRank()*nrows + (int64_t)ss*nbasis] == myGsite[ss]*nbasis ); @@ -320,11 +317,7 @@ public: // The operator contracts out(s,b) = sum_a A[p](s)(a,b) in(nbr,a) // (GeneralCoarsenedMatrix.h Mult kernel): the stored site matrix // acts TRANSPOSED, so element (a,b) lands at dense row (s,b), - // column (nbr,a). BUG LEDGER 2026-08-14: the original mapping - // wrote (s,a),(nbr,b) -- caught by the IMPORT CERTIFICATE on its - // FIRST fresh-import exercise (Test_schur_dense_coarse); every - // production slab predated this path (probe-import era), so no - // production output is suspect. + // column (nbr,a). for(int b=0; b fp64 rank-major import -> - // fp64 recursion -> ONE terminal rounding into the fp32 apply - // slab. Everything downstream (device residency, split-K apply, - // VERIFY) is untouched. + // 3c. The inverse: distributed recursive Schur, END-TO-END fp64. + // stencil (ComplexD) -> fp64 rank-major import -> fp64 recursion -> + // ONE terminal rounding into the fp32 apply slab. Everything + // downstream (device residency, split-K apply, VERIFY) is fp32. //////////////////////////////////////////////////////////////////// template void InvertDense(CoarseOp &Op) diff --git a/Grid/algorithms/multigrid/PVdagMMultiGridParams.h b/Grid/algorithms/multigrid/PVdagMMultiGridParams.h index aa4bbfaa7..adfabd925 100644 --- a/Grid/algorithms/multigrid/PVdagMMultiGridParams.h +++ b/Grid/algorithms/multigrid/PVdagMMultiGridParams.h @@ -29,17 +29,13 @@ NAMESPACE_BEGIN(Grid); // sub-structs for the Hermitian chain. Read from XML (or JSON -- the // serialisation macros give both) via ReadPVdagMMultiGridParams below and // printed at startup by the macro's operator<<, so every log identifies its -// own run. The constructor defaults ARE the banked optimum; an -// unconfigured run reproduces the best recorded point. Update them when a -// better point is banked, and date the change. +// own run. The constructor defaults ARE a tuned operating point, not +// arbitrary: an unconfigured run reproduces it, so changing them changes what +// an unconfigured run does. Smoother mmax == nstep (full GCR history). // -// Current optimum: 2026-08-24, slurm-5335492 F4, 48^3x96 Ls=24 on 288 GCDs, -// 17.2 s/RHS at Nrhs=4, 32.2 s at Nrhs=1 (exact-halo FINAL ~1e-8). -// Smoother mmax == nstep (full GCR history). -// -// There are NO environment-variable controls anywhere in this subsystem -// (ruling 2026-09-05): parameters live here, library-internal constants are -// hard defaults in their classes, verbosity is the --log channels. +// There are NO environment-variable controls anywhere in this subsystem: +// parameters live here, library-internal constants are hard defaults in their +// classes, verbosity is the --log channels. // Research instruments (GCR coefficient recording/replay, Chebyshev and // stationary smoother variants, power iteration) are deliberately NOT part // of this interface: they remain programmatic, for algorithmic studies, @@ -77,7 +73,7 @@ struct MGOuterParams : Serializable { struct MGDenseParams : Serializable { GRID_SERIALIZABLE_CLASS_MEMBERS(MGDenseParams, - int, LeafSpan); // big-leaf span in blocks (banked 9); see BlockCyclicSchurInverse + int, LeafSpan); // big-leaf span in blocks; see BlockCyclicSchurInverse MGDenseParams() : LeafSpan(9) {}; }; diff --git a/Grid/algorithms/multigrid/PVdagMOperators.h b/Grid/algorithms/multigrid/PVdagMOperators.h index 6e2bef4ed..dcdde3d0b 100644 --- a/Grid/algorithms/multigrid/PVdagMOperators.h +++ b/Grid/algorithms/multigrid/PVdagMOperators.h @@ -63,11 +63,6 @@ public: void OpDir (const Field &in, Field &out,int dir,int disp) { GRID_ASSERT(0); } void OpDirAll (const Field &in, std::vector &out){ GRID_ASSERT(0); }; void Op (const Field &in, Field &out){ Field tmp(in.Grid()); _Mat.M(in,tmp); _PV.Mdag(tmp,out); out = out + shift*in; } - // BUG LEDGER 2026-09-06: the example-local original read - // _PV.M(tmp,out); _Mat.Mdag(in,tmp); - // consuming tmp before writing it -- latent, never exercised (nothing in - // the chain calls AdjOp on the shifted fine operator). Corrected here to - // the adjoint of Op, matching PVdagMLinearOperator::AdjOp plus the shift. void AdjOp (const Field &in, Field &out){ Field tmp(in.Grid()); _PV.M(in,tmp); _Mat.Mdag(tmp,out); out = out + shift*in; } void HermOpAndNorm(const Field &in, Field &out,RealD &n1,RealD &n2){ GRID_ASSERT(0); } void HermOp(const Field &in, Field &out){ Field tmp(in.Grid()); Op(in,tmp); AdjOp(tmp,out); } @@ -90,8 +85,4 @@ public: void HermOp (const Field &in, Field &out) { Field tmp(in.Grid()); Op(in,tmp); AdjOp(tmp,out); } }; -// The non-Hermitian spectral-edge diagnostic that used to live here as -// PowerIteration() now lives beside its Hermitian sibling as -// Grid::NonHermitianPowerMethod (Grid/algorithms/iterative/PowerMethod.h). - NAMESPACE_END(Grid); diff --git a/Grid/stencil/Stencil.h b/Grid/stencil/Stencil.h index 97c7b04c8..c1243657e 100644 --- a/Grid/stencil/Stencil.h +++ b/Grid/stencil/Stencil.h @@ -453,15 +453,11 @@ public: double *dbuf =(double *) packet.recv_buf; float *fbuf =(float *) packet.compressed_recv_buf; - // BUG FIX 2026-09-07: the lane structure is a carefully designed - // GPU optimisation -- lane = threadIdx.y puts adjacent threads on - // adjacent words, fully coalesced. On CPU builds there are no - // SIMT lanes and acceleratorSIMTlane() is identically zero, so - // this loop converted only 1/nsimd of the words and the rest of - // the halo silently kept STALE buffer content (deterministic - // wrong answers whenever the previous occupant differed; caught - // on the 2-rank laptop build via Benchmark_dwf's sloppy Cshift - // check). GRID_SIMT keeps the optimisation; CPU loops the lanes. + // The lane index is a GPU coalescing optimisation: lane = threadIdx.y + // puts adjacent threads on adjacent words. On a CPU build there are no + // SIMT lanes and acceleratorSIMTlane() is identically zero, so the plain + // form would convert only 1 word of every nsimd; the #else covers all + // lanes explicitly. accelerator_forNB(ss,outer,nsimd,{ #ifdef GRID_SIMT int lane = acceleratorSIMTlane(nsimd); @@ -534,7 +530,7 @@ public: double *dbuf =(double *) packet.send_buf; float *fbuf =(float *) packet.compressed_send_buf; - // BUG FIX 2026-09-07: CPU lane coverage -- see DecompressPacket. + // CPU lane coverage as in DecompressPacket. accelerator_forNB(ss,outer,nsimd,{ #ifdef GRID_SIMT int lane = acceleratorSIMTlane(nsimd); diff --git a/examples/Example_pvdagm_multigrid.cc b/examples/Example_pvdagm_multigrid.cc index 8c651fc7d..06bf00cb8 100644 --- a/examples/Example_pvdagm_multigrid.cc +++ b/examples/Example_pvdagm_multigrid.cc @@ -27,7 +27,7 @@ Author: Peter Boyle // ./Example_pvdagm_multigrid --grid 48.48.48.96 --mpi ... \ // --pvdagm-params params.xml // -// With no --pvdagm-params the banked defaults run (hot-start gauge field +// With no --pvdagm-params the built-in defaults run (hot-start gauge field // unless Config is set in the file). A missing file gets a template // written next to it and the program exits: the template documents every // parameter. The effective parameters are always printed, so the log