Better commenting

This commit is contained in:
Peter Boyle committed 2026-09-09 14:39:56 -04:00
1 parent ec682693c1
commit 7f6b0f409c
6 files changed
+29 -56

No files matched your search

@@ -88,16 +88,14 @@ public:
double telLeafMaxInv; double telLeafMaxInv;
uint64_t nLeaf; uint64_t nLeaf;
// BIG LEAVES. Below span s blocks a sub-block lives on <= s of the Pr // 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 // process rows / s of the Pc columns; when s is small relative to the grid
// GRID the SUMMA rings run on a few ranks while the rest block in their // the SUMMA rings run on a few ranks while the rest block in their next
// next SendToRecvFrom (histogram 2026-08-27: 93% of ring time in the // SendToRecvFrom. Instead: gather the (s*nb)^2 sub-block to one rank,
// 3.7 MB single-block panels of exactly these levels). Instead: gather // invert locally (inverseLU), scatter back. Fires only while
// the (s*nb)^2 sub-block to one rank, invert locally (inverseLU), scatter // span < min(Pr,Pc) -- when the whole grid participates the rings are not
// back. Fires only while span < min(Pr,Pc) -- when the whole grid // degenerate and gathering would only concentrate memory (and at the top
// participates the rings are not degenerate and gathering would only // of the tree would gather the whole matrix). Larger s concentrates more
// concentrate memory (and at the top of the tree, gather the whole // memory on the root; smaller s loses ring parallelism.
// 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).
int leafSpan = 9; int leafSpan = 9;
uint64_t nBigLeaf = 0; int64_t maxBigW = 0; uint64_t nBigLeaf = 0; int64_t maxBigW = 0;
double tBigGather = 0, tBigInv = 0, tBigScatter = 0; double tBigGather = 0, tBigInv = 0, tBigScatter = 0;
+7 -15
View File
@@ -51,9 +51,7 @@ NAMESPACE_BEGIN(Grid);
// complement (BlockCyclicSchurInverse): fp64 rank-major import -> // complement (BlockCyclicSchurInverse): fp64 rank-major import ->
// RowsToCyclic -> in-place recursion (pure point-to-point SUMMA rings and // RowsToCyclic -> in-place recursion (pure point-to-point SUMMA rings and
// local leaves; bitwise reproducible) -> CyclicToRows -> ONE terminal // local leaves; bitwise reproducible) -> CyclicToRows -> ONE terminal
// rounding into the fp32 apply slab. Distributed at every N and P; banked // rounding into the fp32 apply slab. Distributed at every N and P.
// 3.87x against its retired 1D predecessor and >10x against SLATE at
// N=138240 on 288 GCDs.
// //
// - Split-K apply through GridBLAS.gemmBatched with EXPLICIT leading dimensions // - 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 // (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() ); GRID_ASSERT( l2r[myLex] == grid->ThisRank() );
// x and the slab columns are in GLOBAL-SITE order (hX[myGsite*nbasis+b]); // 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). // 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 // The two coincide only on one rank, so scatter through the inverse map.
// Frontier VERIFYed 0.9965 (2026-08-27). Scatter through the inverse map.
std::vector<int64_t> g2rm; BuildRankMajorMap(g2rm); std::vector<int64_t> g2rm; BuildRankMajorMap(g2rm);
std::vector<int64_t> rm2g(N); for(int64_t g=0; g<N; g++) rm2g[g2rm[g]] = g; std::vector<int64_t> rm2g(N); for(int64_t g=0; g<N; g++) rm2g[g2rm[g]] = g;
for(int ss=0; ss<lsites; ss++) GRID_ASSERT( rm2g[(int64_t)grid->ThisRank()*nrows + (int64_t)ss*nbasis] == myGsite[ss]*nbasis ); for(int ss=0; ss<lsites; ss++) GRID_ASSERT( rm2g[(int64_t)grid->ThisRank()*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) // The operator contracts out(s,b) = sum_a A[p](s)(a,b) in(nbr,a)
// (GeneralCoarsenedMatrix.h Mult kernel): the stored site matrix // (GeneralCoarsenedMatrix.h Mult kernel): the stored site matrix
// acts TRANSPOSED, so element (a,b) lands at dense row (s,b), // acts TRANSPOSED, so element (a,b) lands at dense row (s,b),
// column (nbr,a). BUG LEDGER 2026-08-14: the original mapping // column (nbr,a).
// 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.
for(int b=0; b<nbasis; b++){ for(int b=0; b<nbasis; b++){
ComplexF *row = &slab[(uint64_t)(ss*nbasis+b)*N + nsite*nbasis]; ComplexF *row = &slab[(uint64_t)(ss*nbasis+b)*N + nsite*nbasis];
for(int a=0; a<nbasis; a++) for(int a=0; a<nbasis; a++)
@@ -552,11 +545,10 @@ public:
} }
//////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////
// 3c. The inverse: distributed recursive Schur, END-TO-END fp64 (decision // 3c. The inverse: distributed recursive Schur, END-TO-END fp64.
// 2026-08-14): stencil (ComplexD) -> fp64 rank-major import -> // stencil (ComplexD) -> fp64 rank-major import -> fp64 recursion ->
// fp64 recursion -> ONE terminal rounding into the fp32 apply // ONE terminal rounding into the fp32 apply slab. Everything
// slab. Everything downstream (device residency, split-K apply, // downstream (device residency, split-K apply, VERIFY) is fp32.
// VERIFY) is untouched.
//////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////
template<class CoarseOp> template<class CoarseOp>
void InvertDense(CoarseOp &Op) void InvertDense(CoarseOp &Op)
@@ -29,17 +29,13 @@ NAMESPACE_BEGIN(Grid);
// sub-structs for the Hermitian chain. Read from XML (or JSON -- the // sub-structs for the Hermitian chain. Read from XML (or JSON -- the
// serialisation macros give both) via ReadPVdagMMultiGridParams below and // serialisation macros give both) via ReadPVdagMMultiGridParams below and
// printed at startup by the macro's operator<<, so every log identifies its // printed at startup by the macro's operator<<, so every log identifies its
// own run. The constructor defaults ARE the banked optimum; an // own run. The constructor defaults ARE a tuned operating point, not
// unconfigured run reproduces the best recorded point. Update them when a // arbitrary: an unconfigured run reproduces it, so changing them changes what
// better point is banked, and date the change. // 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, // There are NO environment-variable controls anywhere in this subsystem:
// 17.2 s/RHS at Nrhs=4, 32.2 s at Nrhs=1 (exact-halo FINAL ~1e-8). // parameters live here, library-internal constants are hard defaults in their
// Smoother mmax == nstep (full GCR history). // classes, verbosity is the --log channels.
//
// 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.
// Research instruments (GCR coefficient recording/replay, Chebyshev and // Research instruments (GCR coefficient recording/replay, Chebyshev and
// stationary smoother variants, power iteration) are deliberately NOT part // stationary smoother variants, power iteration) are deliberately NOT part
// of this interface: they remain programmatic, for algorithmic studies, // of this interface: they remain programmatic, for algorithmic studies,
@@ -77,7 +73,7 @@ struct MGOuterParams : Serializable {
struct MGDenseParams : Serializable { struct MGDenseParams : Serializable {
GRID_SERIALIZABLE_CLASS_MEMBERS(MGDenseParams, 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) {}; MGDenseParams() : LeafSpan(9) {};
}; };
@@ -63,11 +63,6 @@ public:
void OpDir (const Field &in, Field &out,int dir,int disp) { GRID_ASSERT(0); } void OpDir (const Field &in, Field &out,int dir,int disp) { GRID_ASSERT(0); }
void OpDirAll (const Field &in, std::vector<Field> &out){ GRID_ASSERT(0); }; void OpDirAll (const Field &in, std::vector<Field> &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; } 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 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 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); } 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); } 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); NAMESPACE_END(Grid);
+6 -10
View File
@@ -453,15 +453,11 @@ public:
double *dbuf =(double *) packet.recv_buf; double *dbuf =(double *) packet.recv_buf;
float *fbuf =(float *) packet.compressed_recv_buf; float *fbuf =(float *) packet.compressed_recv_buf;
// BUG FIX 2026-09-07: the lane structure is a carefully designed // The lane index is a GPU coalescing optimisation: lane = threadIdx.y
// GPU optimisation -- lane = threadIdx.y puts adjacent threads on // puts adjacent threads on adjacent words. On a CPU build there are no
// adjacent words, fully coalesced. On CPU builds there are no // SIMT lanes and acceleratorSIMTlane() is identically zero, so the plain
// SIMT lanes and acceleratorSIMTlane() is identically zero, so // form would convert only 1 word of every nsimd; the #else covers all
// this loop converted only 1/nsimd of the words and the rest of // lanes explicitly.
// 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.
accelerator_forNB(ss,outer,nsimd,{ accelerator_forNB(ss,outer,nsimd,{
#ifdef GRID_SIMT #ifdef GRID_SIMT
int lane = acceleratorSIMTlane(nsimd); int lane = acceleratorSIMTlane(nsimd);
@@ -534,7 +530,7 @@ public:
double *dbuf =(double *) packet.send_buf; double *dbuf =(double *) packet.send_buf;
float *fbuf =(float *) packet.compressed_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,{ accelerator_forNB(ss,outer,nsimd,{
#ifdef GRID_SIMT #ifdef GRID_SIMT
int lane = acceleratorSIMTlane(nsimd); int lane = acceleratorSIMTlane(nsimd);
+1 -1
View File
@@ -27,7 +27,7 @@ Author: Peter Boyle <pboyle@bnl.gov>
// ./Example_pvdagm_multigrid --grid 48.48.48.96 --mpi ... \ // ./Example_pvdagm_multigrid --grid 48.48.48.96 --mpi ... \
// --pvdagm-params params.xml // --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 // 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 // written next to it and the program exits: the template documents every
// parameter. The effective parameters are always printed, so the log // parameter. The effective parameters are always printed, so the log