Index-based mlebabeclap_gsrb kernel - #5620
Merged
WeiqunZhang merged 1 commit intoSep 2, 2026
Merged
Conversation
WeiqunZhang
approved these changes
Sep 2, 2026
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
Rewrites
mlebabeclap_gsrb(2D and 3D) from aBox-based kernel that loops internallyinto a per-cell
(i,j,k,n)kernel, and calls it fromMLEBABecLap::FsmooththroughAMREX_HOST_DEVICE_PARALLEL_FOR_4Dinstead ofAMREX_LAUNCH_HOST_DEVICE_LAMBDA.In 3D this also removes the hand-written quadruple
forloop and the// amrex::Loop here causes gcc 8 to crash.workaround that went with it, since there isno longer a loop inside the kernel at all.
The arithmetic is untouched: only the loop wrapper is removed and the body re-indented.
The kernel already took its EB geometry as an
EBData, so no argument rework was neededbeyond dropping
Box const& box/int ncompin favour of(i,j,k,n). Red/blackGauss-Seidel only writes cells of one parity and only reads neighbours of the other, so
the switch from a sequential
amrex::Loopto a concurrentParallelForintroduces noordering dependence. Results are bitwise unchanged.
Additional background
This is part of the ongoing split of the stale WIP #4922 ("GPU specific kernels for
MLEBABecLap"), as requested there; it supersedes the
mlebabeclap_gsrbpart of #4922.The templating on
Tthat #4922 also introduced has been dropped, per the review commenton that PR, and #4922's change of
dhxfromm_b_scalar/(h[0]*h[0])tom_b_scalar*dxinv[0]*dxinv[0]has deliberately been left out so that the answers staybit-for-bit the same.
This PR is independent of the
mlebabeclap_adotx/mlebabeclap_adotx_centroidPRs(different kernels, different function in
AMReX_MLEBABecLap_F.cpp) and can be merged inany order relative to them. The
(i,j,k,n)form is the prerequisite for a later PR thatadds a fused (
MultiFab-wideParallelFor) GPU path forFsmooth; that follow-up is notincluded here.
Testing (all commands run from the repo root unless noted):
cmake -S . -B build-split -DAMReX_SPACEDIM=3 -DAMReX_EB=ON -DAMReX_LINEAR_SOLVERS_EM=OFF -DAMReX_ENABLE_TESTS=ON -DAMReX_TEST_TYPE=Small -DAMReX_MPI=ON -DCMAKE_BUILD_TYPE=Releasethen
cmake --build build-split -j8andctest --test-dir build-split --output-on-failure→ builds clean, 9/9 tests pass. A 2D library build
(
-DAMReX_SPACEDIM=2 -DAMReX_EB=ON) also builds clean.Tests/LinearSolvers/CellEB,make -j8 COMP=llvm USE_MPI=FALSE DIM=3andDIM=2,run over seven configurations (
sphere,sphere+eb_is_dirichlet=1,rotated_box,two_spheres,flower, two-levelsphere, periodicsphere;n_cell=64,verbose=2). This test drives the V-cycle, soFsmoothis called on every level.The MLMG/BiCGStab residual histories and the initial/final max, 1- and 2-norm residuals
are bitwise identical to
developmentin both 2D and 3D.CUDA build on an NVIDIA RTX A5000 (CUDA 13.2, gcc 11.4):
Tests/LinearSolvers/CellEB,make -j8 COMP=gnu USE_MPI=FALSE USE_CUDA=TRUE CUDA_ARCH=86 DIM=3,same seven configurations. MLMG iteration counts are identical to
development(9, 12, 11, 3, 11, 28, 11); the residual values differ only at the level of the
run-to-run nondeterminism of the unmodified binary (GPU reductions are not
bit-reproducible), and all solves converge to
resid/resid0 ~ 1e-13.The same CellEB matrix was re-run on Linux/gcc 11.4 (
make -j8 COMP=gnu USE_MPI=FALSE,DIM=3andDIM=2) against adevelopmentbuild in a sibling worktree: residualhistories again bitwise identical in both dimensions.
Performance
Neither timing changed measurably; this PR is a refactor, not an optimisation.
CPU,
Tests/LinearSolvers/CellEBmain3d.gnu.TEST.ex(
make -j8 COMP=gnu USE_MPI=FALSE DIM=3, gcc 11.4, single rank),inputs n_cell=128 eb2.geom_type=sphere eb_is_dirichlet=1 verbose=1.All binaries were built first and then timed interleaved in the same session
(3 reps each) on a shared machine, so the absolute numbers are inflated but the
comparison is fair. Best of 3, MLMG
Timers: Solve[s]:developmentmax_grid_size=32max_grid_size=64GPU (NVIDIA RTX A5000, CUDA 13.2,
make -j8 COMP=gnu USE_MPI=FALSE USE_CUDA=TRUE CUDA_ARCH=86 DIM=3), same test atn_cell=256, again built up front and timed interleaved, best of 3:developmentmax_grid_size=32(512 boxes)max_grid_size=64(64 boxes)MLMG converged in the same number of iterations in every one of these runs.
Also worth noting for this one: the CPU path changes from a sequential
amrex::Loopinside the kernel toAMREX_HOST_DEVICE_PARALLEL_FOR_4D, andvlo/vhiare now recomputed per cell rather than per box. Neither shows up in the timings above.Checklist
The proposed changes:
P.S Generated using Claude Code