2 * Browser validation suite: fp32 WebGPU transforms vs the f64 reference.
3 * Results are written to #log, console, and window.__RESULTS__ (read by
4 * scripts/test-gpu.mjs).
5 */
6import { ShtPlan, requestShtDevice } from '../src/sht.ts';
7import { ShtReference, randomSpectrum } from '../src/reference.ts';
8import { lmIndex, type ShtConfig } from '../src/layout.ts';
9import type { FourierMode } from '../src/sht.ts';
11interface CaseResult {
12 name: string;
13 pass: boolean;
14 detail: string;
15}
16const results: CaseResult[] = [];
17const logEl = document.getElementById('log')!;
18function log(line: string) {
19 console.log(line);
20 logEl.textContent += '\n' + line;
21}
23function relL2(a: ArrayLike<number>, b: ArrayLike<number>): number {
24 let num = 0,
25 den = 0;
26 for (let i = 0; i < a.length; i++) {
27 num += (a[i] - b[i]) ** 2;
28 den += b[i] ** 2;
29 }
30 return Math.sqrt(num / (den || 1));
31}
33async function runCase(
34 device: GPUDevice,
35 name: string,
36 cfg: ShtConfig,
37 fourier: FourierMode,
38 tolSynth: number,
39 tolAnalys: number,
40 tolRound: number,
41) {
42 const t0 = performance.now();
43 const plan = await ShtPlan.create(device, cfg, { fourier });
44 const ref = new ShtReference(cfg);
45 const q0 = randomSpectrum(cfg, 42 + cfg.lmax);
47 // synthesis vs f64 reference
48 const spatGpu = await plan.synth(q0);
49 const spatRef = ref.synth(q0);
50 const eSynth = relL2(spatGpu, spatRef);
52 // analysis of the reference field vs f64 reference
53 const spatRef32 = Float32Array.from(spatRef);
54 const qGpu = await plan.analys(spatRef32);
55 const qRef = ref.analys(spatRef32);
56 const eAnalys = relL2(qGpu, qRef);
58 // GPU round trip: synth -> analys, compare to original spectrum
59 const qRound = await plan.analys(spatGpu);
60 const eRound = relL2(qRound, q0);
62 const dt = (performance.now() - t0).toFixed(0);
63 const pass = eSynth < tolSynth && eAnalys < tolAnalys && eRound < tolRound;
64 const detail =
65 `mode=${plan.fourierMode} synth=${eSynth.toExponential(2)}/${tolSynth} ` +
66 `analys=${eAnalys.toExponential(2)}/${tolAnalys} round=${eRound.toExponential(2)}/${tolRound} (${dt}ms)`;
67 results.push({ name, pass, detail });
68 log(`${pass ? 'PASS' : 'FAIL'} ${name} ${detail}`);
69 plan.destroy();
70}
72async function runSingleModeCase(device: GPUDevice, name: string, cfg: ShtConfig, l: number, m: number) {
73 // synthesize a single (l, m) mode and analyze it back: spectrum should
74 // come back as the unit vector, and the field should match Y_lm exactly.
75 const plan = await ShtPlan.create(device, cfg, {});
76 const q = new Float32Array(2 * plan.nlm);
77 q[2 * lmIndex(cfg.lmax, l, m)] = 1.0;
78 const spat = await plan.synth(q);
79 const qBack = await plan.analys(spat);
80 const e = relL2(qBack, q);
81 const pass = e < 2e-5;
82 results.push({ name, pass, detail: `round=${e.toExponential(2)}` });
83 log(`${pass ? 'PASS' : 'FAIL'} ${name} round=${e.toExponential(2)}`);
84 plan.destroy();
85}
87async function main() {
88 if (!navigator.gpu) {
89 (window as any).__RESULTS__ = { fatal: 'navigator.gpu undefined (WebGPU unavailable)' };
90 log('FATAL: WebGPU unavailable');
91 return;
92 }
93 const adapter = await navigator.gpu.requestAdapter();
94 const info = adapter ? `${adapter.info?.vendor ?? '?'} / ${adapter.info?.architecture ?? '?'}` : 'none';
95 log(`adapter: ${info}`);
96 const device = await requestShtDevice();
97 device.addEventListener('uncapturederror', (ev) => {
98 log(`UNCAPTURED GPU ERROR: ${(ev as GPUUncapturedErrorEvent).error.message}`);
99 });
101 try {
102 // small, FFT path
103 await runCase(device, 'lmax=15 fft', { lmax: 15, mmax: 15, nlat: 32, nphi: 32 }, 'fft', 2e-6, 2e-6, 2e-6);
104 // small, DFT path (cross-check of the Fourier stages)
105 await runCase(device, 'lmax=15 dft', { lmax: 15, mmax: 15, nlat: 32, nphi: 36 }, 'dft', 2e-6, 2e-6, 2e-6);
106 // moderate
107 await runCase(device, 'lmax=63 fft', { lmax: 63, mmax: 63, nlat: 64, nphi: 128 }, 'auto', 5e-6, 5e-6, 5e-6);
108 // the SHTNS fp32 comfort zone boundary (SHT_L_RESCALE_FLY_FLOAT = 128)
109 await runCase(device, 'lmax=127 fft', { lmax: 127, mmax: 127, nlat: 128, nphi: 256 }, 'auto', 1e-5, 1e-5, 1e-5);
110 // reduced mmax
111 await runCase(device, 'lmax=127 mmax=40', { lmax: 127, mmax: 40, nlat: 144, nphi: 128 }, 'auto', 1e-5, 1e-5, 1e-5);
112 // beyond the comfort zone: rescaling must kick in (sin^m underflows f32
113 // around m ~ 90 at mid-latitudes); accuracy degrades gracefully
114 await runCase(device, 'lmax=255', { lmax: 255, mmax: 255, nlat: 256, nphi: 512 }, 'auto', 5e-5, 5e-5, 5e-5);
115 await runCase(device, 'lmax=399', { lmax: 399, mmax: 399, nlat: 400, nphi: 1024 }, 'auto', 2e-4, 2e-4, 2e-4);
116 // single-mode checks incl. a high-m sectoral mode (pure rescale territory)
117 await runSingleModeCase(device, 'mode (l=3,m=2)', { lmax: 15, mmax: 15, nlat: 32, nphi: 32 }, 3, 2);
118 await runSingleModeCase(device, 'mode (l=200,m=200)', { lmax: 200, mmax: 200, nlat: 224, nphi: 512 }, 200, 200);
119 } catch (e) {
120 results.push({ name: 'exception', pass: false, detail: String(e) });
121 log(`EXCEPTION: ${e instanceof Error ? e.stack ?? e.message : e}`);
122 }
124 const failed = results.filter((r) => !r.pass);
125 log(failed.length === 0 ? 'ALL GPU TESTS PASSED' : `${failed.length} GPU TEST(S) FAILED`);
126 (window as any).__RESULTS__ = { results, ok: failed.length === 0 };
127}
129main();