pith. sign in

arxiv: 2604.20766 · v1 · submitted 2026-04-22 · 🧮 math.NA · cs.NA

A provably stable numerical method for the anisotropic diffusion equation in confined magnetic fields: Curvilinear coordinates and multi-block domains

Pith reviewed 2026-05-09 23:25 UTC · model grok-4.3

classification 🧮 math.NA cs.NA
keywords anisotropic diffusioncurvilinear coordinatessummation by partsmulti-block domainsnumerical stabilityfinite difference methodsmagnetic fields
0
0 comments X

The pith

A summation-by-parts penalty scheme for anisotropic diffusion stays stable under curvilinear coordinate mappings and multi-block decompositions.

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

The paper constructs a numerical method that solves the anisotropic diffusion equation on domains described by curvilinear coordinates and split into multiple structured blocks. Diffusion perpendicular to magnetic field lines is discretized with summation-by-parts operators and simultaneous approximation terms, while diffusion parallel to the lines uses a volume penalty. The authors derive discrete energy estimates showing that the semi-discrete solution remains bounded when the curvilinear operators and interface penalties obey the summation-by-parts property. This stability result extends an earlier Cartesian version of the method and is checked on a sequence of manufactured solutions plus a qualitative test drawn from a magnetic equilibrium code.

Core claim

The semi-discrete scheme based on curvilinear SBP-SAT operators, volume penalties for parallel diffusion, and interface penalties for block coupling produces a discrete energy that is non-increasing, thereby establishing stability for the anisotropic diffusion problem in confined magnetic fields.

What carries the argument

Curvilinear summation-by-parts simultaneous-approximation-term operators together with volume and interface penalties that preserve the summation-by-parts property.

Load-bearing premise

The chosen coordinate mappings and block interfaces must admit SBP operators and penalties whose discrete energy contribution is negative semi-definite.

What would settle it

A numerical run on a smooth curved multi-block mesh in which the computed discrete energy increases for a manufactured solution would disprove the stability claim.

Figures

Figures reproduced from arXiv: 2604.20766 by Dean Muir, Kenneth Duru, Matthew Hole, Stuart Hudson.

Figure 1
Figure 1. Figure 1: Diagram of a multi-block curvilinear domain with two sub-blocks. The physical sub-domains (left) are [PITH_FULL_IMAGE:figures/full_fig_p005_1.png] view at source ↗
Figure 2
Figure 2. Figure 2: Effect of the dilation parameter on the mesh geometry. The value for the dilation parameter is listed to [PITH_FULL_IMAGE:figures/full_fig_p009_2.png] view at source ↗
Figure 3
Figure 3. Figure 3: Convergence for the perpendicular only manufactured solution on the circular domains with dilation of [PITH_FULL_IMAGE:figures/full_fig_p010_3.png] view at source ↗
Figure 4
Figure 4. Figure 4: Self convergence test in the five block circle test with magnetic field given by ( [PITH_FULL_IMAGE:figures/full_fig_p011_4.png] view at source ↗
Figure 5
Figure 5. Figure 5: Solution along chords extending from −0.9 ≤ x ≤ −0.5, y = 0 for varying resolutions for the second order case. The legend indicates the grid resolution such that n = nx = ny. Top: κ∥ = 106 . Bottom: κ∥ = 109 . We plot the reference solution with a partial Poincaré section in the top of [PITH_FULL_IMAGE:figures/full_fig_p012_5.png] view at source ↗
Figure 6
Figure 6. Figure 6: Top: Reference solution in upper half plane with partial Poincaré section overlayed to illustrate the solution conforming to the island shape. Bottom: The point-wise relative error. This can be seen to be maximised at the barrier between the interior grid and the outer domains, likely because the grid is not aligned with the magnetic field. The separatrix can also be seen in the error field as two curves o… view at source ↗
Figure 7
Figure 7. Figure 7: SPEC Poincaré section and low resolution multiblock grid. Note that the grid resolution here is purely [PITH_FULL_IMAGE:figures/full_fig_p014_7.png] view at source ↗
Figure 8
Figure 8. Figure 8: Contour plot showing the solution contours between [PITH_FULL_IMAGE:figures/full_fig_p015_8.png] view at source ↗
read the original abstract

We present a robust and accurate numerical method for the anisotropic diffusion equation in curvilinear coordinates. This study extends the recent work [Muir et al., Computer Physics Communications, 2025] for solving the anisotropic diffusion equation in magnetic fields from Cartesian meshes to to curvilinear coordinates and complex geometries. The method uses summation by parts with simultaneous approximation terms for computing the diffusion perpendicular to field lines. The diffusion along field lines is computed using a penalty approach, similar to a simultaneous approximation term, but applied across the volume. To extend the method to complex geometry we use a multi-block approach with piecewise smooth structured meshes. That is, the domain is split into sub-grids, with locally adjacent boundaries coupled weakly using penalties. We prove the semi-discrete stability for the curvilinear implementation by deriving discrete energy estimates. The approach is verified though a number of numerical tests, which demonstrate the convergence properties of the method in multi-domain approach. Finally, we present a qualitative result generated in complex geometry and magnetic field, which is generated by the Stepped Pressure Equilibrium Code.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Referee Report

1 major / 3 minor

Summary. The paper presents a numerical method for the anisotropic diffusion equation using summation-by-parts (SBP) operators with simultaneous approximation terms (SAT) in curvilinear coordinates and multi-block domains. It extends prior Cartesian-mesh work by introducing volume penalties for parallel diffusion and interface penalties for multi-block coupling, derives discrete energy estimates to prove semi-discrete stability, demonstrates convergence on numerical tests, and shows a qualitative application in complex geometry generated by the Stepped Pressure Equilibrium Code.

Significance. If the stability proof holds, the work supplies a robust, provably stable discretization for anisotropic diffusion in confined magnetic fields on complex geometries. This is valuable for plasma-physics and fusion applications where Cartesian meshes are insufficient. The explicit derivation of discrete energy estimates, the multi-block extension, and the reported convergence tests constitute concrete strengths that support reliability beyond heuristic validation.

major comments (1)
  1. The semi-discrete stability proof relies on the curvilinear SBP operators and the volume/interface penalties producing a non-positive contribution to the discrete energy. The manuscript should explicitly verify (or cite the precise conditions on the coordinate mappings) that the summation-by-parts property and negative semi-definiteness continue to hold after the curvilinear transformation and multi-block decomposition; without this, the energy estimate does not close for arbitrary mappings.
minor comments (3)
  1. Abstract: repeated word 'to to' in the phrase 'from Cartesian meshes to to curvilinear coordinates'.
  2. Abstract: 'verified though a number of numerical tests' should read 'verified through'.
  3. The description of the parallel-diffusion penalty as 'similar to a simultaneous approximation term, but applied across the volume' would benefit from a brief equation or reference to the precise penalty form to avoid ambiguity with standard SAT terms.

Simulated Author's Rebuttal

1 responses · 0 unresolved

We thank the referee for the positive summary, recognition of the work's significance for fusion applications, and recommendation of minor revision. We address the single major comment below and have incorporated a clarification to strengthen the presentation of the stability result.

read point-by-point responses
  1. Referee: The semi-discrete stability proof relies on the curvilinear SBP operators and the volume/interface penalties producing a non-positive contribution to the discrete energy. The manuscript should explicitly verify (or cite the precise conditions on the coordinate mappings) that the summation-by-parts property and negative semi-definiteness continue to hold after the curvilinear transformation and multi-block decomposition; without this, the energy estimate does not close for arbitrary mappings.

    Authors: We agree that the closure of the energy estimate requires the underlying SBP property to be preserved under the coordinate transformation. Our derivation in Section 3 proceeds by substituting the curvilinear SBP operators (defined via the standard transformation of the reference-element operators) into the discrete inner products and showing that all volume and interface penalty contributions are non-positive. This step relies on the known fact that the curvilinear SBP property holds when the mapping is sufficiently smooth and the Jacobian is positive (conditions stated, for example, in the foundational works on curvilinear SBP operators). In the revised manuscript we have added an explicit paragraph immediately preceding the energy estimate that cites these mapping requirements and notes that the multi-block SAT penalties are constructed identically to the single-block case, so the cancellation at interfaces remains unchanged. The proof therefore holds for the class of piecewise-smooth mappings employed in the Stepped Pressure Equilibrium Code geometries we consider; we do not claim validity for completely arbitrary (non-smooth) mappings. revision: yes

Circularity Check

0 steps flagged

Minor self-citation to prior Cartesian work; core stability derivation is independent

full rationale

The paper claims semi-discrete stability via explicit derivation of discrete energy estimates for the curvilinear SBP-SAT scheme, relying on the summation-by-parts property of the operators and negative semi-definiteness of the penalty terms. This is a direct mathematical derivation from stated operator properties rather than a fit or self-referential definition. The citation to Muir et al. (2025) provides context for the Cartesian precursor but is not load-bearing for the new curvilinear/multi-block proof, which is presented as self-contained. No step equates the stability result to its inputs by construction.

Axiom & Free-Parameter Ledger

0 free parameters · 2 axioms · 0 invented entities

The method rests on the existence of summation-by-parts operators in curvilinear coordinates and the ability to decompose the domain into piecewise smooth structured blocks that can be coupled by penalties while preserving the energy estimate.

axioms (2)
  • domain assumption Summation-by-parts finite-difference operators exist and satisfy the required telescoping property for the chosen curvilinear coordinate mappings.
    Invoked to obtain the discrete energy estimate for perpendicular diffusion.
  • domain assumption The domain admits a decomposition into locally adjacent piecewise smooth structured blocks whose interfaces can be treated with penalty terms.
    Required for the multi-block extension and weak coupling.

pith-pipeline@v0.9.0 · 5500 in / 1394 out tokens · 31957 ms · 2026-05-09T23:25:50.491902+00:00 · methodology

discussion (0)

Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.

Reference graph

Works this paper leans on

21 extracted references · 17 canonical work pages

  1. [1]

    and Koike, K

    Mark H. Carpenter, David Gottlieb, and Saul Abarbanel. “Time-Stable Boundary Conditions for Finite- DifferenceSchemesSolvingHyperbolicSystems:MethodologyandApplicationtoHigh-OrderCompactSchemes”. en. In:Journal of Computational Physics111.2 (Apr. 1994), pp. 220–236.issn: 0021-9991.doi:10.1006/ jcph.1994.1057.url:https://www.sciencedirect.com/science/artic...

  2. [2]

    Chacon and G

    L. Chacon and G. Di Giannatale.An asymptotic-preserving semi-Lagrangian algorithm for the anisotropic heat transport equation with arbitrary magnetic fields. arXiv:2408.02829 [math]. Aug. 2024.doi:10.48550/ arXiv.2408.02829.url:http://arxiv.org/abs/2408.02829(visited on 05/26/2025)

  3. [3]

    An asymptotic-preserving semi-Lagrangian algorithm for the time-dependent anisotropic heat transport equation

    L. Chacón, D. del-Castillo-Negrete, and C.D. Hauck. “An asymptotic-preserving semi-Lagrangian algorithm for the time-dependent anisotropic heat transport equation”. en. In:Journal of Computational Physics272 (Sept. 2014), pp. 719–746.issn: 00219991.doi:10.1016/j.jcp.2014.04.049.url:https://linkinghub. elsevier.com/retrieve/pii/S0021999114003258(visited on...

  4. [4]

    An asymptotic-preserving method for highly anisotropic elliptic equations based on a Micro–Macro decomposition

    Pierre Degond et al. “An asymptotic-preserving method for highly anisotropic elliptic equations based on a Micro–Macro decomposition”. en. In:Journal of Computational Physics231.7 (Apr. 2012), pp. 2724–2740. issn: 0021-9991.doi:10.1016/j.jcp.2011.11.040.url:https://www.sciencedirect.com/science/ article/pii/S0021999111006966(visited on 07/07/2021)

  5. [5]

    Review of summation-by-parts operators with simultaneous approximation terms for the numerical solution of partial differential equations

    David C. Del Rey Fernández, Jason E. Hicken, and David W. Zingg. “Review of summation-by-parts operators with simultaneous approximation terms for the numerical solution of partial differential equations”. en. In: Computers & Fluids95 (May 2014), pp. 171–196.issn: 0045-7930.doi:10 . 1016 / j . compfluid . 2014 . 02.016.url:https://www.sciencedirect.com/sc...

  6. [6]

    Stable and high order accurate difference methods for the elastic wave equation in discontinuous media

    Kenneth Duru and Kristoffer Virta. “Stable and high order accurate difference methods for the elastic wave equation in discontinuous media”. en. In:Journal of Computational Physics279 (Dec. 2014), pp. 37–62.issn: 0021-9991.doi:10.1016/j.jcp.2014.08.046.url:https://www.sciencedirect.com/science/article/ pii/S0021999114006226(visited on 11/11/2021)

  7. [7]

    Journal of Computational Physics295, 456–474 (2015) https://doi.org/10.1016/j.jcp.2015

    Bram van Es, Barry Koren, and Hugo J. de Blank. “Finite-volume scheme for anisotropic diffusion”. In: Journal of Computational Physics306 (Feb. 2016), pp. 422–442.issn: 0021-9991.doi:10.1016/j.jcp.2015. 11.041.url:https://www.sciencedirect.com/science/article/pii/S0021999115007810(visited on 08/15/2024)

  8. [8]

    Lawrence Evans.Partial Differential Equations, Second edition. en. Vol. 19. Graduate Studies in Mathematics. ISSN:1065-7339. AmericanMathematicalSociety,Mar.2010.isbn: 978-0-8218-4974-3978-1-4704-1144-2.doi: 10.1090/gsm/019.url:http://www.ams.org/gsm/019(visited on 04/08/2021)

  9. [9]

    Helical temperature perturbations associated with tearing modes in tokamak plasmas

    Richard Fitzpatrick. “Helical temperature perturbations associated with tearing modes in tokamak plasmas”. In:Physics of Plasmas2.3 (Mar. 1995), pp. 825–838.issn: 1070-664X.doi:10.1063/1.871434.url:https: //doi.org/10.1063/1.871434(visited on 05/27/2023)

  10. [10]

    A Skew-Symmetric Discontinuous Galerkin Spectral Element Discretization and Its Re- lation to SBP-SAT Finite Difference Methods

    Gregor J. Gassner. “A Skew-Symmetric Discontinuous Galerkin Spectral Element Discretization and Its Re- lation to SBP-SAT Finite Difference Methods”. In:SIAM Journal on Scientific Computing35.3 (Jan. 2013). Publisher: Society for Industrial and Applied Mathematics, A1233–A1253.issn: 1064-8275.doi:10.1137/ 120890144.url:https://epubs.siam.org/doi/10.1137/1...

  11. [11]

    Solving Maxwell’s Equations Using the Ultra Weak Variational Formulation.Journal of Computational Physics, 223(2):731–758, 2007

    S. Günter et al. “Modelling of heat transport in magnetised plasmas using non-aligned coordinates”. en. In: Journal of Computational Physics209.1 (Oct. 2005), pp. 354–370.issn: 0021-9991.doi:10.1016/j.jcp. 2005.03.021.url:https://www.sciencedirect.com/science/article/pii/S0021999105001373(visited on 03/13/2021)

  12. [12]

    On heat conduction in an irregular magnetic field. Part 1

    Per Helander, Stuart R. Hudson, and Elizabeth J. Paul. “On heat conduction in an irregular magnetic field. Part 1”. en. In:Journal of Plasma Physics88.1 (Feb. 2022), p. 905880122.issn: 0022-3778, 1469-7807. doi:10 . 1017 / S002237782100129X.url:https : / / www . cambridge . org / core / product / identifier / S002237782100129X/type/journal_article(visited...

  13. [13]

    Temperature Contours and Ghost Surfaces for Chaotic Magnetic Fields

    S. R. Hudson and J. Breslau. “Temperature Contours and Ghost Surfaces for Chaotic Magnetic Fields”. en. In:Physical Review Letters100.9 (Mar. 2008), p. 095001.issn: 0031-9007, 1079-7114.doi:10 . 1103 / PhysRevLett.100.095001.url:https://link.aps.org/doi/10.1103/PhysRevLett.100.095001(visited on 08/19/2020)

  14. [14]

    Computation of multi-region relaxed magnetohydrodynamic equilibria

    S. R. Hudson et al. “Computation of multi-region relaxed magnetohydrodynamic equilibria”. In:Physics of Plasmas19.11 (Nov. 2012), p. 112502.issn: 1070-664X.doi:10 . 1063 / 1 . 4765691.url:http : / / aip . scitation.org/doi/10.1063/1.4765691(visited on 12/03/2017). [15]liuyxpp/CubicHermiteSpline.jl: Pure Julia implementation of 1D & 2D cubic Hermite spline...

  15. [15]

    Summation by Parts Operators for Finite Difference Approximations of Second-Derivatives with Variable Coefficients

    Ken Mattsson. “Summation by Parts Operators for Finite Difference Approximations of Second-Derivatives with Variable Coefficients”. en. In:Journal of Scientific Computing51.3 (June 2012), pp. 650–682.issn: 1573- 7691.doi:10.1007/s10915-011-9525-z.url:https://doi.org/10.1007/s10915-011-9525-z(visited on 10/01/2021)

  16. [16]

    A provably stable approach to solving the anisotropic diffusion equation in magnetic fields

    Dean Muir. “A provably stable approach to solving the anisotropic diffusion equation in magnetic fields”. In: (2024).url:https://hdl.handle.net/1885/733728434(visited on 02/16/2025)

  17. [17]

    A provably stable numerical method for the anisotropic diffusion equation in confined magnetic fields

    Dean Muir et al. “A provably stable numerical method for the anisotropic diffusion equation in confined magnetic fields”. In:Computer Physics Communications310 (May 2025), p. 109536.issn: 0010-4655.doi: 10 . 1016 / j . cpc . 2025 . 109536.url:https : / / www . sciencedirect . com / science / article / pii / S0010465525000396(visited on 02/13/2025)

  18. [18]

    arXiv:2303.15447 [cs, math]

    Dean Muir et al.An efficient method for the anisotropic diffusion equation in magnetic fields. arXiv:2303.15447 [cs, math]. Feb. 2023.doi:10.48550/arXiv.2303.15447.url:http://arxiv.org/abs/2303.15447(visited on 06/05/2023)

  19. [19]

    Stepped pressure equilibrium with relaxed flow and applications in reversed-field pinch plasmas

    Z. S. Qu et al. “Stepped pressure equilibrium with relaxed flow and applications in reversed-field pinch plasmas”. en. In:Plasma Physics and Controlled Fusion62.5 (Apr. 2020). tex.ids= quSteppedPressureEqui- librium2020 publisher: IOP Publishing, p. 054002.issn: 0741-3335.doi:10.1088/1361-6587/ab7fc5.url: http://iopscience.iop.org/article/10.1088/1361-658...

  20. [20]

    Code Verification by the Method of Manufactured Solutions

    Patrick J. Roache. “Code Verification by the Method of Manufactured Solutions”. In:Journal of Fluids Engi- neering124.1 (Nov. 2001), pp. 4–10.issn: 0098-2202.doi:10.1115/1.1436090.url:https://doi.org/10. 1115/1.1436090(visited on 10/16/2022)

  21. [21]

    On Numerical Solution of Strongly Anisotropic Diffusion Equation on Misaligned Grids

    M. V. Umansky, M. S. Day, and T. D. Rognlien. “On Numerical Solution of Strongly Anisotropic Diffusion Equation on Misaligned Grids”. In:Numerical Heat Transfer, Part B: Fundamentals47.6 (June 2005). Pub- lisher: Taylor & Francis _eprint: https://doi.org/10.1080/10407790590928946, pp. 533–554.issn: 1040-7790. doi:10 . 1080 / 10407790590928946.url:https : ...