2 * Import a reference HDF5 file — geometry, initial and final spherical-
3 * harmonic coefficients for a run of this repo's solver, in the format
4 * documented alongside the sibling test-data repo's case files (see
5 * ../turing-surface-test-data/cases/) — run this repo's own solver from that
6 * file's exact initial condition, and report the numerical error against
7 * its final state.
8 *
9 * This is the regression check for the surface Laplace-Beltrami correction:
10 * replay a saved-off run and see how far this repo's own output has drifted
11 * (or use --niter to probe how much the correction term itself matters).
12 *
13 * npm run ref -- --in data/schnak-spots.h5
14 * npm run ref -- --in data/schnak-spots.h5 --niter 0
15 */
16import { requestShtDevice, describeAdapter } from '../src/sht/sht.ts';
17import { ModelSession } from '../src/mgpu/session.ts';
c90d0e2Check reference files in the browser's compare modeJeremy Magland 18import { extractReferenceCase, type H5Node } from '../src/compare/referenceCase.ts';
66ae13eUpdated command now tracks L_infty error too.Owen Melia 19import { relL2, relLinf } from '../src/mgpu/digest.ts';
90108a6Command for testing against a reference implementationOwen Melia 20import { installWebGpu, errMsg, NO_ADAPTER_HINT } from './nodeWebGpu.ts';
21import * as h5wasm from 'h5wasm/node';
23const USAGE = `usage: npm run ref -- --in <file> [options]
66ae13eUpdated command now tracks L_infty error too.Owen Melia 25 --in <file> the reference HDF5 file to check against (required)
26 --niter <n> override the solve iteration count (default: the file's own)
27 --tolerance <n> if given, exit 1 when any reported relL2 meets or exceeds it
28 --tolerance-linf <n> if given, exit 1 when any reported relLinf meets or exceeds it
29 --json machine-readable output
32Runs this repo's solver from the file's exact initial spectral state, to the
66ae13eUpdated command now tracks L_infty error too.Owen Melia 33same physical end time, and reports the relative-L2 and relative-L-infinity
34(max-norm) error of the resulting state against the file's final state (and,
35as a sanity check, of the regenerated geometry against the file's own
36geometry coefficients). --tolerance and --tolerance-linf gate independently:
37either can fail the run on its own.`;
39function fail(msg: string, code = 1): never {
40 console.error(`ref: ${msg}`);
41 process.exit(code);
42}
44const argv = process.argv.slice(2);
45if (argv.includes('--help') || argv.includes('-h')) {
46 console.log(USAGE);
47 process.exit(0);
48}
49let inFile: string | null = null;
50let niterOverride: number | null = null;
51let tolerance: number | null = null;
66ae13eUpdated command now tracks L_infty error too.Owen Melia 52let toleranceLinf: number | null = null;
90108a6Command for testing against a reference implementationOwen Melia 53const wantJson = argv.includes('--json');
54for (let i = 0; i < argv.length; i++) {
55 const a = argv[i];
56 if (a === '--json') continue;
57 const valued = (name: string): string | null => {
58 if (a === `--${name}`) return argv[++i];
59 if (a.startsWith(`--${name}=`)) return a.slice(name.length + 3);
60 return null;
61 };
62 const inv = valued('in');
63 if (inv !== null) {
64 inFile = inv;
65 continue;
66 }
67 const niterv = valued('niter');
68 if (niterv !== null) {
69 niterOverride = Number(niterv);
70 if (!Number.isInteger(niterOverride) || niterOverride < 0) {
71 fail(`--niter must be an integer >= 0 (got '${niterv}')`, 2);
72 }
73 continue;
74 }
66ae13eUpdated command now tracks L_infty error too.Owen Melia 75 const tolLinfv = valued('tolerance-linf');
76 if (tolLinfv !== null) {
77 toleranceLinf = Number(tolLinfv);
78 if (!Number.isFinite(toleranceLinf)) fail(`--tolerance-linf must be a number (got '${tolLinfv}')`, 2);
79 continue;
80 }
90108a6Command for testing against a reference implementationOwen Melia 81 const tolv = valued('tolerance');
82 if (tolv !== null) {
83 tolerance = Number(tolv);
84 if (!Number.isFinite(tolerance)) fail(`--tolerance must be a number (got '${tolv}')`, 2);
85 continue;
86 }
87 fail(`unrecognized argument '${a}'\n\n${USAGE}`, 2);
88}
89if (!inFile) fail(`--in <file> is required\n\n${USAGE}`, 2);
91let device: GPUDevice | null = null;
92let session: ModelSession | null = null;
93let h5file: InstanceType<typeof h5wasm.File> | null = null;
95try {
96 await h5wasm.ready;
97 h5file = new h5wasm.File(inFile, 'r');
c90d0e2Check reference files in the browser's compare modeJeremy Magland 98 const rc = extractReferenceCase(h5file as H5Node, inFile);
100 h5file = null;
c90d0e2Check reference files in the browser's compare modeJeremy Magland 102 const { model, geometry: geometryModel, params, geometryParams, lmax, steps } = rc;
103 const niter = niterOverride ?? rc.niter;
104 const fileGeom = rc.geometryCoeffs;
105 const fileInitial = rc.initial;
106 const fileFinal = rc.final;
90108a6Command for testing against a reference implementationOwen Melia 108 const runtime = await installWebGpu();
109 device = await requestShtDevice().catch((e: unknown) => {
110 throw new Error(`${errMsg(e)}\n${NO_ADAPTER_HINT}`);
111 });
112 const adapter = await describeAdapter(device);
114 session = await ModelSession.create({
115 device,
116 model,
117 params,
118 lmax,
119 geometry: geometryModel,
120 geometryParams,
121 niter,
122 });
66ae13eUpdated command now tracks L_infty error too.Owen Melia 124 const errorOf = (a: Float32Array, b: Float32Array) => ({ relL2: relL2(a, b), relLinf: relLinf(a, b) });
66ae13eUpdated command now tracks L_infty error too.Owen Melia 127 Gx: errorOf(session.geometry.X, fileGeom.X),
128 Gy: errorOf(session.geometry.Y, fileGeom.Y),
129 Gz: errorOf(session.geometry.Z, fileGeom.Z),
132 session.loadState(fileInitial);
133 session.step(steps);
66ae13eUpdated command now tracks L_infty error too.Owen Melia 135 const stateError: Record<string, { relL2: number; relLinf: number }> = {};
90108a6Command for testing against a reference implementationOwen Melia 136 for (const name of model.state) {
137 // Sequential: GpuModel.read() shares one staging buffer across calls.
138 const ours = await session.read(name);
66ae13eUpdated command now tracks L_infty error too.Owen Melia 139 stateError[name] = errorOf(ours, fileFinal[name]);
142 const allErrors = [...Object.values(geometryError), ...Object.values(stateError)];
66ae13eUpdated command now tracks L_infty error too.Owen Melia 143 const worstL2 = Math.max(...allErrors.map((e) => e.relL2));
144 const worstLinf = Math.max(...allErrors.map((e) => e.relLinf));
145 const passL2 = tolerance === null ? null : worstL2 < tolerance;
146 const passLinf = toleranceLinf === null ? null : worstLinf < toleranceLinf;
147 const checks = [passL2, passLinf].filter((p): p is boolean => p !== null);
148 const pass = checks.length === 0 ? null : checks.every(Boolean);
150 if (wantJson) {
151 console.log(
152 JSON.stringify(
153 {
154 in: inFile,
155 model: model.key,
156 geometry: geometryModel.key,
157 grid: { lmax, nlm: session.sht.nlm },
158 niter,
159 steps,
160 dt: params.dt,
161 T: steps * (params.dt ?? 0),
162 backend: { adapter, runtime, precision: 'fp32' },
163 geometryError,
164 stateError,
166 worstLinf,
169 passL2,
170 passLinf,
172 },
173 null,
174 2,
175 ),
176 );
177 } else {
178 console.log(`ref: ${inFile}`);
179 console.log(
180 ` model ${model.label} (${model.state.join(', ')})\n` +
181 ` geometry ${geometryModel.label} ` +
182 geometryModel.params.map((p) => `${p.key}=${geometryParams[p.key]}`).join(' ') +
183 `\n grid lmax ${lmax} · nlm ${session.sht.nlm}\n` +
c90d0e2Check reference files in the browser's compare modeJeremy Magland 184 ` niter ${niter}${niterOverride !== null ? ` (file: ${rc.niter})` : ''}\n` +
90108a6Command for testing against a reference implementationOwen Melia 185 ` run ${steps} steps, dt=${params.dt} (T=${(steps * (params.dt ?? 0)).toFixed(2)})\n`,
186 );
66ae13eUpdated command now tracks L_infty error too.Owen Melia 187 const fmtErr = (v: { relL2: number; relLinf: number }) =>
188 `relL2 ${v.relL2.toExponential(3)} relLinf ${v.relLinf.toExponential(3)}`;
189 console.log(` geometry check (regenerated vs file):`);
190 for (const [k, v] of Object.entries(geometryError)) console.log(` ${k} ${fmtErr(v)}`);
191 console.log(`\n final state (this run vs file):`);
192 for (const [k, v] of Object.entries(stateError)) console.log(` ${k} ${fmtErr(v)}`);
90108a6Command for testing against a reference implementationOwen Melia 193 if (tolerance !== null) {
194 console.log(
66ae13eUpdated command now tracks L_infty error too.Owen Melia 195 `\n worst relL2 ${worstL2.toExponential(3)} vs tolerance ${tolerance.toExponential(3)}: ` +
196 (passL2 ? 'PASS' : 'FAIL'),
197 );
198 }
199 if (toleranceLinf !== null) {
200 console.log(
201 ` worst relLinf ${worstLinf.toExponential(3)} vs tolerance-linf ${toleranceLinf.toExponential(3)}: ` +
202 (passLinf ? 'PASS' : 'FAIL'),
204 }
205 }
207 session.destroy();
208 device.destroy();
209 process.exit(pass === false ? 1 : 0);
210} catch (e) {
211 h5file?.close();
212 session?.destroy();
213 device?.destroy();
214 fail(errMsg(e));
215}