Repository navigation
Optimize hot kernels and L-BFGS minimizer; fix correctness issues - #2
Merged
Merged
Conversation
Code review of the hot kernels and minimizer drivers: latent red-black race for odd n in rosso_krauth, redundant objective evaluation and reduction/sync overhead in lbfgs, duplicate edge sampling in gradient(), RNG pool contention in dynamics, plus build-config and convention hazards (line-tension mapping, default build type). https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
- rk_block: reject odd line lengths up front; with periodic boundaries sites 0 and n-1 share red-black parity but are nearest neighbours, so one colour pass would update both concurrently. Add a test. - depinning_ramp: integer ramp counter instead of accumulating a floating-point loop variable (bit-identical for the default step). - CMake: default CMAKE_BUILD_TYPE to Release; rename the placeholder project name to qew. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
phi at the current iterate equals the accepted trial value from the previous line search, so recomputing it each iteration was one full objective evaluation (a noise sample per site + reduction) wasted per step -- in the common first-trial-accepted case, half the objective cost. Only the first iteration evaluates it explicitly. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
objective() is refactored onto a shared body that reads the line position through a device-callable functor, with a new overload evaluating E at x0 + alpha*d directly inside the reduce. The Armijo loop in lbfgs_minimize no longer writes h per trial, so the saved x0 copy and the failed-search undo copy disappear; the accepted step is committed with a single axpy. One analytic definition of the energy is shared by both paths, keeping objective/gradient consistency intact. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
- s/y pair construction and the sy/ss/yy dot products collapse from one parallel_for plus three parallel_reduces into a single fused reduce (one pass over the data, one device sync instead of four). - The accepted pair is std::swap'ed into the history slot and g/g_new are swapped, replacing three full-vector deep_copies per iteration (safe: the FFT plans are bound to q and d only). https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
The per-kernel timings cannot see minimiser-loop overheads (line-search objective evals, two-loop recursion, history updates), which is where the lbfgs optimisations land. Measured on this machine (Serial, n=4096, arclength, FilteredNoise): 5.73s/12893 fevals before the lbfgs changes, 4.77s/10424 fevals after, at an identical 2470 iterations. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
The paper (arXiv:2410.21838, eq. (10)) controls the model through the
dimensionless ratio lambda_p/xi = (Gamma/U0)^{2/3}, varied by scanning
the line tension at fixed disorder. Previously the same config field
name meant that ratio in static_sweep/depinning_ramp but a physical
length (divided by xi internally) in dynamic_anneal.
All drivers now take pinning_length_over_xi and share
Gamma = amplitude * (lambda_p/xi)^{3/2},
pinned by a mapping assertion per driver. Line tensions are numerically
unchanged in every run. HDF5 key strings are untouched; the values carry
the ratio, which already held for static/depinning output (the dynamic
demo's recorded value changes 0.01 -> 0.1 accordingly).
https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
- Replace integer % with compare-and-adjust wraps in the bicubic sample (8 mods per gather), the model stencils, and the Rosso-Krauth sweep; the gx == nx rounding edge case is folded back explicitly. Integer division is expensive on GPUs and these are the hot loops. - DeviceNoise grows sample_value() (objective path: no dw weights, no derivative accumulators) and sample_dy() (gradient/RK path: dV/dx was computed everywhere and consumed nowhere). Kernels route through requires-detected helpers, so analytic test noises that only provide sample() keep working unchanged. Arithmetic is unchanged; the golden noise-interpolation and golden L-BFGS cross-checks pass. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
Each noise edge midpoint was bicubic-gathered twice, once by each adjacent site. gradient() now runs an edge kernel (one dV/dy gather and, for arclength, one normalized slope per edge) into a GradientWorkspace, followed by a cheap per-site combine. Edge values and the combination order are unchanged (FP addition commutes), so results are identical; the golden tests and the unchanged L-BFGS trajectory (2470 iters / 10424 fevals) confirm it. lbfgs/fire/LangevinDynamics own persistent workspaces; a workspace-free overload keeps the old call signature for tests. Serial benchmark (n=262144): grad 33.5 -> 26.5 ms, dyn 52.4 -> 34.4 ms, end-to-end preconditioned L-BFGS (n=4096) 4.77 -> 3.55 s. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
All moves are forward-only, so mean_h = mean_h0 + cumulative_advance/n; the per-sweep mean-h parallel_reduce (one of two reductions per sweep, with sweeps up to 2e5) is replaced by a host-side accumulation of the advance sum that is computed anyway. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
One pool state per contiguous lane (chunk of ~32 sites) instead of a locked get_state/free_state per site, the Kokkos-recommended pattern. Serial CPU timing is unchanged within noise (the kick is Box-Muller math there, with uncontended locks); the motivation is wide-GPU lock serialisation, to be verified on the CUDA runner. The zero-temperature path now skips the RNG branch entirely. Draw order changes, so thermal trajectories differ realization-wise; the equipartition and zero-diffusion determinism tests pass. https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd
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
This PR implements a comprehensive optimization and correctness-fix program for the qEW codebase, addressing both physics consistency and performance bottlenecks in the minimizer and force-kernel hot paths. All changes preserve the documented analytic consistency between
objective()andgradient(), and the golden-file test contracts remain unchanged.Key Changes
Phase 0: Correctness & Build Hygiene
rosso_krauth.h): Reject odd line lengths withstd::invalid_argumentsince periodic red-black sweeps have a data race when n is odd (sites 0 and n-1 share parity but are nearest neighbours). Add death test.qew_runs.h): Replace floating-point loop counter in force ramp with integer index to eliminate accumulated rounding errors in HDF5 keys.CMAKE_BUILD_TYPE=Releasewhen unset (with status message); rename project toqew.Phase 1: Physics Convention (Line Tension)
static_sweep,depinning_ramp,dynamic_anneal). All now use the dimensionless Larkin ratiolambda_p/xi(from arXiv:2410.21838 eq. (10)) and share the formulaGamma = amplitude * (lambda_p/xi)^{3/2}.pinning_length_over_xiand add per-driver unit tests pinning the mapping.Phase 2: L-BFGS Overhead Reduction
phi = objective(h, ...)above the iteration loop; reuse the accepted trial value across iterations, eliminating one full objective evaluation per iteration in the common case (first trial accepted).objective(x0, d, alpha, p, noise)overload that evaluates atx0 + alpha*dwithout materializing it. Commithonly on acceptance, removing one kernel launch andnstores per trial, plus thedeep_copyundo path.parallel_reduces (sy,ss,yy) into one with a 3-component reducer. Replacedeep_copyhistory updates withstd::swap(safe: FFT plans bound toq/donly).Phase 3: Hot-Kernel Work Reduction
filtered_noise.h,qew_model.h,rosso_krauth.h): Replace(i±1) % nstencil wraps with conditional add/subtract; clamp thegx == nxrounding edge case explicitly before computing fractions. Same arithmetic, bit-identical results.DeviceNoise::sample()intosample_value()(no derivatives) andsample_dy()(y-derivative only). Addbspline_weights_value()helper. Route calls via C++20requireshelpers (noise_value,noise_dvdy) that fall back to fullsample()for analytic test stubs. Objective uses value-only path; gradient uses dy-only.dv_dyand (for arclength) normalized slope per edge into scratch views; combine kernel uses these to build per-site gradient. Eliminates 2× bicubic gather duplication (one gather per edge instead of two). IntroduceGradientWorkspacestruct; keep back-compat overload allocating internally.Phase 4: Rosso–Krauth & Dynamics
https://claude.ai/code/session_01QE8Wd5igjty82gh7GvPuVd