/ proof-of-concept / plate-mechanics-sim
Sign in
proof-of-concept / plate-mechanics-sim
plate-mechanics-sim / packages / core / src / methods / implicit.ts
80 lines · 4.4 KBCodeBlameHistory
2546961Plate mechanics lab: Kirchhoff-Love plates solved by numbl in the browserJeremy Magland 1import type { MethodDefinition } from "./types";
78d40e9Free edges for both finite-difference methodsJeremy Magland 4 FREE_EDGE_EQUATIONS,
5 GOVERNING_EQUATIONS,
78d40e9Free edges for both finite-difference methodsJeremy Magland 7 RADIATION_EQUATIONS,
8 STRESS_EQUATIONS,
9} from "./common";
11export const IMPLICIT_METHOD: MethodDefinition = {
12 id: "implicit",
13 name: "Implicit finite differences",
14 shortName: "Implicit FD",
15 summary:
16 "Same stencils as the explicit scheme, stiffness averaged in time; one sparse factorization, no step limit.",
17 solverFunction: "plate_implicit",
18 planFunction: "plate_implicit_plan",
19 supports: {
81f3b30Elliptical plates, a WebGPU time loop for the explicit method, and WAV exportJeremy Magland 20 shape: ["rectangle", "ellipse"],
21 bc: ["simply-supported", "clamped", "free", "pinned-center"],
2546961Plate mechanics lab: Kirchhoff-Love plates solved by numbl in the browserJeremy Magland 22 excitation: ["impulse", "harmonic"],
23 audio: "clip",
24 },
25 readingOrder: [
26 "plate_implicit.m",
27 "plate_implicit_plan.m",
78d40e9Free edges for both finite-difference methodsJeremy Magland 30 "plate_rigid_modes.m",
31 "plate_curvature_free.m",
34 "plate_force_time.m",
35 "plate_force_footprint.m",
36 "plate_stress.m",
37 "plate_pad.m",
38 "plate_curvature.m",
39 "plate_listener.m",
40 "plate_spectrum.m",
41 "plate_result_schema.m",
42 ],
43 equations: [
44 ...GOVERNING_EQUATIONS,
45 {
46 tex: String.raw`\frac{w^{n+1} - 2w^n + w^{n-1}}{\Delta t^2} + \sigma_0\frac{w^{n+1} - w^{n-1}}{\Delta t} - \sigma_1 \mathsf L\frac{w^{n+1} - w^{n-1}}{\Delta t} + \kappa^2 \mathsf B\left(\theta w^{n+1} + (1-2\theta) w^n + \theta w^{n-1}\right) = \frac{F^n s}{\mu}`,
47 caption:
81f3b30Elliptical plates, a WebGPU time loop for the explicit method, and WAV exportJeremy Magland 48 "Theta scheme: the biharmonic term is averaged over three time levels. L and B are the sparse Laplacian and biharmonic matrices built from the same stencils as the explicit scheme, on the plate's footprint.",
50 {
51 tex: String.raw`\mathbf A\, w^{n+1} = \mathbf C\, w^n - \mathbf A_{\mathrm{old}}\, w^{n-1} + \frac{\Delta t^2}{\mu} F^n s, \qquad \mathbf A = (1+\sigma_0\Delta t)\mathbf I - \sigma_1\Delta t\,\mathsf L + \theta\kappa^2\Delta t^2\,\mathsf B`,
52 caption:
53 "One linear system per step with a constant symmetric positive definite matrix: factor A once by sparse Cholesky, then two triangular solves per step.",
54 },
55 {
56 tex: String.raw`\theta \ge \tfrac14 \;\Rightarrow\; \text{stable for every } \Delta t, \qquad \frac{T_{\text{num}}}{T} \approx 1 + \frac{(\omega\Delta t)^2}{12} \;\; (\theta = \tfrac14)`,
57 caption:
58 "No stability limit, so the step is chosen for accuracy. The price is period elongation: modes near half the sample rate run slow by tens of percent unless several steps per sample are used.",
59 },
78d40e9Free edges for both finite-difference methodsJeremy Magland 60 ...FREE_EDGE_EQUATIONS,
62 ...CURVED_EDGE_EQUATIONS,
64 ...RADIATION_EQUATIONS,
65 ],
66 validity: [
67 "Second-order truncation error in dx and dy, identical to the explicit scheme, since the spatial stencils are the same.",
68 "Second order in dt but dispersive: high modes are detuned rather than unstable. Two steps per sample keeps a 10 kHz mode within about 4 percent.",
81f3b30Elliptical plates, a WebGPU time loop for the explicit method, and WAV exportJeremy Magland 69 "The probe is snapped to the nearest grid node with a full stencil inside the plate; stress fields are evaluated only at the saved frames.",
2546961Plate mechanics lab: Kirchhoff-Love plates solved by numbl in the browserJeremy Magland 70 "Cost per step is a sparse triangular solve pair, about a millisecond for a 61 x 41 grid in the browser, independent of the plate stiffness.",
71 "Audio covers only the simulated window.",
78d40e9Free edges for both finite-difference methodsJeremy Magland 72 "Free edges: the moment and shear conditions are natural conditions of a discrete bending energy, met to first order at the boundary and second order inside; the three rigid-body motions are removed after every step, so the displacement is that of a plate on soft supports, not one drifting in space.",
81f3b30Elliptical plates, a WebGPU time loop for the explicit method, and WAV exportJeremy Magland 73 "Pinned centre: the support is one grid node, so the stress right at it is resolution-dependent, and the two rigid rotations about it are removed rather than restrained by a modelled washer.",
74 "Ellipse: the edge is a staircase of grid nodes, so the geometry, the edge conditions and the stress along the edge are first order in the grid spacing. The stress on the outermost ring of nodes is copied from the nearest node with a full stencil, which puts the peak stress of a clamped edge about one spacing inside. For simply supported edges the curvature term that the staircase loses is added back (the plate paradox), so the disk's frequencies converge to the exact ones.",
76 durationOf: (spec) => spec.implicit.duration,
77 gridOf: (spec) => [spec.implicit.nx, spec.implicit.ny],
78 settingsKey: "implicit",
79 usPerCellStep: 1.0,
80};
↑↓ move↵openescclose