2 * One initial condition, on every variant's grid.
3 *
4 * The host's seeded perturbation is one normal deviate per *grid point*
5 * (src/mgpu/noise.ts), so two sessions at different lmax seeded from the same
6 * integer do not start from the same field — they start from unrelated fields
7 * that merely share a random seed. Comparing them would compare two different
8 * problems, and every number the comparison produced would be meaningless.
9 *
10 * So the field is built once, band-limited at the *coarsest* variant's lmax,
11 * and evaluated on each variant's own grid:
12 *
13 * 1. white noise on the coarsest grid
14 * 2. analysed there -> coefficients up to lmax_min
15 * 3. zero-padded into each variant's coefficient layout
16 * 4. synthesized on that variant's grid
17 *
18 * Steps 3 and 4 are exact: the field is band-limited at lmax_min, and every
19 * variant's band contains that, so each one receives the *same function*
20 * sampled where it needs it. Running each model's own `init` on it then leaves
21 * every session holding the identical spectral state (zero-padded), which is
22 * what makes a pointwise comparison at later times mean something.
23 *
24 * The coarsest variant gets the projected field too, not the raw white noise
25 * it was analysed from — otherwise it alone would start somewhere slightly
26 * different from the others.
28 * A model whose `init` calls `randnfun3` (all of the shipped ones do) draws its
29 * perturbation from a Fourier series on the surface's bounding box instead, and
30 * that needs no projection: it is a function of space, evaluated wherever it is
31 * asked, so one coefficient table *is* one field on every variant's grid. It
32 * still has to be drawn once rather than per session — see `sharedModes`.
34import { lmIndex, nlmCalc } from '../sht/layout.ts';
35import { seededNoise } from '../mgpu/noise.ts';
36import type { ModelSession } from '../mgpu/session.ts';
38/**
39 * Re-index coefficients from a band limit into a wider one's layout, zero-
40 * filling the degrees the source does not have. Both layouts are SHTNS
41 * m-major with mmax = lmax, so nothing but the index mapping changes.
42 */
43export function prolongCoeffs(
44 q: Float32Array,
45 lmaxFrom: number,
46 lmaxTo: number,
47): Float32Array {
48 if (lmaxTo === lmaxFrom) return q;
49 if (lmaxTo < lmaxFrom) {
50 throw new Error(`prolongCoeffs: cannot widen ${lmaxFrom} into a smaller ${lmaxTo}`);
51 }
52 const out = new Float32Array(2 * nlmCalc(lmaxTo, lmaxTo));
53 for (let m = 0; m <= lmaxFrom; m++) {
54 for (let l = m; l <= lmaxFrom; l++) {
55 const from = 2 * lmIndex(lmaxFrom, l, m);
56 const to = 2 * lmIndex(lmaxTo, l, m);
57 out[to] = q[from];
58 out[to + 1] = q[from + 1];
59 }
60 }
61 return out;
62}
64/**
65 * The same band-limited perturbation, sampled on each session's grid. Order
66 * follows `sessions`. Nothing may be in flight on any session's transform
67 * plan — the one-off analys/synth here use the plan's own scratch buffers.
68 */
69export async function sharedNoise(
70 sessions: ModelSession[],
71 amp: number,
72 seed: number,
73): Promise<Float32Array[]> {
74 let base = sessions[0];
75 for (const s of sessions) if (s.cfg.lmax < base.cfg.lmax) base = s;
76 const coeffs = await base.sht.analys(seededNoise(base.npts, amp, seed));
77 const out: Float32Array[] = [];
78 for (const s of sessions) {
79 out.push(await s.sht.synth(prolongCoeffs(coeffs, base.cfg.lmax, s.cfg.lmax)));
80 }
81 return out;
82}
beac00aMerge main into random-fieldsJeremy Magland 83
84/**
85 * The random field every variant seeds from, drawn once — from `reference`,
86 * whose numbers the study quotes — or null for a model that does not call
87 * `randnfun3`.
88 *
89 * One table for all of them is not merely an economy (the draw is interpreter
90 * time, and at a fine wavelength seconds of it). Each session would otherwise
91 * draw from *its own* bounding box, and a box comes from grid samples of the
92 * surface: at different lmax those differ in the last digits, and the draw is
93 * sensitive to the box — a different mode count consumes the RNG differently
94 * and the fields stop being the same one. Drawing once removes the question.
95 */
96export function sharedModes(
97 reference: ModelSession,
98 seed: number,
99): Promise<Float32Array | null> {
100 return reference.drawSeedModes(seed);
101}