|
MPI-AMRVAC 3.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
|
mod_magnetic_helicity is the MPI post-processing backend for relative magnetic helicity in a finite Cartesian volume. It is intentionally separate from mod_magnetic_topology: the latter traces field lines with OpenMP and must run on one MPI rank, whereas this module distributes planes with MPI and does not call the field-line tracer.
For a magnetic field (\mathbf B) in a finite volume, the reference field (\mathbf B_p) is the unique current-free field with the same normal component on every physical face. With vector potentials satisfying the same tangential-boundary convention, the gauge-invariant relative helicity is
[ H_m=\int_V(\mathbf A+\mathbf A_p)\cdot(\mathbf B-\mathbf B_p)\,dV . ]
The implementation also reports the gauge-invariant current-carrying and volume-threading decomposition attributed to Berger (2003) and written in a finite-volume form by Valori, Démoulin, and Pariat (2012):
[ H_J=\int_V(\mathbf A-\mathbf A_p)\cdot(\mathbf B-\mathbf B_p)\,dV, \qquad H_{PJ}=2\int_V\mathbf A_p\cdot(\mathbf B-\mathbf B_p)\,dV, ]
so that Hm = HJ + HPJ is checked independently. The factor 2 is part of the code's HPJ definition, not an additional post-processing multiplier. HJ measures the current-carrying (non-potential) part of the field; HPJ is the volume-threading/mutual term between the reference field and the current-carrying field. Neither is the field-line twist Tw from the topology backend. The optional abs_HJ_over_abs_Hm is omitted (ratio_is_valid=false) when the total helicity is at the round-off scale.
The energy diagnostics are
[ E=\frac12\int_V|\mathbf B|^2dV,\qquad E_p=\frac12\int_V|\mathbf B_p|^2dV,\qquad E_{free}=E-E_p. ]
In a normalized CGS run, helicities are written in (\mathrm{Mx}^2) and energies in erg; with SI_unit=.true. they are written in (\mathrm{Wb}^2) and joule. Code-unit columns are always retained.
The reference service mod_magnetic_reference_fv solves
[ \nabla^2\phi=0,\qquad \partial_n\phi=\mathbf B\cdot\hat{\mathbf n} \quad\text{on all six faces},\qquad \mathbf B_p=\nabla\phi . ]
The boundary value is the total normal field. For cell-centred storage it is obtained by the adjacent interior/ghost interpolation; for staggered storage the face value ws is used. B0field is added in both paths. A pure-Neumann problem is compatible only when the signed sum of the six face fluxes is zero. The code reports the relative imbalance and stops when it exceeds mh_max_flux_imbalance; it never silently rebalances the magnetogram. The additive constant in (\phi) is removed by the multigrid mean constraint, which does not change (\nabla\phi). The first reference solve uses one FMG pass and subsequent residual reductions use standalone V-cycles. The reusable reference-service default remains 1.d-8; the helicity user entry point defaults to mh_mg_tolerance=1.d-4.
The vector potentials use the DeVore–GV construction. For a selected axis (s) (mh_gauge_axis=1, 2, or 3; the default is 3, the z direction),
[ A_s=A_{p,s}=0,\qquad \mathbf A=\mathbf b+\hat{\mathbf e}_s\times \int_s^{s_{max}}\mathbf B\,ds',\qquad \mathbf A_p=\mathbf b_p+\hat{\mathbf e}_s\times \int_s^{s_{max}}\mathbf B_p\,ds'. ]
The top-face transverse vector (\mathbf b=\mathbf b_p) is built from a one-dimensional integral of the top-face normal field. Along the selected axis, a second-order recurrence matches the centred finite-volume curl, so the interior checks curl_A_error and curl_Ap_error directly test the discrete reconstruction.
The analysis first remaps an AMR snapshot through level_io to one fixed uniform level. It then assembles one plane at a time. Plane owners are distributed round-robin across MPI ranks; neighbouring owners exchange only the accumulated integration state. The full three-dimensional magnetic field is not replicated on each rank. Work arrays scale as (O(N_\perp^2)) per rank in the scalar path, in addition to the distributed reference multigrid data. Debug VTI output is streamed separately and is off by default.
Use the independent conversion type:
An empty mh_output_file uses <base_filename>_helicity.csv. Set mh_write_debug_vti=.true. and optionally mh_debug_vti_file to stream Bp, A, Ap, h_m, h_J, and h_PJ for numerical inspection. The optional mh_write_mg_timing prints the reference-solve multigrid timer table without changing the CSV schema.
The CSV header contains snapshot and time fields, code and physical-unit Hm, HJ, HPJ, abs_HJ_over_abs_Hm, E, Ep, Efree, net-flux and boundary-normal errors, epsilon_div_B, curl_A_error, curl_Ap_error, decomposition_error, mg_residual, mg_cycles, and ratio_is_valid.
Run with MPI and keep OpenMP/BLAS at one thread when benchmarking:
The directly usable two-step TDm/RBSL example is `tests/demo4/MagneticHelicity_Cart/README.md`. The topology documentation remains at `magnetic_topology_qsl.md`; it documents the single-rank OpenMP ‘convert_type='magnetic_topology’` path only.
The relative-helicity definition and its boundary condition are anchored in Berger–Field (1984) and Finn–Antonsen (1985). The HJ/HPJ split is the Berger (2003) current-carrying/volume-threading decomposition, while Valori, Démoulin, and Pariat (2012) provide the finite-volume derivation and the implementation/benchmark comparison used here. Fu/Yu (2023) is not used as a foundational definition citation and is not needed in this module documentation.
Version 1 accepts only three-dimensional, uniformly spaced, non-stretched Cartesian analysis meshes. Native AMR, spherical coordinates, and online time-step diagnostics are outside this conversion path; AMR data must be remapped with level_io. Incompatible six-face net flux is a hard error. Debug VTI is a numerical audit product, not the primary science output.
1.9.8