// WebGPU time loop: packs the set-up into storage buffers, compiles the three kernels and runs // the sub-steps of a chunk of audio samples as one command buffer. Works in a browser and under // Node with the `webgpu` (Dawn) package, which the tests use. import { HOLD_HEIGHT } from "../model/exciter.ts"; import { strainEnergies } from "../model/forces.ts"; import { DECIMATION_SPAN } from "../model/radiation.ts"; import type { Setup } from "../model/setup.ts"; import type { ChunkResult, Diagnostics, SimulatorBackend } from "../model/simulator.ts"; import { ELEM_STRIDE, HIT_HEADER, HITU_HEADER, NODE_STRIDE, SYNC_HEADER, WORKGROUP_NODES, shaderSource } from "./shaders.ts"; // evaluated lazily: under Node the WebGPU globals are installed after this module is imported const storageUsage = () => GPUBufferUsage.STORAGE | GPUBufferUsage.COPY_DST | GPUBufferUsage.COPY_SRC; export class GpuSimulator implements SimulatorBackend { readonly kind = "gpu" as const; readonly setup: Setup; readonly device: GPUDevice; samplesDone = 0; private step = 0; private readonly G: number; private readonly HF: number; private readonly HU: number; private readonly buffers: Record = {}; private readonly pipelines: Record = {}; private readonly bindGroups: Record = {}; private readonly chunkCapacity: number; readonly movieEvery: number; private readonly framesCapacity: number; private readonly stagingPressure: GPUBuffer; private readonly stagingState: GPUBuffer; private readonly stagingHits: GPUBuffer; private readonly stagingHitsU: GPUBuffer; private readonly stagingMovie: GPUBuffer | null; /** ``movieEvery`` > 0 captures w every that many audio samples; the chunk size is reduced so a chunk's frames stay under 64 MB. */ constructor(device: GPUDevice, setup: Setup, chunkCapacity = 8192, movieEvery = 0) { this.device = device; this.setup = setup; this.movieEvery = movieEvery; const Nn = setup.mesh.n; if (movieEvery > 0) chunkCapacity = Math.max(128, Math.min(chunkCapacity, Math.floor((64e6 / (4 * Nn)) * movieEvery))); this.chunkCapacity = chunkCapacity; this.framesCapacity = movieEvery > 0 ? Math.ceil(chunkCapacity / movieEvery) + 2 : 0; const { mesh, shell: sh, schedule: sc, radiation: rad, substeps: M, dt } = setup; const N = mesh.n, E = mesh.e, K = sc.K, P = sc.P, W = DECIMATION_SPAN, D = rad.ringLength; this.G = Math.ceil(N / WORKGROUP_NODES); this.HF = HIT_HEADER + 5 * P; this.HU = HITU_HEADER + 3 * P; const G = this.G, HF = this.HF, HU = this.HU; // ---- pack the set-up const elemF = new Float32Array(ELEM_STRIDE * E + W * M); for (let e = 0; e < E; e++) { const b = ELEM_STRIDE * e; for (let i = 0; i < 3; i++) { elemF[b + i] = mesh.Gx[3 * e + i]; elemF[b + 3 + i] = mesh.Gy[3 * e + i]; } elemF[b + 6] = mesh.area[e]; elemF[b + 7] = sh.z0x[e]; elemF[b + 8] = sh.z0y[e]; elemF[b + 9] = sh.Cm[3 * e]; elemF[b + 10] = sh.Cm[3 * e + 1]; elemF[b + 11] = sh.Cm[3 * e + 2]; elemF[b + 12] = 2 * setup.config.material.sigma1 * sh.rhoHA[e]; for (let i = 0; i < 81; i++) elemF[b + 13 + i] = sh.Ke[81 * e + i]; } elemF.set(rad.taps, ELEM_STRIDE * E); const elemU = new Uint32Array(3 * E); elemU.set(mesh.tri); const nodeF = new Float32Array(NODE_STRIDE * N); for (let n = 0; n < N; n++) { const b = NODE_STRIDE * n; for (let d = 0; d < 5; d++) { nodeF[b + d] = sh.Minv[5 * n + d]; nodeF[b + 5 + d] = sh.M[5 * n + d]; } for (let k = 0; k < 6; k++) for (let d = 0; d < 5; d++) nodeF[b + 10 + 5 * k + d] = sh.Rm[(k * N + n) * 5 + d]; } const adj = new Uint32Array(N + 1 + sh.adjEntries.length); adj.set(sh.adjOffsets); adj.set(sh.adjEntries, N + 1); const hitsF = new Float32Array(HF * K), hitsU = new Uint32Array(HU * K); for (let k = 0; k < K; k++) { const hf = HF * k, hu = HU * k; hitsF[hf] = sc.kPoint[k]; hitsF[hf + 1] = sc.alpha[k]; hitsF[hf + 2] = sc.mass[k]; hitsF[hf + 3] = sc.gap[k]; hitsF[hf + 4] = sc.velocity[k]; hitsF[hf + 5] = HOLD_HEIGHT; hitsF[hf + 6] = 0; hitsF[hf + 7] = 0; hitsF[hf + 8] = 0; for (let p = 0; p < P; p++) { for (let i = 0; i < 3; i++) { hitsF[hf + HIT_HEADER + 3 * p + i] = sc.weights[(k * P + p) * 3 + i]; hitsU[hu + HITU_HEADER + 3 * p + i] = sc.nodes[(k * P + p) * 3 + i]; } hitsF[hf + HIT_HEADER + 3 * P + p] = sc.offsets[k * P + p]; } hitsU[hu] = sc.releaseStep[k]; hitsU[hu + 1] = 0; } const radCoef = new Float32Array(6 * E); for (let e = 0; e < E; e++) { radCoef[6 * e] = rad.coefF[e]; radCoef[6 * e + 1] = rad.coefB[e]; radCoef[6 * e + 2] = rad.i0F[e]; radCoef[6 * e + 3] = rad.frF[e]; radCoef[6 * e + 4] = rad.i0B[e]; radCoef[6 * e + 5] = rad.frB[e]; } const mk = (name: string, data: Float32Array | Uint32Array, usage = storageUsage()) => { const buf = device.createBuffer({ label: name, size: Math.max(16, data.byteLength), usage }); device.queue.writeBuffer(buf, 0, data.buffer, data.byteOffset, data.byteLength); this.buffers[name] = buf; return buf; }; mk("elemF", elemF); mk("elemU", elemU); mk("nodeF", nodeF); mk("adj", adj); mk("radCoef", radCoef); mk("hitsF", hitsF); mk("hitsU", hitsU); mk("state", new Float32Array(12 * N)); mk("elemForce", new Float32Array(15 * E)); mk("sync", new Float32Array(SYNC_HEADER + 2 * G * 12)); mk("rad", new Float32Array(E * W + E * D)); mk("pressure", new Float32Array(2 * chunkCapacity)); const chunkBuf = device.createBuffer({ label: "chunk", size: 16, usage: GPUBufferUsage.UNIFORM | GPUBufferUsage.COPY_DST }); this.buffers.chunk = chunkBuf; const staging = (size: number) => device.createBuffer({ size, usage: GPUBufferUsage.MAP_READ | GPUBufferUsage.COPY_DST }); this.stagingPressure = staging(8 * chunkCapacity); this.stagingState = staging(12 * N * 4); this.stagingHits = staging(Math.max(16, HF * K * 4)); this.stagingHitsU = staging(Math.max(16, HU * K * 4)); if (this.framesCapacity > 0) { mk("movie", new Float32Array(this.framesCapacity * N)); this.stagingMovie = staging(this.framesCapacity * N * 4); } else this.stagingMovie = null; // ---- pipelines const module = device.createShaderModule({ code: shaderSource({ N, E, K, P, M, W, D, G, dt, sigma0: setup.config.material.sigma0, nl: sh.nl }) }); const bgl = (entries: [number, "read-only-storage" | "storage" | "uniform"][]) => device.createBindGroupLayout({ entries: entries.map(([binding, type]) => ({ binding, visibility: GPUShaderStage.COMPUTE, buffer: { type } })), }); const L0 = bgl([[0, "read-only-storage"], [1, "read-only-storage"], [2, "read-only-storage"], [3, "storage"], [4, "storage"], [5, "read-only-storage"], [6, "storage"], [7, "storage"]]); const L1 = bgl([[0, "storage"], [1, "read-only-storage"], [2, "read-only-storage"], [3, "read-only-storage"], [4, "storage"], [5, "storage"], [6, "storage"]]); const L2 = bgl([[0, "read-only-storage"], [1, "read-only-storage"], [2, "read-only-storage"], [3, "storage"], [4, "uniform"]]); const L3 = bgl([[0, "read-only-storage"], [1, "storage"], [2, "uniform"], [3, "read-only-storage"]]); const empty = device.createBindGroupLayout({ entries: [] }); const mkPipe = (entryPoint: string, layouts: GPUBindGroupLayout[]) => device.createComputePipeline({ label: entryPoint, layout: device.createPipelineLayout({ bindGroupLayouts: layouts }), compute: { module, entryPoint }, }); this.pipelines.element = mkPipe("elementPass", [L0]); this.pipelines.node = mkPipe("nodePass", [empty, L1]); this.pipelines.gather = mkPipe("gatherPass", [empty, empty, L2]); if (this.framesCapacity > 0) this.pipelines.capture = mkPipe("capturePass", [empty, empty, empty, L3]); const B = this.buffers; const bg = (layout: GPUBindGroupLayout, bufs: GPUBuffer[]) => device.createBindGroup({ layout, entries: bufs.map((buffer, binding) => ({ binding, resource: { buffer } })) }); this.bindGroups.element = bg(L0, [B.elemF, B.elemU, B.state, B.elemForce, B.hitsF, B.hitsU, B.sync, B.rad]); this.bindGroups.node = bg(L1, [B.state, B.elemForce, B.nodeF, B.adj, B.hitsF, B.hitsU, B.sync]); this.bindGroups.gather = bg(L2, [B.sync, B.rad, B.radCoef, B.pressure, B.chunk]); if (this.framesCapacity > 0) this.bindGroups.capture = bg(L3, [B.state, B.movie, B.chunk, B.sync]); this.bindGroups.empty = device.createBindGroup({ layout: empty, entries: [] }); } /** Audio sample completed by the gather pass at sub-step n, or -1 when none is due. */ private sampleDueAt(n: number): number { const M = this.setup.substeps, W = DECIMATION_SPAN; if (n < 1 || (n - 1) % M !== 0) return -1; const j0 = (n - 1) / M; if (j0 === 0) return -1; const j = j0 - W / 2; return j >= 0 ? j : -1; } get maxChunk() { return this.chunkCapacity; } async runChunk(nSamples: number): Promise { if (nSamples > this.chunkCapacity) throw new Error("chunk larger than the pressure buffer"); const { device, setup } = this, M = setup.substeps, E = setup.mesh.e, K = setup.schedule.K, N = setup.mesh.n; const n0 = this.step, n1 = n0 + nSamples * M; let base = -1, count = 0; for (let n = n0; n < n1; n++) { const j = this.sampleDueAt(n); if (j >= 0) { if (base < 0) base = j; count++; } } // frames captured in this chunk: the state at the start of sub-step n where n / M is a multiple of movieEvery const frameSamples: number[] = []; const every = this.movieEvery; if (every > 0) for (let n = n0; n < n1; n++) if (n % (M * every) === 0) frameSamples.push(n / M); const frameBase = frameSamples.length ? frameSamples[0] / every : 0; device.queue.writeBuffer(this.buffers.chunk, 0, new Int32Array([base < 0 ? 0 : base, Math.max(count, 1), frameBase, Math.max(every, 1)])); const enc = device.createCommandEncoder(); enc.clearBuffer(this.buffers.pressure); const pass = enc.beginComputePass(); const wgE = Math.ceil((E + K) / 64), wgN = this.G; for (let n = n0; n < n1; n++) { if (every > 0 && n % (M * every) === 0) { pass.setPipeline(this.pipelines.capture); pass.setBindGroup(0, this.bindGroups.empty); pass.setBindGroup(1, this.bindGroups.empty); pass.setBindGroup(2, this.bindGroups.empty); pass.setBindGroup(3, this.bindGroups.capture); pass.dispatchWorkgroups(wgN); } pass.setPipeline(this.pipelines.element); pass.setBindGroup(0, this.bindGroups.element); pass.dispatchWorkgroups(wgE); if (this.sampleDueAt(n) >= 0) { pass.setPipeline(this.pipelines.gather); pass.setBindGroup(0, this.bindGroups.empty); pass.setBindGroup(1, this.bindGroups.empty); pass.setBindGroup(2, this.bindGroups.gather); pass.dispatchWorkgroups(1); } pass.setPipeline(this.pipelines.node); pass.setBindGroup(0, this.bindGroups.empty); pass.setBindGroup(1, this.bindGroups.node); pass.dispatchWorkgroups(wgN); } pass.end(); enc.copyBufferToBuffer(this.buffers.pressure, 0, this.stagingPressure, 0, 8 * this.chunkCapacity); const nFrames = frameSamples.length; if (nFrames > 0) enc.copyBufferToBuffer(this.buffers.movie, 0, this.stagingMovie!, 0, nFrames * N * 4); device.queue.submit([enc.finish()]); this.step = n1; this.samplesDone += nSamples; await this.stagingPressure.mapAsync(GPUMapMode.READ); const all = new Float32Array(this.stagingPressure.getMappedRange().slice(0)); this.stagingPressure.unmap(); const frames: Float32Array[] = []; if (nFrames > 0) { await this.stagingMovie!.mapAsync(GPUMapMode.READ, 0, nFrames * N * 4); const data = new Float32Array(this.stagingMovie!.getMappedRange(0, nFrames * N * 4).slice(0)); this.stagingMovie!.unmap(); for (let f = 0; f < nFrames; f++) frames.push(data.slice(f * N, (f + 1) * N)); } const size = Math.max(count, 1); return { firstSample: base < 0 ? this.samplesDone : base, pressureFront: all.slice(0, count), pressureBack: all.slice(size, size + count), frames, frameSamples }; } private async read(buffer: GPUBuffer, staging: GPUBuffer, bytes: number): Promise { const enc = this.device.createCommandEncoder(); enc.copyBufferToBuffer(buffer, 0, staging, 0, bytes); this.device.queue.submit([enc.finish()]); await staging.mapAsync(GPUMapMode.READ); const out = staging.getMappedRange().slice(0); staging.unmap(); return out; } async diagnostics(): Promise { const { setup } = this, sh = setup.shell, N = setup.mesh.n, sc = setup.schedule, HF = this.HF, HU = this.HU; const st = new Float32Array(await this.read(this.buffers.state, this.stagingState, 12 * N * 4)); const hf = new Float32Array(await this.read(this.buffers.hitsF, this.stagingHits, Math.max(16, HF * sc.K * 4))); const hu = new Uint32Array(await this.read(this.buffers.hitsU, this.stagingHitsU, Math.max(16, HU * sc.K * 4))); let ke = 0, woh = 0; for (let i = 0; i < 5 * N; i++) ke += 0.5 * sh.M[i] * st[5 * N + i] * st[5 * N + i]; for (let n = 0; n < N; n++) woh = Math.max(woh, Math.abs(st[5 * n + 2] / sh.hNode[n])); const { membrane, bending } = strainEnergies(sh, st.subarray(0, 5 * N)); let ex = 0; const hits = []; for (let k = 0; k < sc.K; k++) { ex += 0.5 * hf[HF * k + 2] * hf[HF * k + 6] ** 2; hits.push({ peakForce: hf[HF * k + 8], contactTime: hu[HU * k + 1] * setup.dt, velocity: hf[HF * k + 6], height: hf[HF * k + 5] }); } return { time: this.step * setup.dt, kinetic: ke, membrane, bending, exciters: ex, maxWoverH: woh, hits }; } async displacement(): Promise { const N = this.setup.mesh.n; const st = new Float32Array(await this.read(this.buffers.state, this.stagingState, 12 * N * 4)); const w = new Float32Array(N); for (let n = 0; n < N; n++) w[n] = st[5 * n + 2]; return w; } /** The state buffer (q at offset 0, 5 floats per node), for a renderer that reads w directly. */ get stateBuffer(): GPUBuffer { return this.buffers.state; } dispose() { for (const b of Object.values(this.buffers)) b.destroy(); this.stagingPressure.destroy(); this.stagingState.destroy(); this.stagingHits.destroy(); this.stagingHitsU.destroy(); this.stagingMovie?.destroy(); } }