From 6cea81bc2c0073e3fa103f78565227c657b6fe92 Mon Sep 17 00:00:00 2001 From: Peter Boyle Date: Thu, 27 Aug 2026 13:42:18 -0400 Subject: [PATCH] More probing why test harness vs. slate was a little faster running Grid's dense inverse than normal environment --- Grid/algorithms/blas/BatchedInverse.h | 4 + .../multigrid/BlockCyclicSchurInverse.h | 44 +++++++++ systems/Frontier/schur2d_env.job | 91 +++++++++++++++++++ tests/debug/Test_schur2d_scale.cc | 32 +++++++ 4 files changed, 171 insertions(+) create mode 100644 systems/Frontier/schur2d_env.job diff --git a/Grid/algorithms/blas/BatchedInverse.h b/Grid/algorithms/blas/BatchedInverse.h index 40aaa107d..96b42a2ed 100644 --- a/Grid/algorithms/blas/BatchedInverse.h +++ b/Grid/algorithms/blas/BatchedInverse.h @@ -114,15 +114,18 @@ public: // NB: written to rocSOLVER's documented z*_64 signatures; not compiled on // HIP in the air-gapped loop -- verify on first hipcc build. /////////////////////////////////////////////////////////////////////////// + double lastGetrfUs = 0.0, lastGetrsUs = 0.0; // split timing of the last inverseLU (HIP path) void inverseLU(int64_t N, ComplexD *A) { #ifdef GRID_HIP rocblas_handle handle = Handle(); deviceVector ipiv((uint64_t)N); deviceVector info(1); + double t0 = usecond(); auto st1 = rocsolver_zgetrf_64(handle, N, N, (rocblas_double_complex *)A, N, &ipiv[0], &info[0]); GRID_ASSERT(st1 == rocblas_status_success); accelerator_barrier(); + lastGetrfUs = usecond()-t0; int64_t info_h = -1; acceleratorCopyFromDevice(&info[0], &info_h, sizeof(int64_t)); GRID_ASSERT(info_h == 0); deviceVector X((uint64_t)N*N); @@ -134,6 +137,7 @@ public: (rocblas_double_complex *)&X[0], N); GRID_ASSERT(st2 == rocblas_status_success); accelerator_barrier(); + lastGetrsUs = usecond()-t0-lastGetrfUs; acceleratorCopyDeviceToDevice((void *)&X[0], (void *)A, (uint64_t)N*N*sizeof(ComplexD)); #else deviceVector bp(1); std::vector ptr(1); ptr[0] = A; diff --git a/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h b/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h index fbe2a4fd0..5c1a96382 100644 --- a/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h +++ b/Grid/algorithms/multigrid/BlockCyclicSchurInverse.h @@ -439,6 +439,50 @@ public: } } /////////////////////////////////////////////////////////////////////// + // ONE-GCD LEAF MICROBENCHMARK (boss only). The big-leaf inverse at + // W=4320 measured 0.53 s with BOTH getri_batched and getrf_64+getrs_64 + // (2026-08-27) -- ~0.4 TF/s on a GCD that runs zgemm at ~15. Time the + // three primitives in isolation on a well-conditioned dense matrix so the + // leaf's rate can be compared with the machine's, and getrf split from + // getrs. Sizes: the W of span 4 / 9 / 18 leaves on accelerator builds; + // tiny on CPU builds (Eigen would take minutes at 4320). + /////////////////////////////////////////////////////////////////////// + if ( me == 0 ) { +#if defined(GRID_HIP) || defined(GRID_CUDA) || defined(GRID_SYCL) + std::vector Ws({1920, 4320, 8640}); +#else + std::vector Ws({240, 480}); +#endif + for(int64_t W : Ws){ + deviceVector M((uint64_t)W*W), C((uint64_t)W*W); + { ComplexD *m = &M[0]; const int64_t WW = W; // diagonally dominant: (i==j ? W : 0) + cos/sin noise + accelerator_for(idx,(uint64_t)W*W,1,{ int64_t j=idx/WW, i=idx-j*WW; double x=0.37*i+0.61*j; + m[idx] = ComplexD((i==j)?(double)WW:0.0,0.0) + ComplexD(std::cos(x),std::sin(1.3*x)); }); + accelerator_barrier(); } + double flopLU = 8.0/3.0*(double)W*W*W; // complex LU ~ (4 real flops per complex mult-add) * (2/3 n^3) + double flopGEMM = 8.0*(double)W*W*W; // complex GEMM + // 1. getri_batched (batch 1) + double tb; + { deviceVector bp(1); std::vector ptr(1); ptr[0]=&M[0]; + acceleratorCopyToDevice(&ptr[0],&bp[0],sizeof(ComplexD*)); + double t0=usecond(); INV.inverseBatched(W,bp); tb=usecond()-t0; } + // 2. inverseLU (getrf_64 + identity getrs_64), timed inside + double tl; { double t0=usecond(); INV.inverseLU(W,&M[0]); tl=usecond()-t0; } + // 3. one zgemm W x W x W for the machine rate + double tg; + { deviceVector ap(1),bp(1),cp(1); std::vector ptr(1); + ptr[0]=&M[0]; acceleratorCopyToDevice(&ptr[0],&ap[0],sizeof(ComplexD*)); acceleratorCopyToDevice(&ptr[0],&bp[0],sizeof(ComplexD*)); + ptr[0]=&C[0]; acceleratorCopyToDevice(&ptr[0],&cp[0],sizeof(ComplexD*)); + double t0=usecond(); + SUMMA.BLAS.gemmBatched(GridBLAS_OP_N,GridBLAS_OP_N,(int)W,(int)W,(int)W,ComplexD(1.0,0.0),ap,(int)W,bp,(int)W,ComplexD(0.0,0.0),cp,(int)W); + SUMMA.BLAS.synchronise(); tg=usecond()-t0; } + std::cout << GridLogMessage << "Schur2D PROBE leaf W=" << W + << ": getri_batched " << tb/1.0e6 << " s (" << flopLU/tb/1.0e6 << " TF/s)" + << " getrf_64+getrs_64 " << tl/1.0e6 << " s (getrf " << INV.lastGetrfUs/1.0e6 << " getrs " << INV.lastGetrsUs/1.0e6 << ")" + << " zgemm " << tg/1.0e6 << " s (" << flopGEMM/tg/1.0e6 << " TF/s)" << std::endl; + } + } + /////////////////////////////////////////////////////////////////////// // The SUMMA's conditions, one at a time, at 8 MB on ring B: // (a) LARGE persistent buffers (the rings use ~0.5 GB Abuf/Bbuf), sending // from offset 0 and from deep inside the region; diff --git a/systems/Frontier/schur2d_env.job b/systems/Frontier/schur2d_env.job new file mode 100644 index 000000000..fd974576a --- /dev/null +++ b/systems/Frontier/schur2d_env.job @@ -0,0 +1,91 @@ +#!/bin/bash -l +#SBATCH --job-name=schur2d-env +#SBATCH --nodes=36 +#SBATCH --ntasks-per-node=8 +#SBATCH --cpus-per-task=7 +#SBATCH --gpus-per-node=8 +#SBATCH --time=0:40:00 +#SBATCH --account=phy157_dwf +#SBATCH --gpu-bind=none +#SBATCH --exclusive +#SBATCH --mem=0 +#SBATCH -S 0 + +############################################################################## +# Environment walk for the 2D Schur inverse (2026-08-27). +# +# The SAME inverse (N=138240, 18x16, nb 480, SCHUR2D_LEAF_SPAN=9, LU leaf) +# runs 20.8 s in Test_schur2d_vs_slate and 26.5 s in the production example, +# and the LOCAL GPU work is what differs: GEMM 2.47 vs 3.11 s, leaf inverse +# 0.53 vs 0.76 s for identical calls; rings about equal. The environments +# differ in (a) OMP_NUM_THREADS 1 vs 7, (b) MPI_THREAD_MULTIPLE init (SLATE +# harness) vs SERIALIZED (example), (c) ~10.6 GB of fine-grid state resident +# on the device in the example. Test_schur2d_scale is the SLATE-free harness +# (same library code, plain Grid_init); walk it from the harness environment +# to the example's one knob at a time. +# +# Readout per cell (identical work in every cell): +# "gemm min/max" vs 2.45/2.61 (harness) / 3.01/3.24 (example) +# "BIG LEAVES ... inverse (max over ranks)" vs 0.53 / 0.76 +# "SCHUR fp64 ... invert" / "GRID : ... invert" vs 20.8 / 26.5 +# "Schur2D PROBE leaf W=" one-GCD getri / getrf+getrs / zgemm TF/s +############################################################################## + +cat << EOF > select_gpu +#!/bin/bash +export GPU_MAP=(0 1 2 3 7 6 5 4) +export NUMA_MAP=(3 3 1 1 2 2 0 0) +export GPU=\${GPU_MAP[\$SLURM_LOCALID]} +export NUMA=\${NUMA_MAP[\$SLURM_LOCALID]} +export HIP_VISIBLE_DEVICES=\$GPU +unset ROCR_VISIBLE_DEVICES +if [ \$SLURM_PROCID = "0" ]; then echo \$*; fi +exec numactl -m \$NUMA -N \$NUMA \$* +EOF +chmod +x ./select_gpu + +root=$HOME/ParallelIO/systems/Frontier +source $root/sourceme-rocm7.2.sh + +export MPICH_GPU_SUPPORT_ENABLED=1 +export MPICH_SMP_SINGLE_COPY_MODE=CMA +export MPICH_OFI_NIC_POLICY=GPU +module load libfabric + +BIN=$root/tests/debug/Test_schur2d_scale +OPTS1="--accelerator-threads 8 --shm 4096 --shm-mpi 1 --device-mem 32000" +vol=48.48.48.96 +MPI_GEOM=3.6.4.4 + +export S2D_N=138240 +export S2D_NB=480 +export SCHUR2D_LEAF_SPAN=9 +export SCHUR2D_LEAF_LU=1 +export SUMMA_HANDSHAKE=1 +# SCHUR2D_PROBE unset: probe on, including the one-GCD leaf microbenchmark + +run_cell () { + name=$1; shift + echo "----- $name : $* -----" + fname=log.env.$name + env "$@" srun -N36 -n288 --kill-on-bad-exit=1 ./select_gpu $BIN --mpi ${MPI_GEOM} --grid $vol $OPTS1 > $fname 2>&1 + echo " exit $?"; sleep 30 + grep -h "PROBE banner\|device ballast\|PROBE leaf\|BIG LEAVES\|gemm min/max\|invert\|certificate" $fname \ + | sed 's/^Grid : Message : [0-9.]* s : //' | cut -c1-170 +} + +# name environment +run_cell E1_omp1 OMP_NUM_THREADS=1 +run_cell E2_omp1_multiple OMP_NUM_THREADS=1 GRID_MPI_THREAD_MULTIPLE=1 MPICH_MAX_THREAD_SAFETY=multiple # = the SLATE harness environment +run_cell E3_omp7 OMP_NUM_THREADS=7 # the example's thread count +run_cell E4_omp7_ballast OMP_NUM_THREADS=7 S2D_BALLAST_GB=10.6 # + the example's device residency +run_cell E5_omp1_ballast OMP_NUM_THREADS=1 S2D_BALLAST_GB=10.6 # ballast alone +run_cell E6_omp1_again OMP_NUM_THREADS=1 # repeat E1: the noise floor + +echo "=========================================================" +echo "summary" +for f in log.env.E*; do + echo "$f: $(grep -h "invert" $f | grep -h "TOTAL\|took" | sed 's/^Grid : Message : [0-9.]* s : //' | cut -c1-80 | head -1) gemm $(grep -oh "gemm min/max [0-9./]*" $f) leaf $(grep -oh "inverse (max over ranks) [0-9.]*" $f)" +done +echo "harness (SLATE build, OMP=1, MULTIPLE): invert 20.8 gemm 2.45/2.61 leaf 0.53; example (OMP=7, SERIALIZED, 10.6 GB resident): 26.5 3.01/3.24 0.76" +echo "=========================================================" diff --git a/tests/debug/Test_schur2d_scale.cc b/tests/debug/Test_schur2d_scale.cc index 91fda5650..609e76c88 100644 --- a/tests/debug/Test_schur2d_scale.cc +++ b/tests/debug/Test_schur2d_scale.cc @@ -76,6 +76,20 @@ static ComplexD Fill(int64_t i, int64_t j, int64_t N) int main(int argc, char **argv) { + // Environment walk (2026-08-27): the SAME inverse runs 20.8 s in the SLATE + // harness and 26.5 s in the example, with the local GPU work (GEMM, leaf + // inverse) 25-40% slower in the example. Knobs to reproduce the example's + // environment here, one at a time (systems/Frontier/schur2d_env.job): + // GRID_MPI_THREAD_MULTIPLE=1 MPI_Init_thread(MULTIPLE) before Grid_init (SLATE harness does this) + // S2D_BALLAST_GB=x x GB of Lattice fields made device-resident before the invert + // (the example carries ~10.6 GB of fine-grid state in the MemoryManager) + // OMP_NUM_THREADS set in the job, read by nothing here but the runtime + if ( getenv("GRID_MPI_THREAD_MULTIPLE") ) { + int provided = 0; + MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &provided); + std::cout << "GRID_MPI_THREAD_MULTIPLE: requested MPI_THREAD_MULTIPLE, provided " << provided + << (provided==MPI_THREAD_MULTIPLE ? " (MULTIPLE)" : " (NOT multiple)") << std::endl; + } Grid_init(&argc, &argv); GridCartesian *grid = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), @@ -125,6 +139,24 @@ int main(int argc, char **argv) BlockCyclicMatrix A0(grid,N,nb,Pr,Pc); BlockCyclicSchurInverse RSI2; + // Device ballast: Lattice fields written on the accelerator so they sit in + // the MemoryManager's device LRU exactly as the example's fine-grid state does. + typedef Lattice,Ns> > BallastField; + std::vector ballast; + if ( getenv("S2D_BALLAST_GB") ) { + double gb = atof(getenv("S2D_BALLAST_GB")); + uint64_t fbytes = (uint64_t)grid->oSites()*sizeof(BallastField::vector_object); + int nf = (int)(gb*1.0e9/(double)fbytes + 0.5); + ballast.reserve(nf); + for(int i=0;ioSites(), 1, { v[ss] = Zero(); }); + } + std::cout << GridLogMessage << "Test_schur2d_scale: device ballast " << nf << " fields x " << fbytes/1.0e6 + << " MB = " << nf*fbytes/1.0e9 << " GB resident (S2D_BALLAST_GB=" << gb << ")" << std::endl; + } + BlockCyclicRedistribute::RowsToCyclic(grid,rowStart,&rows1d[0],myrows,A); double t2 = usecond(); if ( A.data.size() )