3A reference solution for validating reaction-diffusion solvers that operate
4on curved surfaces, specifically the discretization of the surface
5Laplace-Beltrami operator. Any solver capable of solving reaction-diffusion
6systems on a closed surface can be checked against this case by starting
7from the same initial condition, running to the same physical end time, and
8comparing the resulting field in spherical-harmonic coefficient space.
10## Geometry
12A triaxial ellipsoid, parametrized by colatitude `θ ∈ [0, π]` and azimuth
13`φ ∈ [0, 2π)`:
15```
16x = ax · sin(θ) · cos(φ)
17y = ay · sin(θ) · sin(φ)
18z = az · cos(θ)
19```
21with semi-axes `ax`, `ay`, `az` (reference defaults: `1.5, 1.0, 0.6`).
23**The surface is band-limited.** The reference solution was not computed on
24the exact analytic ellipsoid above, but on its spherical-harmonic
25reconstruction truncated at a finite degree `lmax`: the analytic `(x, y, z)`
26was expanded in spherical harmonics and then re-synthesized from that
27truncated expansion, which very slightly rounds off the surface relative to
28the ideal ellipsoid. The truncated surface's own coefficients are included
29in the data file (see below). **A solver being validated against this case
30should reconstruct its working surface (and induced metric) from those
31coefficients, not from the analytic formula**, so that any disagreement in
32the final field reflects the reaction-diffusion solve and not a difference
33in which surface was actually used.
35## Equations
37The Allen–Cahn equation, with `Δ_g` the Laplace-Beltrami operator of the
38surface above. Unlike the Schnakenberg cases in this data set, Allen–Cahn is
39a single-species phase-field equation: an interface forms between `u ≈ +1`
40and `u ≈ −1` regions and then coarsens over time.
42```
43∂u/∂t = eps2·Δ_g(u) + u − u³
44```
46| parameter | meaning | reference value |
47|---|---|---|
48| `eps2` | interface-width parameter (`ε²`) | 0.01 |
50A solver is free to choose its own numerical method for `Δ_g` and its own
51timestep — those are not part of what is being tested. What matters is
52reaching the same physical end time `T`, computed as `steps × dt` from the
53values stored in the data file's metadata (see below), starting from the
54supplied initial condition.
56## Initial condition
58The initial `u` field is supplied as data — spherical-harmonic coefficients
59at `t = 0` — rather than as a formula to regenerate. A solver should
60synthesize its starting field directly from these coefficients. This
61sidesteps any need to reproduce whatever pseudorandom process originally
62generated the perturbation: the state at `t = 0` is simply data to load.
64## Comparison convention
66To compare against this reference, a solver's own final `u` field must be
67expressed in the same spherical-harmonic convention used throughout this
68data set:
70- **Orthonormal** spherical harmonics, **including the Condon-Shortley
71 phase**.
72- Because `u` is a real-valued field, only coefficients for `m ≥ 0` are
73 stored, as complex numbers. The `m < 0` coefficients are implied by
74 `Q(l, −m) = (−1)^m · conj(Q(l, m))`, and every `m = 0` coefficient has a
75 zero imaginary part.
76- **Ordering**: "m-major" — for `m = 0` up to `lmax`, and within each `m`,
77 for `l = m` up to `lmax`. The flat index of coefficient `(l, m)` is:
79 ```
80 index(l, m) = m·(lmax + 1) − m·(m − 1)/2 + (l − m)
81 ```
83- The total coefficient count is `nlm = (lmax + 1)(lmax + 2)/2`. Coefficients
84 are stored as a flat array of length `2·nlm`, with the real and imaginary
85 part of each coefficient adjacent: `[re(0), im(0), re(1), im(1), ...]`.
87A solver's comparison output should be a flat array of the same length,
88`2·nlm`, in the same ordering, for direct numerical comparison against the
89reference.
91## Data file format
93Each reference run is stored as a single HDF5 file. All coefficient arrays
94are 32-bit floats (the reference solution itself was computed in single
95precision, so no comparison should expect agreement finer than that).
96Scalar metadata is stored as HDF5 attributes, grouped by what they describe:
98```
99/ (attrs: a short description of the run, and the
100 name of the species stored under initial/ and
101 final/, e.g. "U")
102├─ spec/ run parameters as attributes: lmax, seed, steps,
103│ warmup, and the iteration count used internally
104│ by the reference solver's own Δ_g approximation
105│ ├─ params/ the model's own parameters: eps2, dt
106│ └─ geometry_params/ the geometry's own parameters: ax, ay, az
107├─ grid/ the reference solver's own working grid, as
108│ attributes: lmax, mmax, nlat, nphi, nlm
109├─ geometry/
110│ ├─ Gx dataset, float32[2·nlm] — the band-limited
111│ ├─ Gy surface's own spherical-harmonic coefficients,
112│ └─ Gz one per Cartesian component
113├─ initial/ one dataset (U), float32[2·nlm] — the initial
114│ condition at t = 0
115└─ final/ one dataset (U), float32[2·nlm] — the solution
116 after evolving to t = steps·dt
117```
119`lmax` (under `spec/` and `grid/`) is the spherical-harmonic truncation
120degree that defines `nlm` and the coefficient ordering above; it is the same
121value everywhere it appears in a given file. `seed` identifies which
122pseudorandom perturbation produced the initial condition, but is provenance
123information only — the initial condition itself is fully specified by the
124`initial/` coefficients, and does not need to be regenerated from the seed.
126Every dataset in `geometry/`, `initial/`, and `final/` has length `2·nlm`
127and follows the ordering convention described above.