10729f3WebGPU time loop for the cymbal2 shell model, with numbl setup and output stagesJeremy Magland 1// Driver: numbl setup stage -> WebGPU time loop -> numbl output stage.
2//
3// node src/run.mjs full default render, writes cymbal2-gpu.wav
4// node src/run.mjs --pinned the validate suite's reduced configuration
5// node src/run.mjs --compare GPU vs CPU (numbl cymbal2_sim) on the pinned config
6// node src/run.mjs --mode mega|multi|persistent backend (default multi)
7// node src/run.mjs --velocity 4 ... parameter overrides
9import { writeFileSync } from 'node:fs';
10import { runSetup, runOutput, runCpuSim, loadMFiles } from './setup.mjs';
11import { buildKernel, buildLayout } from './kernel.mjs';
12import { buildStages } from './kernels2.mjs';
13import { buildPersistent } from './kernels3.mjs';
14import { installNodeWebGpu, requestDevice, runMegakernel, runMultiDispatch, runPersistent } from './gpu.mjs';
15import { encodeWav } from './wav.mjs';
17function parseArgs(argv) {
18 const a = { params: {}, chunk: 2048, ws: 256, wav: 'cymbal2-gpu.wav', compare: false, bench: false, mode: 'multi' };
19 for (let i = 2; i < argv.length; i++) {
20 const t = argv[i];
21 if (t === '--pinned') Object.assign(a.params, { sample_rate: 44100, oversample: 1, tail_s: 1.5 });
22 else if (t === '--compare') { a.compare = true; Object.assign(a.params, { sample_rate: 44100, oversample: 1, tail_s: 1.5 }); }
23 else if (t === '--bench') a.bench = true;
24 else if (t === '--velocity') a.params.strike_velocity = Number(argv[++i]);
25 else if (t === '--count') a.params.strike_count = Number(argv[++i]);
26 else if (t === '--rate') a.params.strike_rate_hz = Number(argv[++i]);
27 else if (t === '--oversample') a.params.oversample = Number(argv[++i]);
28 else if (t === '--tail') a.params.tail_s = Number(argv[++i]);
29 else if (t === '--chunk') a.chunk = Number(argv[++i]);
30 else if (t === '--ws') a.ws = Number(argv[++i]);
31 else if (t === '--wav') a.wav = argv[++i];
32 else if (t === '--mode') a.mode = argv[++i];
33 else if (t === '--set') { const [k, v] = argv[++i].split('='); a.params[k] = Number(v); }
34 else throw new Error(`unknown argument ${t}`);
35 }
36 return a;
37}
39function hostDiagnostics(fhist, state, kernel, setup) {
40 let peakF = 0, contactSteps = 0, events = 0, prev = 0;
41 for (let i = 0; i < fhist.length; i++) {
42 const f = fhist[i];
43 if (f > peakF) peakF = f;
44 if (f > 0) { contactSteps++; if (prev === 0) events++; }
45 prev = f;
46 }
47 const L = kernel.layout;
48 return {
49 peakForceN: peakF,
50 contactMs: (1000 * contactSteps) / setup.simRate,
51 contactEvents: events,
52 peakXiMm: state[L.SC + 12] * 1000,
53 peakDispThickness: state[L.SC + 13] / setup.thickness,
54 contactVel: Array.from(state.subarray(L.ST_CVEL, L.ST_CVEL + setup.K)),
55 strokesPlayed: Math.min(state[L.SC + 4] + (state[L.SC + 3] > 0 ? 1 : 0), setup.K),
56 };
57}
59const rms = (x) => {
60 let s = 0;
61 for (let i = 0; i < x.length; i++) s += x[i] * x[i];
62 return Math.sqrt(s / x.length);
63};
65async function main() {
66 const args = parseArgs(process.argv);
68 console.log('setup (numbl): grid, per-wavenumber inverses, mallet, schedule ...');
69 let t0 = performance.now();
70 const setup = runSetup(args.params);
71 console.log(` ${(performance.now() - t0).toFixed(0)} ms; grid ${setup.nrt}x${setup.nth}, nfree ${setup.nfree}, nsteps ${setup.nsteps}, simRate ${setup.simRate}`);
73 await installNodeWebGpu();
74 const { device, info } = await requestDevice();
75 console.log(`device: ${info}`);
77 let kernel;
78 let gpu;
79 t0 = performance.now();
80 const progress = (k, n, ms) => {
81 if (k === n || k % (1024 * 64) === 0) {
82 process.stdout.write(`\r ${k}/${n} steps, ${(1000 * ms / k).toFixed(1)} us/step`);
83 }
84 };
85 if (args.mode === 'persistent') {
86 const L = buildLayout(setup);
87 kernel = L;
88 const pk = buildPersistent(L, { chunk: args.chunk >= 4096 ? 512 : Math.min(args.chunk, 2048) });
89 console.log(`persistent kernel: ${pk.nwg} workgroups of ${pk.wg}, ${pk.chunk} steps/dispatch, ops ${(L.opsData.byteLength / 1024).toFixed(0)} KB`);
90 gpu = await runPersistent(device, L, pk, progress);
91 } else if (args.mode === 'mega') {
92 kernel = buildKernel(setup, { ws: args.ws, chunk: args.chunk });
93 console.log(`megakernel: ws ${args.ws}, ${args.chunk} steps/dispatch, ops ${(kernel.opsData.byteLength / 1024).toFixed(0)} KB`);
94 gpu = await runMegakernel(device, kernel, setup.nsteps, progress);
95 } else {
96 const L = buildLayout(setup);
97 kernel = L; // shares .layout for diagnostics
98 const stages = buildStages(L, {});
99 console.log(`multi-dispatch: ${stages.order.length} dispatches/step, wg ${stages.wg} x ${stages.nwg}, ops ${(L.opsData.byteLength / 1024).toFixed(0)} KB`);
100 gpu = await runMultiDispatch(device, L, stages, progress);
101 }
102 console.log('');
103 const usPerStep = (1000 * gpu.gpuMs) / setup.nsteps;
104 console.log(`gpu integration: ${(gpu.gpuMs / 1000).toFixed(2)} s total, ${usPerStep.toFixed(1)} us/step`);
106 let bad = 0;
107 let firstBad = -1;
108 for (let i = 0; i < gpu.acc.length; i++) if (!Number.isFinite(gpu.acc[i])) { bad++; if (firstBad < 0) firstBad = i; }
109 if (process.env.DUMP_ACC) {
110 const { writeFileSync: wf } = await import('node:fs');
111 wf(process.env.DUMP_ACC, Buffer.from(gpu.acc.buffer, gpu.acc.byteOffset, gpu.acc.byteLength));
112 console.log(`acc dumped; nonfinite ${bad}${bad ? ' first at ' + firstBad : ''}`);
113 }
114 if (bad > 0 && !process.env.DUMP_ACC) throw new Error(`simulation diverged: ${bad} non-finite acceleration samples (first at ${firstBad})`);
115 if (bad > 0 && process.env.DUMP_ACC) { console.log('skipping output stage (diverged)'); return; }
117 const d = hostDiagnostics(gpu.fhist, gpu.state, kernel, setup);
118 console.log(`contact: peakF ${d.peakForceN.toFixed(2)} N, cvel ${d.contactVel[0].toFixed(3)} m/s, ` +
119 `${d.contactMs.toFixed(2)} ms over ${d.contactEvents} events, peakXi ${d.peakXiMm.toFixed(3)} mm, ` +
120 `peakDisp ${d.peakDispThickness.toFixed(1)} thicknesses`);
122 if (args.compare) {
123 console.log('cpu reference (numbl, unmodified cymbal2_sim) ...');
124 const cpu = runCpuSim(setup.mfiles, args.params);
125 console.log(` ${(cpu.ms / 1000).toFixed(1)} s total, ${(1000 * cpu.ms / cpu.nsteps).toFixed(1)} us/step`);
126 console.log(` cpu: peakF ${cpu.peakForceN.toFixed(2)} N, cvel ${cpu.contactVel.toFixed(3)}, contactMs ${cpu.contactMs.toFixed(2)}, events ${cpu.contactEvents}, peakXi ${cpu.peakXiMm.toFixed(3)} mm`);
128 // Contact window alignment and force agreement
129 const firstC = (f) => f.findIndex((x) => x > 0);
130 const kg = firstC(gpu.fhist), kc = firstC(cpu.fhist);
131 console.log(` first contact step: gpu ${kg}, cpu ${kc}${kg === kc ? '' : ' << MISALIGNED'}`);
132 let femax = 0;
133 for (let i = 0; i < cpu.fhist.length; i++) {
134 const e = Math.abs(gpu.fhist[i] - cpu.fhist[i]);
135 if (e > femax) femax = e;
136 }
137 console.log(` max |force diff| ${femax.toExponential(2)} N (peak ${cpu.peakForceN.toFixed(1)} N)`);
139 // Short-horizon acceleration agreement (before chaotic divergence)
140 for (const w of [256, 1024, 4096]) {
141 let num = 0, den = 0;
142 for (let i = kc; i < Math.min(kc + w, cpu.acc.length); i++) {
143 num += (gpu.acc[i] - cpu.acc[i]) ** 2;
144 den += cpu.acc[i] ** 2;
145 }
146 console.log(` acc rel-L2 over ${w} steps after contact: ${Math.sqrt(num / Math.max(den, 1e-30)).toExponential(2)}`);
147 }
148 console.log(` acc rms: gpu ${rms(gpu.acc).toExponential(4)}, cpu ${rms(cpu.acc).toExponential(4)}`);
149 }
151 console.log('output stage (numbl cymbal2_output) ...');
152 const out = runOutput(loadMFiles(), args.params, gpu.acc);
153 console.log(` y: ${out.y.length} samples, rms ${rms(out.y).toFixed(4)}, peakAccel ${out.peakAccel.toExponential(3)} m/s^2`);
154 writeFileSync(args.wav, encodeWav(out.y, args.params.sample_rate ?? 48000));
155 console.log(`wrote ${args.wav}`);
156}
158main().then(() => process.exit(0)).catch((e) => {
159 console.error(e);
160 process.exit(1);
161});