Pith. sign in

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 →

arxiv 1908.02387 v2 pith:HGGN6PLX submitted 2019-08-06 physics.plasm-ph physics.comp-ph

classification physics.plasm-phphysics.comp-ph PACS 52.30.-q52.55.Fa52.65.Kj
keywords verticaldisplacementeventVDEbenchmarknonlinearMHDsimulationresistivewallhalocurrentreducedslipmotionconditiontokamak
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

This paper establishes a benchmark case for vertical displacement events (VDEs), the vertical loss of control of a tokamak plasma toward the vessel wall, and uses it to compare three nonlinear MHD codes: M3D-C1, NIMROD, and JOREK. It claims that the linear VDE growth rates recovered from all three codes agree within about 10% across a range of wall resistivities, and that axisymmetric nonlinear simulations of a full VDE show close agreement in plasma position, plasma current, and eddy and halo currents in the wall. The practical payoff is a reusable verification case for any code aimed at VDE design studies. The paper also explains why JOREK, running a reduced MHD model, matches the full-MHD codes: the slip motion condition holds when the total toroidal field stays within 5% of the vacuum toroidal field.

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.

Watch

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

Editorial extensions of the paper, not claims the author makes directly.

  • 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.
Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

3 major / 7 minor

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)
  1. [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.
  2. [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.
  3. [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)
  1. [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.
  2. [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.
  3. [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.
  4. [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.
  5. [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.
  6. [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.
  7. [References] Reference [16] is an arXiv preprint; if a published version exists, please update it.

Circularity Check

0 steps flagged · score 0.0 of 10

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 2 free parameters · 5 assumptions · 0 invented entities

No quantity is fitted to force agreement between the codes; all three codes solve independently implemented resistive MHD models. The listed axioms are the modeling assumptions that make the benchmark physically meaningful, especially the reduced-MHD simplification in JOREK and the artificial quench used for nonlinear evolution.

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
    Chosen to set the effective edge temperature to 1 eV and avoid numerical problems at low edge temperature. It is not fitted to the code agreement, but it defines the edge-resistivity regime used for the benchmark.
  • Artificial thermal quench multipliers = kappa_perp x 500, Dn x 20 at first wall contact
    Ad hoc model choices introduced in Section 6 to mimic the enhanced transport that would come from 3D instabilities in an axisymmetric simulation. They affect the nonlinear comparison and are not derived from physics.
assumptions (5)
  • domain assumption Reduced MHD with the Slip Motion Condition represents full-MHD n=0 VDE dynamics in this configuration.
    Section 7 invokes F=RBphi approximately equal to F_vacuum, within 5%, to justify JOREK's reduced MHD model. This is plausible but not formally proved.
  • domain assumption Spitzer resistivity with fixed effective charge and Coulomb logarithm models the plasma response.
    Section 3 specifies eta(Te) = 1.03e-4 Z ln Lambda Te^-3/2 with Z=1 and ln Lambda=17. Benchmark consistency relies on all codes using the same resistivity model.
  • domain assumption Thin-wall and thick-wall resistive wall representations give equivalent VDE response for this geometry.
    NIMROD uses a thin-wall/vacuum coupling while M3D-C1 and JOREK use thick-wall or Green's function wall treatments (Sections 4-5). Agreement suggests equivalence, but no separate proof is given.
  • ad hoc to paper The artificial thermal quench and time-shifting preserve the comparability of nonlinear VDE simulations.
    Section 6 introduces factor-500 and factor-20 enhancements to transport coefficients at wall contact, and traces are shifted to align first wall contact. These choices are not physics-derived.
  • domain assumption The NSTX geqdsk equilibrium is an adequate common initial condition for all three codes.
    Section 5 notes small inconsistencies in the geqdsk equilibrium that relax after a few nonlinear steps; the exact cause is left for future work.

how reviews work

0 comments
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 reproduced from arXiv: 1908.02387 by the authors.

Figure 1
Figure 1. Mesh used for VDE benchmark case with M3D￾C 1 (black). Shown are the ideal wall domain boundary (purple), the coils (orange), the thick resis￾tive wall (between green and blue line) and the separatrix (red). The mesh has approximately 35000 ele￾ments. The M3D-C1 code is a high-order finite el￾ement code that solves the nonlinear time￾dependent extended MHD equations. It uses a split-implicit time advance in or￾der t… view at source ↗
Figure 2
Figure 2. Electron temperature at the outboard midplane versus the distance from the separatrix in a set of 2D nonlinear VDE simulations with different values of κk/κ⊥. VDE growth rates (γ) during the early drift phase have been obtained via an exponential fit to the time traces of Zaxis (during the initial 0.1 m of the displacement). Simulations are based on a DIII-D like equilibrium and have been performed with M3D-C1 . dom… view at source ↗
Figure 3
Figure 3. a) Linear VDE growth rates obtained from M3D-C1 simulations for different values of the wall resistivity (ηw) and of the resistivity in the open field line region (ηedge). Contour plots show the toroidal current density eigenfunctions of a case with ηedge = 3.1 × 10−5 Ω m, ηw = 3.0 × 10−7 Ω m where response currents form in the wall (b) and a case with ηedge = 3.1 × 10−5 Ω m, ηw = 1.0 × 10−1 Ω m where response curre… view at source ↗
Figures from the paper (5 more)
Figure 4
Figure 4. Figure 4: Equilibrium poloidal magnetic flux of the VDE benchmark case (M3D￾C 1 ). Also shown are the separatrix (red line) and the resistive wall (green and blue lines). 3. Benchmark set up The equilibrium used for this benchmark case is loosely based on the NSTX discharge #139…
Figure 5
Figure 5. Figure 5: Comparison of VDE growth rates from linear M3D-C1 and NIMROD simulations. The growth rates deviate by between 4% and 13%. for these computations, bicubic and biquartic elements have been applied. As described in Ref. [6], the NIMROD computations presented here use a th…
Figure 6
Figure 6. Figure 6: Comparison of VDE growth rates from the linear phase of 2D nonlinear M3D-C1 , NIMROD and JOREK simulations. They deviate between 0.3% and 15%. Also shown are the results of linear M3D-C1 calculations. corners at R = 0.02 m are not rounded. The computations also used fi…
Figure 7
Figure 7. Figure 7: Contour plots show the poloidal magnetic flux in the 2D nonlinear M3D￾C 1 simulation at the point in time when the plasma first becomes limited by the wall (a) and close to the end of the VDE (b). The time traces of the thermal energy (c) in the M3D-C1 , NIMROD and JOR…
Figure 8
Figure 8. Figure 8: Comparison of time traces from a 2D nonlinear simulation performed with JOREK, NIMROD and M3D-C1 : a) vertical position of magnetic axis, b) radial position of magnetic axis, c) toroidal current inside the LCFS and the open field line region, d) toroidal current inside…

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

19 extracted references · 18 canonical work pages

  1. [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

  2. [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

  3. [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

  4. [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

  5. [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)

  6. [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

  7. [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

  8. [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

Show all 19 references
  1. [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

  2. [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

  3. [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

  4. [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

  5. [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

  6. [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

  7. [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

  8. [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

  9. [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

  10. [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

  11. [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

Pith tools

Reviewed August 14, 2026 · model on record in the stance chip above.