[spossann]

Summary

Profiling examples that show what polar splines can do: Poisson, CurlCurl and GyrokineticPoisson on a cylinder and on a torus. The manufactured solution has its largest gradient on the axis r = 0. Each case compares:

• C1 polar splines on IGAPolarCylinder / IGAPolarTorus (no hole);
• tensor-product splines on HollowCylinder / HollowTorus with a 1e-2 hole, the geometry of the existing profiling examples.

Stacked on #758 (polar StiffnessPreconditioner). The base is polar-stiffness, so this diff shows only the new commits.

Layout

• One params file per model and geometry: profiling/examples/{Poisson,CurlCurl}/polar_{cylinder,torus}_strong_scaling/params_*.py. Command-line flags select the variant:
  • --splines polar|tensor
  • --precond
  • --poloidal-solver schur|lu|ilu
  • --hole-bc free|dirichlet
  • --num-elements
• One submit script per model and geometry, with --variant: submit_{poisson,curl_curl}_polar_{cylinder,torus}_strong_scaling.py. The variant table sets the flags through ProfilingCase.launch(param_flags=...), and the script is run once per variant. Several ProfilingCases can't share one process: without --upload, setup_run deletes the shared output folder.
• Variants:
  • Poisson: tensor_{none,stiffness,multigrid}, polar_{none,tensor_mass,mass,mass_ilu,stiffness,stiffness_ilu}.
  • CurlCurl: tensor_{none,stiffness,hiptmair_xu,multigrid}, polar_{none,tensor_mass,mass,stiffness,stiffness_lu}.
  • There are no ILU variants for CurlCurl, because ILU is not robust there at σ = 1e-2 (see Polar splines: StiffnessPreconditioner (grad, curl) and stab_average for GyrokineticPoisson #758).
• Grid: (64, 66, 64), degree (3, 3, 3), for both spline spaces. Polar splines need a multiple of 3 elements in η2.
• Exact solution: a Cartesian Gaussian (w = 0.2) centred w/√2 off the axis, where its gradient is largest. The right-hand side is derived with sympy; for CurlCurl, E = g c. On top of the existing outputs, the error close to the axis (r < 0.1) and variant.json are saved.
• BC at the hole (tensor-product variants):
  • Poisson uses the natural BC.
  • CurlCurl uses Dirichlet. With the natural BC, the BC at the hole contradicts the solution, whose curl is largest on the axis. With σ = 1e-2 the error then grows to O(10³).

Fix: post-processing of polar fields

Output._load_femfields wrote the saved coefficients into a PolarVector. That either failed, or the field evaluated to zero, because SplineFunction.__call__ re-extracts from the polar vector. Now the saved tensor-product (extracted) coefficients are loaded into vector_stencil, and the polar vector is recovered with the restart extraction operator: once per snapshot, cached per space. New test test_pproc_polar (serial and parallel pproc) fails without the fix.

GyrokineticPoisson

• profiling/examples/GyrokineticPoisson/polar_{cylinder,torus}_strong_scaling/params_gyrokinetic_poisson.py and their submit scripts, with the same flags.
• Equilibria, as in the existing hollow GKP examples:
  • ScrewPinch (a=1, R0=1, Lz=2π), which_geometry="cylindrical";
  • AdhocTorus (a=1, R0=3, q_kind=0, p_kind=1), which_geometry="toroidal".
  • Defaults stab_mat="M0ad_withT" and diffusion_mat="M1gyro".
• Manufactured solution: the same off-axis Gaussian, times a modulation along η3, so the average ⟨φ⟩ is non-zero.
  • Cylinder: ⟨φ⟩ is exact.
  • Torus: ⟨φ⟩(r) = ∫g R dθ / ∫R dθ, tabulated by trapezoidal quadrature in θ.
  • ρ is computed with sympy in Cartesian coordinates. At startup, the symbolic n0, p and B0 are checked against the struphy equilibrium (n_xyz, p_xyz, b_xyz).
  • Leaving ⟨φ⟩ out of ρ raises the error from 6e-4 to 36%, so the averaged term is really tested.
• Variants: tensor_{none,stiffness,multigrid}, polar_{none,tensor_mass,mass,mass_ilu,stiffness,stiffness_ilu}.

Small grid (16 × 18 × 16, serial):

variant	cylinder: its	torus: its
polar, stiffness (Schur / ILU)	7 / 7	5 / 5
polar, none	55	76
polar, polar mass	419	311
polar, tensor mass	738	521
tensor-product, none	63	82
tensor-product, stiffness	82	72
tensor-product, multigrid	264	226

The error in φ depends only on the spline space: polar 4.9e-4 (cylinder) and 6.2e-4 (torus), tensor-product 1.5e-2 and 1.4e-2.

Torus at 32 × 33 × 32, 4 ranks, stiffness preconditioner:

splines	its	solve time	error in φ	peak memory (4 ranks)
polar	4	3.1 s	4.0e-5	7.0 GB
tensor-product	109	31 s	2.2e-2	5.6 GB

• Full-size runs (64 × 66 × 64) are for the cluster. They did not fit into 15 GB on a workstation; already at 32³ the tensor-product run takes 5.6 GB.
• The polar stiffness preconditioner is not exact for GKP (4–7 iterations rather than 2), because it drops the C^23 (η2–η3) coupling of M1gyro, see Polar splines: StiffnessPreconditioner (grad, curl) and stab_average for GyrokineticPoisson #758.
• The tensor-product preconditioners are worse than none on the hollow cylinder with the free hole BC: stiffness 82, multigrid 264, none 63.

Compacted layout of the existing cube/hollow examples

• The existing Poisson, CurlCurl and GyrokineticPoisson cube/hollow strong scaling cases now use the same layout as the polar ones:
  • one params file per model and geometry ({geo}_strong_scaling/params_*.py, flags --precond and --num-elements);
  • one submit script per model and geometry (--variant, default none).
• Removed: 19 per-variant params folders and 19 submit scripts (about −6300 lines).
• The labels of the former variants are unchanged (e.g. poisson_hollow_cylinder_strong_scaling_multigrid), so the result histories in the profiling-data repo continue. Names and descriptions are also byte-identical.
• The flags of a run are recorded in the packaged job script (slurm_script_str).
• Verified:
  • All 27 former variants (at 16³) give identical CG iterations and errors with the old and the new files.
  • A review compared each variant's effective parameters, names, descriptions and labels against the old files. The explicit multigrid=MultiGridOptions() was removed because it equals the propagator default.

Results

Small grid (16 × 18 × 16, serial; Poisson cylinder 16 × 18 × 8), Stiffness preconditioner:

case	polar: iterations, error	tensor-product (hole): iterations, error
Poisson cylinder	2, 8.8e-4	81 (no pc), 1.35e-2
Poisson torus	2, 5.4e-4	89 (no pc), 1.36e-2
CurlCurl cylinder	6, 1.4e-3	65, 1.24
CurlCurl torus	8, 1.3e-3	118, 24.6

Full size, Poisson cylinder, 64 × 66 × 64, 4 ranks, Stiffness preconditioner:

splines	iterations	solve time	max error in φ
polar	2	4.0 s	2.1e-6
tensor-product	192	158 s	2.9e-2

The tensor-product error does not converge under refinement because of the hole.

Notes for interpreting the runs:

• At full size, the unpreconditioned CurlCurl variants (and probably polar_tensor_mass) will hit maxiter = 10000. Their timings then measure 10000 iterations, not a solve.
• The mass-matrix preconditioners are worse than none for Poisson (they approximate M0, not the Laplacian). They are kept as reference variants.
• MultiGrid on the hollow torus with the free hole BC is worse than no preconditioner for Poisson. This is not a regression of this PR.

Tests

• test_pproc_polar (new): serial and 2 ranks.
• post_processing/tests/test_output.py and test_pproc.py (including --with-mpi on 2 ranks) pass.
• All six params files run in serial and on 2 ranks. All 17 GKP variants run at small size. The reviewer ran all 37 variants at small size: all exit 0.
• Two separate review passes (Poisson/CurlCurl, and GKP with an independent finite-difference check of the operator) found no bugs; their suggestions are applied.

🤖 Generated with Claude Code

[max-models]

In these submit scripts, there are way more parameters than we had in the old ones. It's not really a problem, but sometimes it may be hard to know exactly what was used when inspecting the results. Looks clear enough to my human eyes though.