{"id":"a025ad32-de80-482e-bdad-e8ddbd03d629","arxiv_id":"2607.08743","paper_version":1,"verdict":"ACCEPT","confidence":"HIGH","novelty_score":4.0,"correctness_risk":"unknown","formal_verification":"none","parameter_count":8,"one_line_summary":"A Kokkos-based C++ port of the SACRA numerical relativity code achieves ~10x speedup on GPU/APU over the Fortran CPU version while preserving waveform accuracy, pi-symmetry, and second-order convergence.","lead":"SACRA-K is a new version of the SACRA numerical relativity code, rewritten in C++ with the Kokkos library so it runs efficiently on GPUs and APUs from different vendors. It reproduces the physics of the original Fortran code while running about 10x faster, which matters for simulating gravitational-wave sources like merging neutron stars.","discovery_kind":"unclear","skeptic_critique":{"model":"glm-5.2","headline":"No significant objection identified. The port is well-validated and the scaling limitations are honestly reported.","rationale":"The reader's verdict of ACCEPT with HIGH confidence is appropriate. The paper is a well-executed engineering contribution: it ports an established NR code to a performance-portable framework, validates it against the parent code across three binary types with quantitative acceptance criteria, demonstrates bitwise symmetry preservation, confirms expected convergence order, and measures speedup on two accelerator architectures. The reader correctly identified the scaling limitation as the weakest point, but the paper reports it transparently and the weak scaling results (77–83% at 256 devices relative to a 32-device baseline) are reasonable for a stencil code with fixed ghost-cell overhead. The absence of a public code release is a legitimate reproducibility concern but does not affect the validity of the claims as presented. The abstract's statement about waveform discrepancies is slightly overgeneralized for the BNS post-merger phase, but the body text qualifies this properly and the inspiral agreement—which carries the physics relevant for current GW detectors—is solidly established. I agree with the reader's assessment and see no reason to adjust the verdict.","tokens_in":22376,"tokens_out":2201,"duration_ms":99845,"concrete_test":"Run a code-to-code waveform comparison (SACRA-K vs SACRA-MPI) for the BNS configuration at N=80 or N=100 instead of N=40, extracting at the same r₀=480M⊙. If the inspiral phase difference |Δφ_GW| remains below 10⁻² rad and the relative amplitude difference stays below 10⁻³, the N=40 validation is confirmed as representative. If the residuals grow substantially at higher resolution, a resolution-dependent porting artifact would be revealed that the current N=40 comparison could mask.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The paper's central claims—correctness of the port, bitwise π-symmetry, second-order convergence, and ~10x speedup—are each supported by appropriate evidence. The reader's concern about the speedup being measured only for the smallest configuration is valid but not load-bearing: the paper transparently reports both strong scaling degradation (to ~10% efficiency at 256 devices) and weak scaling (77–83% relative to a 32-device baseline at 256 devices), and the abstract correctly qualifies the speedup as applying to 'the smallest test configuration.' One minor overgeneralization exists: the abstract states waveform discrepancies 'remain at or below the distinguishability threshold of current gravitational-wave detectors' without qualifying that for BNS post-merger, phase differences grow to O(1) rad (Section 3.1.3, Figure 4). The body text is honest about this, attributing it to chaotic amplification of roundoff-level differences in the differentially rotating remnant, and argues correctly that post-merger signal is not the primary target for current detector waveform modeling. This does not undermine the central claim that SACRA-K is a correct and substantially faster port. The code-to-code comparison at N=40 (Section 3.1) is low resolution, but the independent convergence test at N=92–188 (Section 3.3) provides a complementary check at production-relevant resolutions and shows the expected second-order behavior, matching SACRA-MPI. No single assumption or measurement is both critical to the central claim and insufficiently supported.","agreement_with_reader":"agree"},"referee_report":{"model":"glm-5.2","summary":"The manuscript presents SACRA-K, a C++ port of the Fortran numerical relativity code SACRA-MPI, built on the Kokkos performance-portability library. The port retains the physics modules of the parent code (BSSN with Z4c constraint damping, Berger-Oliger AMR, HRSC hydrodynamics from NANASI) while enabling execution on CPUs, NVIDIA GPUs, and AMD APUs from a single source. The authors validate the port through code-to-code waveform comparisons (BBH, BHNS, BNS), a bitwise pi-symmetry test, and a four-resolution self-convergence test for BNS. Performance is benchmarked on A100 and MI300A hardware, showing roughly an order-of-magnitude speedup over the Fortran code on a 72-core CPU node in the smallest configuration, with strong and weak scaling measured up to 256 devices.","tokens_in":23147,"tokens_out":1157,"duration_ms":171834,"significance":"Performance-portable numerical relativity codes are of growing importance to the field as computing architectures diversify and next-generation detectors demand higher-fidelity waveforms. SACRA-MPI is an established code with a strong track record in BNS and BHNS science; porting it to GPU/APU architectures via Kokkos is a valuable contribution. The validation suite is thorough and the acceptance criteria are well-motivated, grounded in inter-code scatter, resolution-dependent variation, and detector distinguishability thresholds. The honest reporting of scaling limitations (strong scaling efficiency dropping to ~10% at 256 devices due to fixed ghost-buffer overhead) is commendable. The pi-symmetry test, showing exact bitwise preservation to the roundoff floor, is a particularly strong correctness check on the Kokkos parallelization. The work is a solid engineering contribution that enables future science campaigns on heterogeneous architectures.","major_comments":[],"minor_comments":[{"comment":"Abstract: The statement that waveform discrepancies 'remain at or below the distinguishability threshold of current gravitational-wave detectors' is not qualified for the BNS post-merger phase, where phase differences grow to O(1) rad (Section 3.1.3, Figure 4). The body text is transparent about this and correctly argues that post-merger signal is not the primary target for current detector waveform modeling, but the abstract should be qualified to avoid overgeneralization.","section":null},{"comment":"Section 3.1.1: The BBH extraction radius r_0 = 24 M_sun is acknowledged as insufficiently large for scientific runs. The authors argue that the finite-radius bias is 'largely common mode' between the two codes. This is reasonable, but a sentence noting that the measured discrepancy may underestimate the total error with respect to the asymptotic waveform would strengthen the discussion.","section":null},{"comment":"Section 3.1.2: The BHNS consistency test meets the acceptance criteria only marginally, with |ΔA|/A riding the 10^-3 line in late inspiral and |Δφ_GW| touching 10^-2 rad just before merger. The authors attribute this to the sensitivity of finite-volume reconstruction at the stellar surface to floating-point operation ordering. This is plausible, but a brief quantitative comparison to the inter-resolution spread at similar configurations would help contextualize whether this marginality is a concern for production runs at higher N.","section":null},{"comment":"Section 4.1, Table 1: The zone-cycles-per-second metric is defined in the text but the exact zone-cycle count formula (number of level evolutions times active cells per level) could be stated more explicitly, perhaps as an equation, to aid reproducibility of the benchmark.","section":null},{"comment":"Section 4.2: The strong scaling degradation to ~10% efficiency at 256 devices is attributed to the growing buffer-to-active-cell ratio. The 6-cell ghost buffer width is set by the stencil requirements (fourth-order finite differencing plus KO dissipation). It would be useful to note whether a smaller buffer is feasible at lower stencil order or whether this is a fundamental constraint of the current numerical scheme.","section":null},{"comment":"Section 2.4: The statement that the entire evolution is carried out on the device through Kokkos kernels is clear, but the specific optimization for the AMD APU (team scratch for 8^3 tiles) is described only briefly. A sentence or two on whether this optimization is also applicable to NVIDIA GPUs (or why it is APU-specific) would help readers understand the portability of the optimization strategy.","section":null},{"comment":"Figure 2: The y-axis label for the bottom panel uses |Δφ_GW| while the text uses both |Δφ_GW| and |Δφ|. Standardizing the notation (the text uses |Δϕ_GW| with a phi symbol in some places and |Δφ_GW| with a varphi in others) would improve readability.","section":null},{"comment":"Section 2.2: The sentence describing the fallback to HLLE solver 'wherever the intermediate state would be unphysical' could briefly clarify how frequently this fallback occurs in practice for the test configurations, as it affects the effective order of the hydrodynamics scheme (not stated explicitly).","section":null}],"recommendation":"minor_revision","confidential_remarks":"The paper is a straightforward and well-executed code port with honest validation. The reader's concern about the speedup being measured only for the smallest configuration is valid but not load-bearing: the abstract correctly qualifies the claim, and the weak scaling results (77-83% at 256 devices relative to a 32-device baseline) provide reasonable evidence that production-scale runs will benefit from GPU acceleration. No novelty or citation concerns. The work fits well within the journal's scope for computational astrophysics methods."},"author_rebuttal":null,"desk_editor":{"model":"glm-5.2","letter":"This is a solid engineering paper. Han, Kiuchi, and Shibata ported the Fortran SACRA-MPI numerical relativity code to C++/Kokkos, and the validation work is thorough. They compare waveforms directly against the parent Fortran code across BBH, BHNS, and BNS configurations, show bitwise π-symmetry preservation, and demonstrate second-order convergence in a four-resolution BNS test matching prior SACRA-MPI results. The ~10x speedup over a 72-core CPU node is measured on both NVIDIA A100 and AMD MI300A, and the cross-vendor portability is the real selling point — one source code, two accelerator ecosystems, no vendor lock-in. The code-to-code waveform agreement is genuinely impressive: phase differences of 10^-5 rad during BBH inspiral, 10^-3 rad for BNS inspiral. The acceptance criteria are well-motivated, tied to actual detector distinguishability thresholds. The π-symmetry test compiling with -ffp-contract=off to prevent FMA breaking symmetry is a nice touch that shows they understand the subtleties. The soft spots are real but proportionate. The headline 10x speedup is measured only on the smallest configuration (4 accelerators, 2×2×1 decomposition). Strong scaling degrades badly — down to ~10% efficiency at 256 devices — because the fixed 6-cell ghost buffer per subdomain grows as a fraction of active cells. The paper is transparent about this, and the weak scaling results (77–83% at 256 devices relative to a 32-device baseline) are more encouraging, but the abstract's framing of the speedup could be clearer about this scope. One minor overgeneralization: the abstract says waveform discrepancies remain at or below detector thresholds, but BNS post-merger phase differences grow to O(1) rad. The body text is honest about this — chaotic amplification of roundoff differences in the turbulent remnant — and correctly notes post-merger signal isn't the primary target for current detector modeling. Still, the abstract overstates slightly. No public code release, which limits reproducibility. The validation at N=40 for code-to-code comparison is low resolution, but the independent convergence test at N=92–188 compensates. This is for the numerical relativity and computational astrophysics community — people who need GPU-accelerated compact binary simulations and care about cross-vendor portability. It deserves a serious referee. The scaling limitations are honestly reported and don't undermine the central claim that this is a correct, substantially faster port.","headline":"Validated Kokkos port of SACRA-MPI to GPU/APU with ~10x speedup; honest about scaling limits","tokens_in":23414,"tokens_out":596,"would_cite":true,"duration_ms":101426,"reading_group":"no","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":[],"pacs":["04.25.D-","04.30.-w","04.40.Dg","95.30.Lz"],"model":"glm-5.2","headline":"Relativity code ported to GPU runs 10x faster, matching original","keywords":["numerical relativity","performance portability","Kokkos","gravitational waves","binary neutron stars","adaptive mesh refinement","GPU acceleration","BSSN formulation"],"falsifier":"If the waveform discrepancies between SACRA-K and SACRA-MPI were found to exceed the distinguishability threshold of current gravitational-wave detectors during the inspiral phase, or if the pi-symmetry were broken at the bitwise level, or if the BNS convergence order fell below second order, the port would be deemed physically invalid. On the performance side, if the per-device throughput advantage over CPU did not hold at production-relevant resolutions or if weak scaling efficiency collapsed below 50% at moderate device counts, the practical utility claim would be undermined.","tokens_in":22632,"feed_emoji":"🚀","tokens_out":1250,"duration_ms":261385,"temperature":0.7,"pith_summary":"The paper presents SACRA-K, a rewrite of the Fortran numerical-relativity code SACRA-MPI into C++ using the Kokkos library, which allows a single source code to run on CPUs, GPUs, and APUs from different vendors. The authors validate that the port reproduces the original code's gravitational-wave waveforms across binary black hole, black hole–neutron star, and binary neutron star configurations to accuracy at or below the distinguishability threshold of current gravitational-wave detectors, preserves pi-symmetry at the bitwise level, and achieves second-order convergence for binary neutron star mergers. In the smallest test configuration, SACRA-K runs roughly an order of magnitude faster on NVIDIA A100 or AMD MI300A accelerators than the Fortran code on a 72-core CPU node, and the authors measure weak scaling up to 256 accelerator devices with efficiencies of 77–83% relative to a 32-device baseline.","feed_headline":"Relativity simulation code ported to GPU runs 10x faster","feed_subtitle":"A C++ rewrite of an established Fortran code matches its physics across three binary types while gaining an order of magnitude in speed on N","key_machinery":"Kokkos performance-portability library (abstracts CUDA, HIP, OpenMP backends); BSSN formulation with Z4c constraint damping for spacetime evolution; box-in-box Berger-Oliger adaptive mesh refinement with time subcycling; finite-difference (4th-order centered/upwind) for geometry and finite-volume PPM/HLLC for hydrodynamics; ghost-buffer communication via GPU-aware MPI; team-scratch memory optimisation on AMD APU","core_discovery":"The central claim is that an established, physics-rich numerical-relativity code (SACRA-MPI) can be ported to a performance-portable C++/Kokkos framework (SACRA-K) without sacrificing physical fidelity, gaining approximately an order of magnitude in throughput on modern accelerators. The correctness of the port is established by three independent checks: waveform agreement with the parent code at or below detector-distinguishability thresholds, exact bitwise preservation of pi-symmetry, and second-order convergence in binary neutron star merger simulations. The performance gain is demonstrated per-device in the smallest configuration and the scaling behaviour is characterised across up to 26","pith_inferences":["The weak scaling results (77–83% at 256 devices relative to a 32-device baseline) suggest that production-scale simulations using hundreds of devices will achieve effective throughput well below the headline 10x speedup, though still substantially above CPU-only performance.","The strong-scaling bottleneck from the fixed 6-cell ghost buffer implies that reducing buffer width or adopting asynchronous communication overlap could yield disproportionate gains at high device counts, a direction the authors acknowledge but have not yet implemented.","If the second-order convergence ceiling is set by the approximate Riemann solver and stellar surfaces, then higher-order hydrodynamics schemes would be needed to improve convergence beyond what is reported, independent of the hardware port.","The exact bitwise pi-symmetry preservation under Kokkos parallel reductions (with -ffp-contract=off) indicates that reproducibility of symmetry properties on GPUs is achievable but requires compiler-flag-level control, which may not transfer automatically to other codes."],"forward_implications":["If the port fidelity holds at production resolutions, gravitational-wave waveform libraries for binary neutron star and black hole–neutron star systems can be generated at roughly one-tenth the previous computational cost, enabling broader parameter-space surveys.","The Kokkos-based approach demonstrated here provides a template for porting other established Fortran numerical-relativity codes to heterogeneous architectures without rewriting the underlying physics.","The observed strong-scaling degradation (10% efficiency at 256 devices) signals that future production runs will need either larger per-device workloads or algorithmic changes to ghost-buffer communication to maintain throughput.","The post-merger phase divergence between SACRA-K and SACRA-MPI suggests that high-precision post-merger waveforms may be intrinsically sensitive to floating-point-level differences, a limitation that affects all GPU-ported NR codes, not just this one."],"fun_headline_variants":["Numerical relativity code ported to Kokkos runs 10x faster on GPUs","SACRA-K preserves physics fidelity while gaining 10x speed on accelerators","Fortran relativity code rewritten in C++ matches waveforms at 10x speed","GPU-portable numerical relativity code validated across three binary types","Cross-vendor relativity code keeps bitwise symmetry and gains 10x throughput"],"cache_read_input_tokens":0,"weakest_assumption_plain":"The headline 10x speedup is established only for the smallest configuration (4 accelerators), while production-scale simulations requiring many more devices face a strong-scaling efficiency that drops to about 10% at 256 devices, driven by the growing fraction of ghost-buffer cells relative to active cells as subdomains shrink. The claim that GPU acceleration enables large parameter-space surveys depends on weak-scaling behaviour that, while reasonable, is measured on a less ","fun_headline_variants_meta":{"raw":{"variants":["Numerical relativity code ported to Kokkos runs 10x faster on GPUs","SACRA-K preserves physics fidelity while gaining 10x speed on accelerators","Fortran relativity code rewritten in C++ matches waveforms at 10x speed","GPU-portable numerical relativity code validated across three binary types","Cross-vendor relativity code keeps bitwise symmetry and gains 10x throughput"]},"model":"glm-5.2","effort":"low","cost_usd":0.0,"raw_usage":{"total_tokens":725,"prompt_tokens":623,"completion_tokens":102,"prompt_tokens_details":null},"tokens_in":623,"tokens_out":102,"duration_ms":84194,"temperature":1.0,"reasoning_tokens":null,"cache_read_input_tokens":0,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-07-10T01:55:55.426735+00:00","model_set":{"reader":"glm-5.2"},"falsifier":"If the waveform discrepancies between SACRA-K and SACRA-MPI were found to exceed the distinguishability threshold of current gravitational-wave detectors during the inspiral phase, or if the pi-symmetry were broken at the bitwise level, or if the BNS convergence order fell below second order, the port would be deemed physically invalid. On the performance side, if the per-device throughput advantage over CPU did not hold at production-relevant resolutions or if weak scaling efficiency collapsed below 50% at moderate device counts, the practical utility claim would be undermined.","supporting_citations":[],"review_version":1}