80c0b6dPlate dynamics lab: plates in realizable mounts under gravity, explicit finite differences on WebGPUJeremy Magland 1import { describe, expect, it } from "vitest";
2import { NodeRunner } from "../src/node";
3import {
4 FD_METHOD,
5 LIBRARY_FILES,
6 buildGpuStages,
7 buildRunBundle,
8 parsePlateResult,
9 withEdits,
10 type ProblemSpec,
11 type RunResult,
12} from "../src/index";
13import { requestGpuDevice, type GpuHandle } from "../src/gpu/fdLoop";
14import { maxAbs, relL2, testSpec } from "./helpers";
16/**
17 * Desktop WebGPU for Node through the optional `webgpu` package (Google Dawn):
18 * `npm i --no-save webgpu`. Without the package, or on a machine without a GPU
19 * adapter, the suite is skipped rather than failed, so CI can run it
20 * opportunistically.
21 */
22async function dawn(): Promise<GpuHandle | null> {
23 try {
24 const specifier = "webgpu";
25 const mod = (await import(specifier)) as {
26 create: (flags: string[]) => GPU;
27 globals: Record<string, unknown>;
28 };
29 Object.assign(globalThis, mod.globals);
30 Object.defineProperty(globalThis, "navigator", {
31 value: { gpu: mod.create([]) },
32 configurable: true,
33 writable: true,
34 });
35 return await requestGpuDevice();
36 } catch (e) {
37 console.log(
38 `gpu.test: no WebGPU device under Node (${e instanceof Error ? e.message : e}); skipping`,
39 );
40 return null;
41 }
42}
44const handle = await dawn();
45const cpu = new NodeRunner();
46const gpu = new NodeRunner({ gpu: () => Promise.resolve(handle!) });
48async function both(spec: ProblemSpec, edits: Record<string, string> = {}) {
49 const stages = buildGpuStages(spec, FD_METHOD)!;
50 const g = await gpu.run({
51 bundle: withEdits(buildRunBundle(spec, FD_METHOD, { gpu: true }), edits),
52 gpu: stages,
53 });
54 const c = await cpu.run({ bundle: withEdits(buildRunBundle(spec, FD_METHOD), edits) });
55 if (!g.ok) throw new Error(`gpu run failed: ${g.error}\n${g.output.join("")}`);
56 if (!c.ok) throw new Error(`cpu run failed: ${c.error}\n${c.output.join("")}`);
57 return { g, c, rg: parsePlateResult(g.value), rc: parsePlateResult(c.value) };
58}
60const fin = (a: Float64Array) => a.filter(Number.isFinite);
62/** Relative L2 error over several rows at once, normalized by the norm of all reference rows together. */
63function relL2Rows(a: Float64Array[], b: Float64Array[]): number {
64 let num = 0;
65 let den = 0;
66 for (let k = 0; k < a.length; k++)
67 for (let i = 0; i < a[k].length; i++) {
68 num += (a[k][i] - b[k][i]) ** 2;
69 den += a[k][i] ** 2;
70 }
71 return a.length ? Math.sqrt(num / Math.max(den, 1e-300)) : 0;
72}
74describe.skipIf(!handle)("explicit finite differences on WebGPU", () => {
75 it.each([
76 "sheet-on-thread",
77 "gong-on-cords",
78 "sheet-on-foam",
79 "cymbal-on-stand",
80 "cantilever-in-vise",
81 "disk-in-hand",
82 ])("reproduces the double-precision loop for %s", async (scenario) => {
83 const spec = testSpec(scenario);
84 spec.fd.duration = 0.01;
85 const { g, rg, rc } = await both(spec);
86 expect(g.compute?.loop).toBe("webgpu");
87 expect(rg.info.extra.loop).toBe("webgpu");
88 expect(rg.info.stable).toBe(true);
89 expect(rg.t.length).toBe(rc.t.length);
90 expect(rg.frames.t.length).toBe(rc.frames.t.length);
91 const eW = relL2(rc.probe.w, rg.probe.w);
92 const eF = relL2(fin(rc.frames.w.data), fin(rg.frames.w.data));
93 const eVm = relL2(rc.probe.vm, rg.probe.vm);
94 const eA = relL2(rc.audio.p, rg.audio.p);
95 // a coordinate or a support that stays at zero has no scale of its own: compare against the largest
96 const eR = relL2Rows(rc.rigid.series, rg.rigid.series);
97 const eS = relL2Rows(rc.supports.series, rg.supports.series);
98 console.log(
99 `gpu.test ${scenario}: rel L2 probe w ${eW.toExponential(2)}, frames ${eF.toExponential(2)}, probe vm ${eVm.toExponential(2)}, audio ${eA.toExponential(2)}, rigid ${eR.toExponential(2)}, supports ${eS.toExponential(2)}, ${rg.info.nsteps} steps in ${g.compute?.gpuMs?.toFixed(0)} ms on the GPU`,
100 );
101 expect(eW).toBeLessThan(2e-3);
102 expect(eF).toBeLessThan(2e-3);
103 expect(eVm).toBeLessThan(2e-2);
104 expect(eA).toBeLessThan(2e-2);
105 expect(eR).toBeLessThan(2e-3);
106 expect(eS).toBeLessThan(2e-2);
107 expect(maxAbs(rg.probe.w)).toBeGreaterThan(0);
108 });
110 it("keeps the swing of a released plate for many steps", async () => {
111 const spec = testSpec();
112 spec.excitation = { type: "release", x: 0.15, y: 0.02, displacement: 0.002, radius: 0.015 };
113 spec.fd.duration = 0.05;
114 const { rg, rc } = await both(spec);
115 const eR = relL2(rc.rigid.series[1], rg.rigid.series[1]);
116 const eW = relL2(rc.probe.w, rg.probe.w);
117 console.log(
118 `gpu.test release: rel L2 swing ${eR.toExponential(2)}, probe ${eW.toExponential(2)}`,
119 );
120 expect(eR).toBeLessThan(1e-3);
121 expect(eW).toBeLessThan(2e-3);
122 });
124 it("stops at the blow-up where the double-precision loop does", async () => {
125 const spec = testSpec();
126 spec.fd.cfl = 1.03;
127 const { rg, rc } = await both(spec);
128 expect(rg.info.stable).toBe(false);
129 expect(rc.info.stable).toBe(false);
130 expect(rg.info.diverged_step).toBeGreaterThan(0);
131 expect(Math.abs(rg.info.diverged_step - rc.info.diverged_step)).toBeLessThanOrEqual(
132 Math.max(10, 0.1 * rc.info.nsteps),
133 );
134 expect(rg.probe.w.every(Number.isFinite)).toBe(true);
135 expect(rg.frames.t.length).toBe(rc.frames.t.length);
136 });
138 it("falls back to numbl when the loop file is edited, and says so", async () => {
139 const spec = testSpec();
140 const loop = LIBRARY_FILES.find((f) => f.name === "plate_fd_loop.m")!.source;
141 const { g, rg, rc } = await both(spec, {
142 "plate_fd_loop.m": loop + "\n% edited in the code view\n",
143 });
144 expect(g.compute?.loop).toBe("numbl");
145 expect(g.compute?.note).toMatch(/edited/);
146 expect(rg.info.extra.loop).toBe("numbl");
147 expect(relL2(rc.probe.w, rg.probe.w)).toBeLessThan(1e-12);
148 });
149});
151describe("explicit finite differences without a GPU", () => {
152 it("runs the whole bundle in numbl when the GPU is asked for but absent", async () => {
153 const spec = testSpec();
154 const r: RunResult = await cpu.run({
155 bundle: buildRunBundle(spec, FD_METHOD, { gpu: true }),
156 gpu: buildGpuStages(spec, FD_METHOD)!,
157 });
158 expect(r.ok).toBe(true);
159 if (!r.ok) return;
160 expect(r.compute?.loop).toBe("numbl");
161 expect(r.compute?.note).toMatch(/not available/);
162 const res = parsePlateResult(r.value);
163 expect(res.info.extra.loop).toBe("numbl");
164 expect(res.info.stable).toBe(true);
165 });
167 it("writes the staged script and the plain script for the same problem", () => {
168 const spec = testSpec();
169 const plain = buildRunBundle(spec, FD_METHOD).files[0].source;
170 const staged = buildRunBundle(spec, FD_METHOD, { gpu: true }).files[0].source;
171 expect(plain).toMatch(/result = plate_fd\(params\);/);
172 expect(staged).toMatch(/S = plate_fd_setup\(params\);/);
173 expect(staged).toMatch(/Rec = plate_fd_loop\(S\);/);
174 expect(staged).toMatch(/result = plate_fd_output\(S, Rec\);/);
175 const stages = buildGpuStages(spec, FD_METHOD)!;
176 expect(stages.loopFile).toBe("plate_fd_loop.m");
177 expect(stages.setup).toMatch(/find\(S\.KB\)/);
178 expect(stages.output).toMatch(/plate_fd_output\(S, Rec\)/);
179 });
180});