{"id":"16ca57f4-fc90-4745-b276-5ab84db14c63","arxiv_id":"2608.06978","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":7.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":2,"one_line_summary":"A new conservative SPMHD scheme in SWIFT passes standard tests and achieves the first coupling of the EAGLE galaxy formation model to magnetohydrodynamics.","lead":"This paper presents a new smoothed particle magnetohydrodynamics (SPMHD) solver implemented in the open-source SWIFT simulation code, designed to add magnetic fields to galaxy formation simulations. It validates the method on a suite of standard tests and three astrophysical applications, including the first coupling of the EAGLE galaxy formation model to an MHD solver.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The RMS plasma-beta estimator (Eq. 47) can suppress the tensile-instability correction in exactly the low-β clumps where Eq. 43 says it is needed; this production-stability claim needs a targeted disordered-two-phase test.","rationale":"The reader's weakest_assumption is the RMS beta_loc estimator, and after reading the method section this is also the point I would stress-test. The scheme's advertised production robustness hinges on Eq. (47) because the tensile-instability correction is the main guard against particle clumping in disordered, feedback-disturbed particle fields; all other dissipative and cleaning terms are standard or have direct validation in Sections 3.1-3.2. My concern is more specific than the reader's: the RMS average is not merely 'possibly unstable'; it has a concrete failure mode in which a low-beta particle is embedded in high-beta neighbours, so Eq. (47) exceeds 10, lambda=0, and the correction is switched off exactly where Eq. (43) identifies the unstable regime. The paper's own standard tests do not isolate this configuration, and the Alfven-wave test described as beta<1 in Section 3.1.2 actually has beta=20 for the stated P=0.1, B=0.1, mu0=1, so that validation point is weaker than claimed. Sections 3.3-3.5, which would contain the decisive full-physics evidence, are truncated. The proposed two-phase disordered-clump test is cheap, decisive, and directly settles whether Eq. (47) over-suppresses the correction. If it passes, the central stability claim is materially supported and the CONDITIONAL verdict can be upgraded; if it fails, the production-stability claim is not yet established. I therefore keep the reader's CONDITIONAL verdict unchanged.","tokens_in":58157,"tokens_out":6156,"duration_ms":64400,"concrete_test":"Run a 3D periodic box with a high-beta background (beta about 100) containing a spherical low-beta clump (beta about 0.1), using the fixed production configuration of Section 3 and deliberately disordered, slightly jittered particle positions to mimic feedback. Evolve for several sound-crossing times and compare three runs: (i) lambda computed from local beta_i (Eq. 46), (ii) lambda from RMS beta_loc_i (Eq. 47), and (iii) lambda forced to 1 inside the clump. If run (ii) shows tensile particle clumping or spurious voids inside the low-beta clump while run (iii) does not, the RMS estimator has disabled the correction in a regime where Eq. (43) requires it; if run (ii) matches run (iii) in stability and shows no clumping, the concern is retired.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central production-stability claim rests on the RMS beta_loc of Eq. (47) modulating the tensile-instability correction lambda(beta_loc) in Eqs. (44)-(45). Because Eq. (47) is an unweighted RMS over all neighbours, it is not a local estimate of the beta that enters the stability criterion (43). In a magnetically dominated clump (beta_i<1) embedded in a high-beta ambient medium, the RMS over the neighbourhood can exceed 10, setting lambda=0 and switching off the correction precisely where the linear-stability condition says the monopole force must be subtracted. The authors motivate the RMS choice by 'extensive experimentation' and by the need to avoid spurious repulsive forces in feedback-disordered particle arrangements, but they provide no analysis of the regime where the RMS and local beta disagree. The standard validation tests do not isolate this failure mode: the Alfven-wave test quoted as 'beta<1' in Section 3.1.2 has beta=2P/B^2=20 with the stated P=0.1, B=0.1, mu0=1 parameters, so it does not actually exercise the low-beta branch of lambda; the shock tubes reach beta<1 only on one side and with ordered lattices. The full-physics runs that would show whether Eq. (47) suppresses clumping or over-suppresses the correction are in the truncated Sections 3.3-3.5, so the load-bearing evidence is currently an empirical assertion rather than a demonstrated property.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript presents a new smoothed particle magnetohrodynamics (SPMHD) formulation implemented in the open-source SWIFT code. The scheme is built on the SPHENIX density-energy SPH solver and adds a direct-induction ideal MHD discretisation, a tensile-instability correction with a novel RMS-based local plasma-beta estimator (Eq. 47), a constrained hyperbolic/parabolic divergence-cleaning scheme, artificial resistivity with an Alfvén-speed signal, Ohmic diffusion, and a cosmological comoving formulation. The authors validate the method on a comprehensive suite of standard tests: Alfvén wave convergence, MHD shock tubes, Orszag-Tang vortex, magnetic rotor, blast wave, Kelvin-Helmholtz and cloud-wind interaction, plus non-ideal diffusion tests. They announce three astrophysical applications (proto-stellar jet launching, galaxy-cluster dynamo, and a Milky Way-like disk with the EAGLE model) as the final part of the test suite. The full text provided to the referee, however, breaks off in Section 3.2.3 and does not contain Sections 3.3–3.5; the announced production-scale applications are therefore not accessible for assessment.","tokens_in":58374,"tokens_out":3645,"duration_ms":40563,"significance":"If the claims hold, this is a valuable contribution: the scheme is implemented in a modern, massively parallel, open-source code; the hyperparameters are kept fixed across the test suite; the validation set is broad and mostly quantitative; and the Alfvén-wave and non-ideal diffusion tests show second-order or better convergence. A successful coupling of the EAGLE galaxy-formation model to an MHD solver would be a first and would open the door to production magnetised galaxy-formation simulations. The central production-stability claim, however, rests on a new regularisation ingredient whose failure mode is not isolated by the presented tests, and the evidence that would demonstrate the claim in the target regime is in the missing Sections 3.3–3.5. The manuscript is therefore scientifically promising but, in its current provided form, does not yet support the headline application-level claims.","major_comments":[{"comment":"The abstract and introduction announce three production-scale applications, culminating in 'the first reported coupling of the EAGLE galaxy formation model to a magnetohydrodynamics solver', but the provided manuscript text ends in the middle of Section 3.2.3 and never presents Sections 3.3, 3.4, or 3.5. The central claim of production viability, and specifically the stable operation of Eq. (47) in full-physics runs with sub-grid feedback, is therefore currently unsupported by any presented evidence. This is not a stylistic issue: the announced applications are the load-bearing demonstration that the method works in the regime for which it was designed. The revision must include these sections, with quantitative diagnostics (e.g., stability over time, divergence-error statistics, clumping checks) rather than only field maps.","section":"Sections 3.3–3.5"},{"comment":"The Alfvén-wave test text states that the chosen parameters 'translate into a β < 1' and that the test assesses the tensile-instability correction in the strong-field regime. With the stated parameters, P = 0.1, B = 0.1, and μ0 = 1, the plasma beta is β = 2μ0P/B^2 = 20, not β < 1. The test therefore does not exercise the β < 2 branch of the switch λ(β_loc) in Eq. (44), and the claim that the scheme is 'robust to particle clumping in the strong field regime' is not demonstrated by this experiment. Either the parameters must be changed so that the low-β branch is actually probed, or the text must be corrected to state that the test covers only the β > 10 regime of λ.","section":"Section 3.1.2, Eq. (44)–(45)"},{"comment":"The RMS-based estimator β_loc (Eq. 47) is an unweighted average over all neighbours and is not the local plasma beta that appears in the stability criterion (43). In a magnetically dominated clump (β_i < 1) embedded in a high-β ambient medium, the RMS over the neighbourhood can be ≫ 1, setting λ = 0 in Eq. (44) and switching off the tensile-instability correction precisely where the linear-stability condition requires it to be active. The authors motivate Eq. (47) by 'extensive experimentation' and by improved behaviour in disordered, feedback-stirred particle arrangements, but no analysis or targeted test is provided for the regime where the RMS and the particle-local beta disagree. The existing validation tests do not isolate this case: the shock tubes and blast wave use ordered or nearly ordered lattices, and the Alfvén wave (as noted above) actually runs at β = 20. The manuscript needs a targeted disordered two-phase test — e.g., a low-β clump in a high-β ambient medium with a randomised or feedback-like particle distribution — or a quantitative study of the β_loc distribution in the application runs. Without this, the central production-stability claim rests on an empirical assertion.","section":"Section 2.3.3, Eq. (47)"},{"comment":"The new time-step conditions introduced for B and ψ are not written as well-defined mathematical expressions. Eq. (74) divides a vector by an expression ending in a dot, and Eq. (75) contains a similar dangling '·' in the denominator; it is unclear what vector or scalar operation is intended. This matters because the authors state in Section 3 that these conditions are disabled by default in the presented tests, but the conditions are part of the method description and will be used by adopters. The equations need to be made explicit, including the definition of the norm or componentwise operation.","section":"Section 2.3.6, Eqs. (74)–(75)"}],"minor_comments":[{"comment":"In the description of the Kelvin-Helmholtz set-up, the text says the central region's physical attributes 'are denoted by the subscript O'; this should be the subscript C (and similarly in the following sentence).","section":"Section 3.1.8"},{"comment":"The caption lists the strong-field run as 'β = 250' in the fourth and fifth rows; from the body text and the simulation description, this should be β = 25.","section":"Section 3.1.9, Fig. 17 caption"},{"comment":"The text refers repeatedly to 'Karapiperis & Schaller (2025)' and 'Shchutskyi et al. (2025)' for validation of cosmological non-ideal MHD and dynamo benchmarks, but these are not presented in this manuscript; since the provided text also omits the application sections, the reader cannot verify the claimed cross-validation.","section":"Section 3.2.2"},{"comment":"The manuscript is unusually long and some figure captions repeat the full setup parameters verbatim from the text; shortening the captions and moving the fully detailed setup data to a table would improve readability.","section":"General"}],"recommendation":"major_revision","confidential_remarks":"The most serious issue is that the provided manuscript text breaks off before the application sections, so the paper's central novelty — the production-scale EAGLE-MHD coupling — is not assessable. I would ask the editor to confirm whether the submission was accidentally truncated; if the sections exist, they must be included in the revision. In addition, the manuscript cites two companion papers (Karapiperis & Schaller 2025; Shchutskyi et al. 2025) for key validation of cosmological and dynamo behaviour; the editor should check whether those works are in fact forthcoming and whether this paper can stand alone. The beta-misstatement in the Alfvén-wave section is a factual error that should be corrected at proof stage, but it is not, by itself, a reason for rejection."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Two things to know about arXiv:2608.06978. First, it is a serious, well-executed SPMHD methods paper: the scheme is conservative, built on SPHENIX, and validated on a broad suite of standard tests with honest reporting of limitations. Second, the production-stability story leans heavily on a new empirical RMS-based local beta estimator (Eq. 47) that modulates the tensile instability correction, and that estimator has a real failure mode that the paper does not analyze: in a magnetically dominated clump (beta < 1) surrounded by high-beta gas, the RMS over the neighbourhood can exceed 10 and switch the correction off exactly when the linear stability criterion says it is needed. The paper's full-physics evidence for this estimator is in the truncated Sections 3.3–3.5, so the claim currently rests on assertion rather than demonstration.\n\nThe genuinely new pieces are the comoving B parametrization (B_c = a^{3γ/2} B, which makes the Alfvén speed transform like the sound speed), the modified shock indicator, the choice of the Alfvén speed as the artificial resistivity signal velocity, and the first EAGLE-MHD coupling claim. The validation work is strong: second-order convergence on the Alfvén wave, good agreement with ATHENA on Orszag-Tang and the rotor, and analytic matches for the non-ideal diffusion tests. The authors also keep hyperparameters fixed across tests and ship code and parameter files, which earns credit.\n\nThe soft spots, in proportion: (1) the RMS beta_loc is the load-bearing empirical ingredient for production stability, and the paper gives no analysis of the regime where the RMS and the local beta disagree; a targeted disordered-two-phase test, or at least a quantification of the discrepancy, should be added. (2) The paper states the Alfvén wave test has beta < 1, but with P=0.1, B=0.1, μ0=1 the beta is 20, so that test does not actually exercise the low-beta branch of the correction. That is a minor but real misstatement. (3) The astrophysical application sections were truncated in my copy, so I could not verify the EAGLE coupling results; the paper's reach claim is strong, so the referee should look at those figures carefully.\n\nWho is this for? Galaxy formation simulators who want MHD in SWIFT, and SPMHD method developers. It deserves a serious referee. My recommendation: engage with it, but require a revision that addresses the RMS estimator's failure mode and corrects the beta claim.","headline":"Solid SPMHD methods paper with honest validation, but the production-stability claim rests on an under-analyzed RMS beta estimator that can fail in exactly the low-beta clumps it is meant to protect.","tokens_in":58980,"tokens_out":4691,"would_cite":true,"duration_ms":41243,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65M75","76W05","85-08"],"pacs":[],"model":"deepseek-v4-flash","headline":"This paper introduces a smoothed particle magnetohydrodynamics (SPMHD) formulation in the SWIFT code, and reports the first coupling of the EAGLE galaxy formation model to an MHD solver, with tests spanning a proto-stellar jet, a cluster…","keywords":["magnetohydrodynamics","smoothed particle hydrodynamics","galaxy formation","methods: numerical","magnetic fields","SWIFT","EAGLE model","cosmological simulations"],"falsifier":"Reproduce the EAGLE disk-galaxy application twice, once with the RMS plasma beta estimate (equation 47) and once with the naive per-particle beta (equation 46) driving the tensile-instability correction; the paper's claim predicts the naive version produces violent particle accelerations or ejections in feedback-disordered regions while the RMS version stays stable.","tokens_in":57893,"feed_emoji":"🧲","tokens_out":10746,"duration_ms":103076,"temperature":0.7,"pith_summary":"The paper sets out to make magnetic fields a standard component of production galaxy-formation simulations. It introduces a smoothed particle magnetohydrodynamics (SPMHD) formulation in the SWIFT code, built on the SPHENIX density-energy hydrodynamics scheme, with adaptive stabilisation: a tensile-instability correction whose strength is set by a new RMS-based local plasma $\\beta$ estimate, a constrained hyperbolic/parabolic divergence cleaner, and artificial resistivity tuned to limit excess diffusion in cosmological flows. The authors argue the scheme is stable, accurate, and cheap enough to couple to effective sub-grid galaxy-formation recipes, and they demonstrate it on a jet-launching proto-stellar core, a massive galaxy cluster dynamo, and a Milky Way-like disk galaxy. The disk run is the first reported coupling of the EAGLE galaxy formation model to a magnetohydrodynamics solver. If the central claim is right, magnetic physics can be included in large cosmological runs at production cost.","feed_headline":"Magnetic fields join the EAGLE galaxy formation model","feed_subtitle":"The scheme stays stable through feedback and reproduces jets, cluster dynamos, and disk magnetic fields.","key_machinery":"The load-bearing object is the SPMHD evolution system in SWIFT, expressed in a density-energy conservative form and paired with the SPHENIX discontinuity-capturing terms. Three mechanisms carry the stability argument: (i) a tensile-instability correction that subtracts a fraction of the monopole force, with the fraction set by the new estimator $\\beta_{\\rm loc} = \\sqrt{\\sum_j \\beta_j^2 / \\sum_j 1}$ (an unweighted RMS over neighbours, equation 47); (ii) a constrained hyperbolic/parabolic divergence-cleaning scalar whose cleaning speed is $c_h = v_{\\rm sig}/2$ and whose damping rate has $\\sigma_p = 1$ for critical damping; and (iii) artificial resistivity with signal velocity equal to the Alfvén speed and switch $\\alpha_{\\rm AR}=\\min(\\alpha_{\\rm max}, h\\|\\nabla\\mathbf{B}\\|/\\|\\mathbf{B}\\|)$, which is invariant under field rescaling. The design goal is that each term acts only where needed, limiting spurious dissipation while keeping particles stable in the disordered arrangements produced by sub-grid feedback.","core_discovery":"The paper's central claim is that a conservative density-energy SPMHD scheme, augmented with three adaptive regularisation ingredients, is stable, accurate, and computationally efficient enough for production galaxy-formation simulations. The ingredients are a tensile-instability correction modulated by an unweighted root-mean-square of neighbouring plasma betas rather than the particle's own $\\beta$; a constrained mixed hyperbolic/parabolic divergence-cleaning scheme with cleaning speed set to half the pairwise signal velocity and critical damping at $\\sigma_p=1$; and artificial resistivity with the Alfvén speed as signal velocity and a Tricco-Price shock indicator that is invariant under $\\mathbf{B}\\to\\lambda\\mathbf{B}$. With all hyperparameters fixed across the test suite, the method reproduces standard laboratory MHD experiments, launches a jet from a forming proto-stellar core, amplifies a magnetic field in a massive galaxy cluster, and evolves magnetic fields in a Milky Way-like disk galaxy, constituting the first reported EAGLE-MHD simulation.","pith_inferences":["We infer that the RMS plasma beta estimator, being less sensitive to disordered particle neighbourhoods, could also suppress spurious tensile corrections in other meshless MHD solvers that struggle with sub-grid feedback, though the paper only tests it in SWIFT.","A testable extension not explored in the paper is to compare long-term magnetic energy growth in cluster runs with and without the divergence cleaner's energy-conserving terms; the paper's energy accounting implies that growth should come from physical dynamo action rather than cleaning artefacts.","Because the divergence-cleaning speed is chosen as half the local signal speed rather than the fast magnetosonic speed, the scheme may be cheaper in cosmological boxes; a direct benchmark against the more common choice would quantify that gain.","The first EAGLE-MHD disk result suggests that magnetic fields can now be included as a standard module in future large-volume cosmological campaigns, with the main open question being the trade-off in computational cost as resolution and box size grow."],"forward_implications":["A galaxy-formation simulation with full sub-grid physics can now carry a magnetic field, so EAGLE-style runs gain a magnetised baseline for disk, cluster, and proto-stellar studies.","Because hyperparameters stay fixed across the entire test suite, the scheme's published configuration is a directly usable production setup rather than a per-problem tuned one.","The same solver handles laboratory shocks, jets, cluster dynamos, and disks, so a single code path can cover cosmological and object-scale MHD without switching methods.","Second-order convergence on smooth Alfvén waves and beyond-second-order convergence on Ohmic diffusion tests mean resolution studies with this scheme improve accuracy predictably.","The reported cluster and disk runs provide concrete magnetic field strengths and topologies that future cosmological MHD simulations can be compared against."],"supporting_citations":[{"why":"Baseline Direct Induction SPMHD scheme this work extends; supplies conservative equations of motion, divergence cleaning, and artificial resistivity prescriptions.","marker":"Price et al. (2018)"},{"why":"Constrained hyperbolic/parabolic divergence cleaning for SPMHD; provides the conjugate-operator energy-conserving formalism.","marker":"Tricco & Price (2012)"},{"why":"Variable cleaning speeds and the evolved quantity $\\psi/c_h$, which the paper adopts for divergence control.","marker":"Tricco et al. (2016a)"},{"why":"Origin of the tensile-instability correction by subtracting the monopole force, which the paper modulates with the local plasma beta estimate.","marker":"Børve et al. (2001)"},{"why":"Meshless MHD baseline and specific recommendations the paper follows for signal velocities, timestep constraints, and resistivity switch forms.","marker":"Hopkins & Raives (2016)"},{"why":"SPHENIX hydrodynamics scheme on which the SPMHD solver is built; supplies the density-energy equations of motion and discontinuity-capturing terms.","marker":"Borrow et al. (2022)"},{"why":"The EAGLE galaxy formation model that is coupled to MHD for the first time in this work.","marker":"Schaye et al. (2015)"},{"why":"Calibration details of EAGLE that the disk-galaxy application inherits.","marker":"Crain et al. (2015)"},{"why":"Mixed hyperbolic/parabolic divergence-cleaning method that underlies the cleaning approach used here.","marker":"Dedner et al. (2002)"}],"fun_headline_variants":["First EAGLE-MHD simulation: magnetic fields in galaxy formation","New SPMHD scheme brings magnetic fields to EAGLE simulations","Stable MHD for galaxy formation: from jets to cluster dynamos","Magnetic fields join galaxy formation simulations with EAGLE","Production-ready MHD scheme for cosmic structure formation"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that the new RMS-based local plasma beta estimate reliably tells the tensile-instability correction where to act, because if it fires in disordered, feedback-driven particle configurations the correction becomes a spurious repulsive force and the scheme becomes unstable.","fun_headline_variants_meta":{"raw":{"variants":["First EAGLE-MHD simulation: magnetic fields in galaxy formation","New SPMHD scheme brings magnetic fields to EAGLE simulations","Stable MHD for galaxy formation: from jets to cluster dynamos","Magnetic fields join galaxy formation simulations with EAGLE","Production-ready MHD scheme for cosmic structure formation"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000619,"raw_usage":{"total_tokens":2918,"prompt_tokens":1036,"completion_tokens":1882,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":652,"completion_tokens_details":{"reasoning_tokens":1810}},"tokens_in":652,"tokens_out":1882,"duration_ms":11015,"temperature":1.0,"reasoning_tokens":1810,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-10T17:14:01.063756+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Reproduce the EAGLE disk-galaxy application twice, once with the RMS plasma beta estimate (equation 47) and once with the naive per-particle beta (equation 46) driving the tensile-instability correction; the paper's claim predicts the naive version produces violent particle accelerations or ejections in feedback-disordered regions while the RMS version stays stable.","supporting_citations":[],"review_version":1}