REVIEW 3 major objections 7 minor 19 references
Axisymmetric simulations of vertical displacement events in tokamaks: A benchmark of M3D-C$^1$, NIMROD and JOREK
T0 review · 3 major / 7 minor · reviewed 2026-08-14 · deepseek-v4-flash
Pith's one-line read This paper establishes a reusable axisymmetric benchmark showing that three independent nonlinear MHD codes—M3D-C1, NIMROD, and JOREK—agree on vertical displacement event growth rates and wall currents.
desk verdict A solid first VDE benchmark among M3D-C1, NIMROD, and JOREK; the reduced-MHD caveat is real but does not sink the paper. read the letter →
The pith
A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.
The reading
What carries the argument
The load-bearing object is the benchmark case itself: a vertically unstable NSTX-like equilibrium inside an axisymmetric rectangular resistive wall, with defined temperature profiles, Spitzer resistivity, and heat and particle diffusion coefficients. The mechanism that explains the key agreement is the slip motion condition, v = -R²∇u×∇φ, which lets the plasma slide across the large toroidal field without doing work against it; when F=RBφ stays within 5% of the vacuum field, reduced MHD reproduces full MHD for the n=0 vertical instability. The codes' resistive-wall treatments differ—M3D-C1 embeds a thick wall in a single mesh, NIMROD uses a thin-wall vacuum coupling, and JOREK uses the STARWALL Green's function—so the agreement exercises the wall-response physics rather than a shared implementation.
What would settle it
Compute the same benchmark with a lower edge resistivity (higher edge temperature) so response currents form in the open-field-line region, and compare JOREK to M3D-C1 and NIMROD; if the growth rates diverge where halo response dominates, the slip motion justification is not general. Equivalently, a 3D nonlinear VDE run in which F deviates from F_vacuum by more than 5% before first wall contact would test whether reduced MHD remains valid exactly when asymmetric forces matter.
Extended reading notes
Core claim
The central claim is that three independently developed nonlinear MHD codes—M3D-C1, NIMROD, and JOREK—reproduce each other's axisymmetric VDE dynamics well enough to serve as mutual verification. In the low-wall-resistivity, cold-open-field-line regime, all three recover the expected linear dependence of the VDE growth rate on wall resistivity; deviations stay around 3% for most cases and below 15% overall. For the full nonlinear evolution, time-shifted traces of the magnetic-axis position, toroidal currents, and wall and halo currents agree closely between the codes, even though M3D-C1 and NIMROD use full MHD while JOREK uses a reduced MHD model. The authors attribute this agreement to the slip motion condition: for the n=0 instability the plasma moves across the large toroidal field with F≡RBφ within 5% of the vacuum value, which is precisely the regime where reduced MHD is valid.
Load-bearing premise
The claim that JOREK's reduced MHD model captures the vertical instability rests on the toroidal field remaining close to its vacuum value (within 5% here); if the field deviates more during the nonlinear evolution, the JOREK agreement with full-MHD codes could be coincidental.
Editorial extensions
If this is right
- The benchmark case, supplied with equilibrium and coil files, gives other VDE codes a standard test case for validating their linear growth rates and nonlinear evolution.
- In the low-wall-resistivity, cold-edge regime, VDE growth rates should scale linearly with wall resistivity; codes that show a different scaling are likely missing wall-response physics.
- For n=0 vertical stability, reduced MHD can stand in for full MHD when RBφ stays close to the vacuum toroidal field, enabling cheaper nonlinear VDE simulations.
- The three-code agreement on halo currents supports using axisymmetric simulations to estimate vessel forces during the early, symmetric phase of a VDE.
Reading between the lines
- A natural next test is to run the same benchmark without the artificial thermal quench, to see whether the agreement persists through scrape-off and current-profile evolution rather than only at the chosen quench time.
- If the slip motion condition degrades in future 3D simulations, JOREK may need the full-MHD STARWALL coupling for asymmetric VDEs; the 5% F deviation observed here gives a quantitative threshold to monitor.
- The benchmark's rectangular-wall geometry could be adapted to test other codes' treatment of halo current closure through conducting structures, since halo width is self-consistently determined rather than prescribed.
Signed reviews
Editorial analysis
A structured set of objections, weighed in public.
Referee Report
Summary. This paper presents an inter-code benchmark for axisymmetric vertical displacement event (VDE) simulations, comparing the MHD codes M3D-C1, NIMROD, and JOREK on a vertically unstable NSTX-like equilibrium enclosed by a rectangular resistive wall. The authors report linear VDE growth rates from M3D-C1 and NIMROD linear simulations and from the linear phase of axisymmetric nonlinear simulations of all three codes, recovering the expected linear scaling with wall resistivity in the appropriate regime. They then compare full nonlinear VDE evolutions (plasma position, toroidal currents, wall eddy currents, and halo currents) after applying a time shift to align first wall contact and after imposing an artificial thermal quench. The paper finds good (10-15%) agreement in linear growth rates and good qualitative agreement in nonlinear traces, and interprets the JOREK agreement in terms of the reduced MHD 'Slip Motion Condition' justified by a claimed maximum 5% deviation of F=RBphi from the vacuum value. The equilibrium, coil data, and parameter tables are provided as reusable benchmark input.
Significance. If the benchmark is valid, it provides a much-needed verification case for VDE simulations, which are critical for predicting disruption forces in ITER and future tokamaks. The paper is strong in its careful parameter specification, use of three independently developed codes, and candid discussion of caveats (artificial thermal quench, time shifting, and model differences). The supplementary files (geqdsk equilibrium, coil positions) make the case directly reusable by other groups, which is a valuable community contribution. However, the interpretation that JOREK's reduced-MHD results mutually verify the full-MHD codes rests on an unshown 5% F-vacuum proximity claim, and the time-shifting of nonlinear traces raises questions about what precisely is being compared. These issues do not invalidate the raw comparison but need to be addressed for the benchmark's cross-verification claim to be fully supported.
major comments (3)
- [Section 7] The statement that 'the total toroidal field differed from the vacuum toroidal field by a maximum of 5%' is central to the paper's claim that JOREK's reduced MHD model reproduces the full MHD models, yet it is not supported by any figure, table, or detailed calculation. This 5% bound is asserted without a time-resolved diagnostic, and no error bound connects a 5% deviation in F to the observed agreement in growth rates, plasma currents, or halo currents. Please provide a time-resolved diagnostic of F≡RBphi (or Bphi) at representative plasma locations (e.g., magnetic axis, outboard midplane, and edge) covering the entire nonlinear evolution, including the late phase when the plasma contacts the wall and halo currents are largest. In addition, explain why a 5% deviation is sufficient for the Slip Motion Condition to justify reduced MHD, and discuss whether this condition degrades at any point during the simulation.
- [Section 6] The nonlinear time traces are shifted so that the first wall-contact times coincide, but the raw times differ greatly: JOREK makes first wall contact at about 126 ms, NIMROD at about 87.4 ms, and M3D-C1 at about 91.5 ms. The paper attributes this to 'exponential dependency on the initial conditions' but does not state what initial perturbation or numerical noise was used in each code. Please report the initial perturbation amplitudes and show the unshifted traces (in the main text or supplementary material). Without this information, readers cannot assess whether the absolute timing of the VDE is part of the benchmark or whether the 'excellent agreement' in Fig. 8 is largely an artifact of the time alignment.
- [Section 6] The artificial thermal quench is imposed by multiplying the perpendicular heat diffusion coefficient by 500 and the particle diffusion coefficient by 20 when the plasma becomes limited by the wall. These are ad-hoc multipliers, yet no sensitivity study is provided. It would strengthen the benchmark to show, for at least one code, how the results (e.g., halo current or wall current) depend on the choice of these factors, or to state explicitly that the benchmark validates the codes' response to a prescribed thermal quench rather than predicting the quench itself. As it stands, the robustness of the excellent agreement to these choices is unclear.
minor comments (7)
- [Abstract and Section 5] The abstract uses 'excellent agreement' for the nonlinear simulations, but the linear growth rates show deviations up to 15% (Section 5). Consider using 'good' or 'favorable' consistently and quantifying the deviations in the abstract.
- [Figure 8(d)] The panel showing the toroidal current inside the LCFS lacks a NIMROD trace, with no explanation in the text. Please state why NIMROD results are omitted for this quantity or add the trace.
- [Table 1 and Section 5] The definition of Te,eff = Te,edge - Te,off is given in the table caption, but the text should clarify that this offset is used only in the Spitzer resistivity evaluation, not in the temperature evolution equation. This is implied but never explicitly stated.
- [Section 2] The estimate that the outer ideal wall influences growth rates by less than 10% is stated without derivation or reference. Please provide a brief justification or cite a previous study for this estimate.
- [Section 6] For JOREK, the halo current is computed from the equilibrium relation j×B = ∇p. Since the plasma may not be in exact force balance during a VDE, please discuss the potential error introduced by this approximation.
- [Section 7, Eq. (1)] The notation v = ξγ is nonstandard; it would be clearer to write v = γξ or v = dξ/dt in the linear stability context.
- [References] Reference [16] is an arXiv preprint; if a published version exists, please update it.
Circularity Check
No circularity found: this is an independent inter-code benchmark with no fitted target quantity reused as an input.
full rationale
This paper is a code-versus-code benchmark, not a derivation. The central comparisons are emergent outputs of three independently implemented MHD codes, and no predicted quantity (growth rate, plasma current, wall current, halo current) is fed back into any of the codes as an input. The growth rates are obtained by exponential fits to the magnetic-axis displacement time trace, but this is a measurement of a simulation output, not a parameter that is then used to produce the compared quantities. The JOREK reduced-MHD agreement is justified by the cited 'Slip Motion Condition' and a 5% deviation of RB_phi from its vacuum value; this is a physics assumption that could be wrong, but it is not circular because the criterion is not defined in terms of the benchmark outputs and the agreement itself is not assumed. Self-citations to the M3D-C1, NIMROD, and JOREK model papers are normal references to the codes under test, and they do not carry a load-bearing circular argument. The admitted need for 'further investigation' of the small linear/nonlinear mismatch is a stated limitation, not a circular step. The benchmark is self-contained as a verification exercise, so the appropriate circularity score is 0.
Assumptions & free parameters
free parameters (2)
- Effective edge temperature offset Te,off =
13.65 eV for Section 5 growth-rate scans; 0 eV for Section 6 full VDE run
- Artificial thermal quench multipliers =
kappa_perp x 500, Dn x 20 at first wall contact
assumptions (5)
- domain assumption Reduced MHD with the Slip Motion Condition represents full-MHD n=0 VDE dynamics in this configuration.
- domain assumption Spitzer resistivity with fixed effective charge and Coulomb logarithm models the plasma response.
- domain assumption Thin-wall and thick-wall resistive wall representations give equivalent VDE response for this geometry.
- ad hoc to paper The artificial thermal quench and time-shifting preserve the comparability of nonlinear VDE simulations.
- domain assumption The NSTX geqdsk equilibrium is an adequate common initial condition for all three codes.
Cite this review
Pith. "Pith review of Axisymmetric simulations of vertical displacement events in tokamaks: A benchmark of M3D-C$^1$, NIMROD and JOREK." pith.science (2026). https://pith.science/paper/HGGN6PLX
@misc{pith2026190802387,
author = {Pith},
title = {Pith review of: Axisymmetric simulations of vertical displacement events in tokamaks: A benchmark of M3D-C$^1$, NIMROD and JOREK},
year = {2026},
howpublished = {\url{https://pith.science/paper/HGGN6PLX}},
note = {Machine review of arXiv:1908.02387}
}
abstract
A benchmark exercise for the modeling of vertical displacement events (VDEs) is presented and applied to the 3D nonlinear magneto-hydrodynamic codes M3D-C$^1$, JOREK and NIMROD. The simulations are based on a vertically unstable NSTX equilibrium enclosed by an axisymmetric resistive wall with rectangular cross section. A linear dependence of the linear VDE growth rates on the resistivity of the wall is recovered for sufficiently large wall conductivity and small temperatures in the open field line region. The benchmark results show good agreement between the VDE growth rates obtained from linear NIMROD and M3D-C$^1$ simulations as well as from the linear phase of axisymmetric nonlinear JOREK, NIMROD and M3D-C$^1$ simulations. Axisymmetric nonlinear simulations of a full VDE performed with the three codes are compared and excellent agreement is found regarding plasma location and plasma currents as well as eddy and halo currents in the wall.
Figures
Figures from the paper (5 more)
Reference graph
Works this paper leans on
-
[1]
Dynamic response of the ITER tokamak during asymmetric VDEs,
T. Schioler, C. Bachmann, G. Mazzone, and G. Sannazzaro, “Dynamic response of the ITER tokamak during asymmetric VDEs,” Fusion Engineering and Design , vol. 86, no. 9-11, pp. 1963– 1966, 2011
work page 1963
-
[2]
Asymmetric wall force and toroidal rotation in tokamak disruptions,
H. Strauss, “Asymmetric wall force and toroidal rotation in tokamak disruptions,” Physics of Plasmas, vol. 22, no. 8, p. 082509, 2015
work page 2015
-
[3]
Reduction of asymmetric wall force in ITER disruptions with fast current quench,
H. Strauss, “Reduction of asymmetric wall force in ITER disruptions with fast current quench,” Physics of Plasmas , vol. 25, no. 2, p. 020702, 2018
work page 2018
-
[4]
Modelling of NSTX hot vertical displacement events using M3D-C 1,
D. Pfefferl´ e, N. Ferraro, S. C. Jardin, I. Krebs, and A. Bhattacharjee, “Modelling of NSTX hot vertical displacement events using M3D-C 1,” Physics of Plasmas , vol. 25, p. 056106, 2018
work page 2018
-
[5]
F. J. Artola Such, Free-boundary simulations of MHD plasma instabilities in tokamaks . PhD thesis, Universit´ e Aix-Marseille, 2018 (https://hal-amu.archives-ouvertes.fr/tel-02012234v1)
work page 2018
-
[6]
Effects of asymmetries in computations of forced vertical displacement events,
C. R. Sovinec and K. Bunkers, “Effects of asymmetries in computations of forced vertical displacement events,” Plasma Physics and Controlled Fusion , vol. 61, p. 024003, 2019
work page 2019
-
[7]
Nonlinear magnetohydrodynamics simulation using high-order finite elements,
C. R. Sovinec, A. H. Glasser, T. A. Gianakon, D. C. Barnes, R. A. Nebel, S. E. Kruger, D. D. Schnack, S. J. Plimpton, A. Tarditi, M. S. Chu, and the NIMROD Team, “Nonlinear magnetohydrodynamics simulation using high-order finite elements,” Journal of Computational Physics, vol. 195, pp. 355–386, March 2004
work page 2004
-
[8]
MHD stability in x-point geometry: simulation of ELMs,
G. Huysmans and O. Czarny, “MHD stability in x-point geometry: simulation of ELMs,” Nuclear fusion, vol. 47, no. 7, p. 659, 2007
work page 2007
Show all 19 references
-
[9]
Coupling JOREK and STARWALL codes for non-linear resistive- wall simulations,
M. Hoelzl, P. Merkel, G. Huysmans, E. Nardon, E. Strumberger, R. McAdams, I. Chapman, S. G¨ unter, and K. Lackner, “Coupling JOREK and STARWALL codes for non-linear resistive- wall simulations,” in Journal of Physics: Conference Series , vol. 401, p. 012010, IOP Publishing, 2012
2012
-
[10]
Multiple timescale calculations of sawteeth and other global macroscopic dynamics of tokamak plasmas,
S. C. Jardin, N. Ferraro, J. Breslau, and J. Chen, “Multiple timescale calculations of sawteeth and other global macroscopic dynamics of tokamak plasmas,” Computational Science & Discovery , vol. 5, p. 014002, May 2012
2012
-
[11]
Multi-region approach to free-boundary three-dimensional tokamak equilibria and resistive wall instabilities,
N. M. Ferraro, S. C. Jardin, L. L. Lao, M. S. Shephard, and F. Zhang, “Multi-region approach to free-boundary three-dimensional tokamak equilibria and resistive wall instabilities,” Physics of Plasmas, vol. 23, p. 056114, May 2016
2016
-
[12]
Jardin, Computational methods in plasma physics
S. Jardin, Computational methods in plasma physics . Chapman & Hall/CRC Computational Science, CRC Press, Taylor & Francis Group, 2010. Krebs et al. – VDE benchmark of nonlinear MHD codes 14
2010
-
[13]
3D two-temperature magnetohydrodynamic modeling of fast thermal quenches due to injected impurities in tokamaks,
N. Ferraro, B. C. Lyons, C. C. Kim, Y.-Q. Liu, and S. C. Jardin, “3D two-temperature magnetohydrodynamic modeling of fast thermal quenches due to injected impurities in tokamaks,” Nuclear Fusion, vol. 59, no. 1, p. 016001, 2018
2018
-
[14]
Studies of plasma equilibrium and transport in a tokamak fusion device with the inverse-variable technique,
R. Khayrutdinov and V. Lukash, “Studies of plasma equilibrium and transport in a tokamak fusion device with the inverse-variable technique,” Journal of Computational Physics , vol. 109, no. 2, pp. 193–201, 1993
1993
-
[15]
Dynamic modeling of transport and positional control of tokamaks,
S. C. Jardin, N. Pomphrey, and J. Delucia, “Dynamic modeling of transport and positional control of tokamaks,” Journal of computational Physics , vol. 66, no. 2, pp. 481–507, 1986
1986
-
[16]
Linear MHD stability studies with the STARWALL code,
P. Merkel and E. Strumberger, “Linear MHD stability studies with the STARWALL code,” arXiv:1508.04911, 2015
2015 arXiv
-
[17]
Non-linear MHD simulations of edge localized modes (ELMs),
G. Huysmans, S. Pamela, E. Van Der Plas, and P. Ramet, “Non-linear MHD simulations of edge localized modes (ELMs),” Plasma Physics and Controlled Fusion , vol. 51, no. 12, p. 124012, 2009
2009
-
[18]
Wigger, Development and Application of a Nonlinear Axisymmetric Resistive MHD-Code
C. Wigger, Development and Application of a Nonlinear Axisymmetric Resistive MHD-Code . PhD thesis, Technische Universit¨ at M¨ unchen, 2011
2011
-
[19]
Stability of tokamaks with respect to slip motions,
E. Rebhan and A. Salat, “Stability of tokamaks with respect to slip motions,” Nuclear Fusion, vol. 16, no. 5, p. 805, 1976
1976
Reviewed August 14, 2026 · model on record in the stance chip above.
Discussion (0). Continue with ORCID to comment.