/ listening-for-chaos / chaos-apps
Sign in
listening-for-chaos / chaos-apps
Driven plate: optional loss, per-mode decay concentrated at the top of the band (rate ~ omega^k, k 1, 2 or 4, tau top in driven periods), switchable during a run; drive and loss power readout; mode map colour bar on its own line
Jeremy Magland <jmagland@flatironinstitute.org> committed commit b22731fd7636 parent bedc4d4 Browse files
5 changed files+108−11
src/driven-plate/README.mdmodified+7−1View file
@@ -1,6 +1,6 @@
11 # driven-plate
22
3-A rectangular plate driven in one mode, without damping, its energy raised slowly, to watch the energy move from the driven mode into the others. Step 1 (2026-10-09): simply supported Nx × Ny plate with von Kármán stretching, modal drive with an energy-target controller, log ramp, height map, energy chart. Step 2 (2026-10-09): mode map (energy of every mode against time).
3+A rectangular plate driven in one mode, without damping, its energy raised slowly, to watch the energy move from the driven mode into the others. Step 1 (2026-10-09): simply supported Nx × Ny plate with von Kármán stretching, modal drive with an energy-target controller, log ramp, height map, energy chart. Step 2 (2026-10-09): mode map (energy of every mode against time). Loss option (2026-10-09): a per-mode sink at the top of the band.
44
55 Plate (`plate.ts`): the square plate of `../plate` (simply supported, flat) generalized to a rectangle, with spacings ax and ay. Nx follows from the masses across the short side and the aspect so the cells are nearly square (default 8 across, aspect 1.3: 11 × 8). Units scaled away (D = ρh = h = 1, short side 1): with energy given against the nonlinear scale, displacement against h and time in driven periods, only the aspect and ν = 0.3 remain. Bending on the five-point Laplacian, exact sine modes ω_mn = l_m / ax² + l_n / ay²; stretching by the stress function, L⁻² by two sine transforms (plain matrix products, no mirror trick yet), force the exact gradient.
66
@@ -15,3 +15,9 @@ First look (default settings): the other modes' energy rises steadily with the d
1515 Mode map (step 2, `history.ts`, `modemap.ts`): `drive.ts` samples every mode's linear energy (sine transforms of u and p) 16 times per driven period and averages them; the page closes a time bin every 0.5 periods (`History`), merging neighbours in pairs and doubling the bin length when there are 1200, so a run of any length fits. The energy chart now reads the same bins (smooth: the per-frame sampling of the swinging off-resonance energies made the other-modes curve a band). The map has one row per mode, ordered by pitch, labelled (m, n) and pitch relative to the driven mode, colour log energy (inferno-like) from 12 decades below the largest so far; new bins are painted incrementally, the whole image only on a merge or a change of the colour scale's top decade. Hover reads a mode and its energy.
1616
1717 Seen (default settings, 2026-10-09; the map, and a bun script listing the leading modes every 100 periods of the default ramp): up to about 0.2 E_nl the other modes that hold energy are the (odd, odd) family of (1, 1), (3, 1) first, then (1, 3), (3, 3), (5, 1), (1, 5), (5, 3), their energies growing as about E³ (the driven mode's own nonlinear response, bound to it); the rest stay at the perturbation's level. Near 0.3 E_nl (1, 3) takes off and passes (3, 1) (6e-3 E_nl at 0.31, 0.2 at 0.56), then (1, 5) (0.1 at 1 E_nl). With the perturbation, modes of other symmetry ((2, 5), (4, 3), (4, 5)) lead by 1.8 E_nl; without it the leaders at 1.8 E_nl are still all (odd, odd) ((1, 5), (3, 7), (5, 3), (5, 7)), up to then identical to the perturbed run. On the map the whole band brightens within a few tens of periods around 1100 periods (about 0.5 E_nl) and stays filled to the grid's top mode (no sink: undamped).
18+
19+Loss (`drive.ts`, Plate and drive card; off by default, switchable during a run): each mode's momentum decays at γ_k = (2 / τ top)(ω_k / ω_top)^k, so its amplitude decays in τ top (ω_top / ω_k)^k (τ top in driven periods, default 10; k = 1, 2 (default) or 4). Applied exactly in the sine basis every 8 steps (two transforms, about 10% extra per step), kinetic energy removed summed in Wl. Not the dashpot loss of `../plate` (force −c L v): on a plate that gives k = 1, which also damps the driven mode noticeably (decay 40 τ top on the default plate, against 1700 τ top with k = 2), so the per-mode form, which can concentrate the loss at the top, was chosen. The readout shows drive and loss power (E_nl per driven period, smoothed over about 20 periods); the energy error is now relative to the larger of E and the energy that has flowed through (Estart + |Wd| + Wl), since with a loss the error grows with the work done.
20+
21+Checks (2026-10-09): decay times of (11, 8) and (6, 4) at tiny amplitude, drive off, against prediction: 9.97 / 10.00, 21.27 / 21.34 (k = 1), 45.39 / 45.53 (k = 2), 207.0 / 207.3 (k = 4) periods (the first try gave twice the prediction: damping only the momentum halves the amplitude's rate, now in γ). Default ramp with the loss (τ top 10, k = 2): energy error 4e-3 at 100 E_nl relative to E (where ~150 E_nl had been dissipated in the last 100 periods), drive and loss power balance at the top (3.7 E_nl per period at 100 E_nl).
22+
23+With and without the loss (2026-10-09; default plate and ramp, then held 100 + 200 periods; mean energy per mode in octaves of pitch above the driven mode, 2-4, 4-8, 8-16, 16-32, 32-64): at 10 E_nl without loss nearly flat, 0.06, 0.08, 0.11, 0.09, 0.09 E_nl (close to equipartition, the pile-up), the driven mode 0.46; with the loss (τ top 10, k = 2) falling, 0.30, 0.19, 0.067, 0.035, 0.022 (about ω⁻¹ per mode), the driven mode 3.0. At 1 E_nl with the loss the fall is steep: 1.4e-2, 5.5e-3, 1.0e-3, 1.4e-4, 1.6e-6 (without: 9.6e-3, 1.6e-2, 1.3e-2, 1.6e-3, 9.9e-4, not yet spread to the top). The onset is the same with or without the loss at this τ top: (1, 3) takes off near 0.3 E_nl in both.
src/driven-plate/drive.tsmodified+59−2View file
@@ -15,6 +15,14 @@
1515 // gentle when the driven mode holds little of the energy. G is set every
1616 // control interval from E; Q is held fixed over each step, and its work Q dq
1717 // is added to Wd, so E - Wd stays constant up to the integrator's error.
18+// Loss (optional, off by default): each mode's momentum decays at its own
19+// rate, gamma_k = (2 / tauTop) (omega_k / omega_top)^power, so the mode's
20+// amplitude decays at gamma_k / 2 (half its time is kinetic): decay time
21+// tauTop for the top mode, longer below. Applied exactly in the sine
22+// basis (p -> S diag(exp(-gamma_k dt)) S p) every LOSS_EVERY steps, as a
23+// sink at the top of the band. power = 1 is what the dashpot loss of
24+// ../plate gives (rate proportional to omega for a plate); 2 and 4 leave the
25+// low modes nearly untouched. The kinetic energy removed goes into Wl.
1826 // Time: velocity Verlet at the step of plate.ts, recomputed from the energy.
1927 // Mode energies for the displays are averaged over time (sampled every
2028 // 1 / SAMPLES_PER_PERIOD of a driven period): a mode's linear energy swings
@@ -26,6 +34,12 @@ const TAU_C = 2; // controller time constant, in driven periods
2634 const G_MAX = 0.2; // cap on G, times the driven mode's omega
2735 const PERTURB = 1e-6; // energy of the perturbation, relative to the starting energy
2836 const SAMPLES_PER_PERIOD = 16;
37+const LOSS_EVERY = 8;
38+
39+export interface Loss {
40+ tauTop: number; // amplitude decay time of the top mode, in driven periods; 0 = off
41+ power: number; // gamma ~ omega^power
42+}
2943
3044 export interface DriveOptions {
3145 m: number; // driven mode (m, n)
@@ -49,6 +63,7 @@ export class DrivenPlate {
4963 E = 0; // energy now (kinetic + potential)
5064 Et: number; // target energy (absolute units)
5165 Wd = 0; // work done by the drive
66+ Wl = 0; // energy removed by the loss
5267 Estart: number; // energy at the start, for the energy error
5368 G = 0; // controller gain
5469 driveOn = true; // false: no drive (G = 0)
@@ -150,6 +165,7 @@ export class DrivenPlate {
150165 const { u, p, F, S } = this, n = this.g.n, M = this.g.mass, sqM = this.sqM;
151166 // Gain update every control interval (a few per driven period).
152167 const every = Math.max(1, Math.round(this.period / 32 / h));
168+ let pending = 0; // steps since the loss was last applied
153169 for (let s = 0; s < steps; s++) {
154170 if (s % every === 0) {
155171 this.E = this.kinetic() + this.V;
@@ -162,6 +178,10 @@ export class DrivenPlate {
162178 for (let k = 0; k < n; k++) p[k] += 0.5 * h * (F[k] + f * S[k]);
163179 this.readMode();
164180 this.Wd += Q * (this.q - qBefore);
181+ if (this.gamma && ++pending === LOSS_EVERY) {
182+ this.applyLoss(pending * h);
183+ pending = 0;
184+ }
165185 if (qBefore < 0 && this.q >= 0) {
166186 const tc = this.t + (h * -qBefore) / (this.q - qBefore);
167187 if (Number.isFinite(this.lastCross)) this.measuredPeriod = tc - this.lastCross;
@@ -175,13 +195,50 @@ export class DrivenPlate {
175195 this.nextSample = this.t + this.period / SAMPLES_PER_PERIOD;
176196 }
177197 }
198+ if (this.gamma && pending > 0) this.applyLoss(pending * h);
178199 this.E = this.kinetic() + this.V;
179200 return steps;
180201 }
181202
182- // Relative energy error: (E - Estart - Wd) / E.
203+ // The loss: decay rates per mode (grid index), null when off.
204+ private gamma: Float64Array | null = null;
205+ loss: Loss = { tauTop: 0, power: 2 };
206+ setLoss(loss: Loss): void {
207+ this.loss = loss;
208+ if (!(loss.tauTop > 0)) {
209+ this.gamma = null;
210+ return;
211+ }
212+ const g = this.g, top = this.modes.omegaMax, rate = 2 / (loss.tauTop * this.period);
213+ this.gamma = new Float64Array(g.n);
214+ for (let b = 0; b < g.Ny; b++) for (let a = 0; a < g.Nx; a++) this.gamma[a + g.Nx * b] = rate * (omegaOf(g, a + 1, b + 1) / top) ** loss.power;
215+ }
216+ // Amplitude decay time of mode (m, n) from the loss, in driven periods
217+ // (Infinity when off).
218+ decayTime(m: number, n: number): number {
219+ if (!(this.loss.tauTop > 0)) return Infinity;
220+ return this.loss.tauTop * (this.modes.omegaMax / omegaOf(this.g, m, n)) ** this.loss.power;
221+ }
222+ private lossHat: Float64Array | null = null;
223+ private applyLoss(dt: number): void {
224+ const gam = this.gamma!, n = this.g.n, P = (this.lossHat ??= new Float64Array(n));
225+ this.plate.sine2(this.p, P);
226+ let dT = 0;
227+ for (let k = 0; k < n; k++) {
228+ const f = Math.exp(-gam[k] * dt), p1 = P[k] * f;
229+ dT += P[k] * P[k] - p1 * p1;
230+ P[k] = p1;
231+ }
232+ this.plate.sine2(P, this.p);
233+ this.Wl += dT / (2 * this.g.mass);
234+ this.readMode();
235+ }
236+
237+ // Relative energy error: E - Estart - Wd + Wl, relative to the energy
238+ // handled (E, or with the loss on, the larger energy that has flowed
239+ // through: the error grows with the work done, not with E).
183240 get energyError(): number {
184- return (this.E - this.Estart - this.Wd) / this.E;
241+ return (this.E - this.Estart - this.Wd + this.Wl) / Math.max(this.E, this.Estart + Math.abs(this.Wd) + this.Wl);
185242 }
186243 get step(): number {
187244 return this.h;
src/driven-plate/index.htmlmodified+12−0View file
@@ -47,6 +47,17 @@
4747 <select id="mode"></select>
4848 <label class="chip" title="Start with a tiny random motion (a millionth of the starting energy) in every other mode; off, the plate starts in the driven mode exactly"><input id="perturb" type="checkbox" checked /><span>Perturbation</span></label>
4949 </div>
50+ <div class="row">
51+ <label class="chip" title="A loss that drains the highest modes, as a sink for the energy that reaches the top of the band"><input id="loss" type="checkbox" /><span>Loss</span></label>
52+ <span class="muted">τ top</span>
53+ <input id="loss-tau" type="number" min="0.1" step="any" value="10" />
54+ <span class="muted">periods, rate ∝</span>
55+ <select id="loss-power">
56+ <option value="1">ω</option>
57+ <option value="2" selected>ω²</option>
58+ <option value="4">ω⁴</option>
59+ </select>
60+ </div>
5061 <div id="params"></div>
5162 </section>
5263 <section class="card">
@@ -75,6 +86,7 @@
7586 <h2>Plate and drive</h2>
7687 <p>Masses across is the number of masses across the short side; the long side gets as many as keep the cells nearly square. Aspect is the ratio of the sides. A square plate has pairs of modes at the same pitch ((1, 2) and (2, 1)), which trade energy easily; a non-square plate splits them. The default 1.3 keeps the lowest modes at least 1% apart in pitch.</p>
7788 <p>The modes are products of sines, numbered (m, n) by their half waves along x and y. The list gives each one's pitch relative to the lowest. Perturbation starts the plate with a tiny random motion in every other mode, a millionth of the starting energy in all; without it the plate starts in the driven mode exactly, and by symmetry only modes with the driven mode's symmetry can then grow (for (1, 1): odd m and odd n), at first; round-off seeds the rest. Any change here starts again.</p>
89+ <p>Loss drains energy from the modes, fastest at the top of the band, as a sink for the energy the cascade brings there (without it the energy has nowhere to go and spreads over every mode up to the grid's highest). Each mode's amplitude decays on its own, the top mode in τ top driven periods and a mode of pitch ω in τ top (ω<sub>top</sub> / ω)<sup>k</sup>, k from the rate's power of ω: with ω (k = 1, the dashpot loss of the other apps) the driven mode decays in about 40 τ top on the default plate, with ω² in about 1700 τ top, with ω⁴ practically never. The drive makes up what the loss takes. The loss can be switched and changed during a run.</p>
7890 <p>Units are scaled away: for the von Kármán plate only the aspect ratio and Poisson's ratio (0.3) matter once energy is given against E<sub>nl</sub>, displacement against the thickness h and time in periods.</p>
7991 </template>
8092 <template id="help-energy">
src/driven-plate/main.tsmodified+22−1View file
@@ -30,6 +30,9 @@ const gridInfo = $("grid-info");
3030 const modeSelect = $<HTMLSelectElement>("mode");
3131 const perturbBox = $<HTMLInputElement>("perturb");
3232 const paramsDiv = $("params");
33+const lossBox = $<HTMLInputElement>("loss");
34+const lossTauInput = $<HTMLInputElement>("loss-tau");
35+const lossPowerSelect = $<HTMLSelectElement>("loss-power");
3336
3437 let P: Params = { Ny: 8, aspect: 1.3 };
3538 let sim: DrivenPlate;
@@ -38,6 +41,8 @@ let logTarget = LOG_START;
3841 let actualSpeed = 0; // periods per second, smoothed
3942 let ampl = 0; // colour scale: recent largest |u|
4043 let hist: History;
44+// Drive and loss power (E_nl per driven period), smoothed.
45+let drivePower = 0, lossPower = 0;
4146
4247 const fmt = (x: number, d = 3) => (x === 0 ? "0" : Math.abs(x) >= 1e-2 && Math.abs(x) < 1e4 ? x.toPrecision(d) : x.toExponential(d - 1));
4348 const speed = () => 10 ** Number(speedSlider.value);
@@ -62,6 +67,8 @@ function reset(): void {
6267 sim = new DrivenPlate(P, { m, n, E0: 10 ** LOG_START, perturb: perturbBox.checked });
6368 logTarget = LOG_START;
6469 targetSlider.value = String(LOG_START);
70+ applyLoss();
71+ drivePower = lossPower = 0;
6572 hist = new History(sim);
6673 modeMap.reset(sim);
6774 ampl = 0;
@@ -74,6 +81,9 @@ function showParams(): void {
7481 `grid ${g.Nx} × ${g.Ny}, sides ${g.Lx.toFixed(2)} × ${g.Ly.toFixed(2)} (short side 1), ν 0.3`,
7582 `driven (${sim.opts.m}, ${sim.opts.n}), pitch ${(sim.omega / M.omega[0]).toFixed(3)} × lowest`,
7683 `top mode ${(M.omegaMax / sim.omega).toFixed(1)} × driven`,
84+ sim.loss.tauTop > 0
85+ ? `loss: decay ${fmt(sim.loss.tauTop)} periods at the top, ${fmt(sim.decayTime(M.m[M.count >> 1], M.n[M.count >> 1]))} mid-band, ${fmt(sim.decayTime(sim.opts.m, sim.opts.n))} driven mode`
86+ : "loss off",
7787 `lowest: ${Array.from(M.omega.slice(0, 8), (w, k) => `(${M.m[k]},${M.n[k]}) ${(w / M.omega[0]).toFixed(2)}`).join(" ")}`,
7888 ];
7989 paramsDiv.textContent = lines.join("\n");
@@ -109,8 +119,11 @@ function frame(now: number): void {
109119 }
110120 setTarget();
111121 if (periods > 0) {
112- const t0 = performance.now();
122+ const t0 = performance.now(), Wd = sim.Wd, Wl = sim.Wl;
113123 sim.advance(periods * sim.period);
124+ const a = Math.min(1, periods / 20); // smooth over ~20 periods
125+ drivePower += a * ((sim.Wd - Wd) / sim.Enl / periods - drivePower);
126+ lossPower += a * ((sim.Wl - Wl) / sim.Enl / periods - lossPower);
114127 const ms = performance.now() - t0;
115128 msPerPeriod = 0.8 * msPerPeriod + 0.2 * Math.max(ms / periods, 1e-3);
116129 }
@@ -131,6 +144,7 @@ function draw(): void {
131144 readout.textContent =
132145 `t ${(sim.t / sim.period).toFixed(0)} periods · E ${fmt(E)} E_nl (target ${fmt(sim.Et / sim.Enl)}) · driven mode ${fmt((100 * driven) / E, 3)}% · other modes ${fmt((100 * others) / E, 2)}%\n` +
133146 `driven pitch ${Number.isFinite(T) ? fmt(1 / T) : "–"} × linear · max |w| ${fmt(sim.maxAbsU(), 2)} h · ` +
147+ (sim.loss.tauTop > 0 ? `drive ${fmt(drivePower, 2)}, loss ${fmt(lossPower, 2)} E_nl/period · ` : "") +
134148 `speed ${running ? fmt(actualSpeed, 2) : "0"} periods/s · energy error ${Math.abs(sim.energyError).toExponential(1)}`;
135149 }
136150
@@ -293,6 +307,13 @@ function applyGrid(): void {
293307 nyInput.addEventListener("change", applyGrid);
294308 aspectInput.addEventListener("change", applyGrid);
295309 modeSelect.addEventListener("change", reset);
310+// Loss: applies to the running plate.
311+function applyLoss(): void {
312+ const tau = Number(lossTauInput.value);
313+ sim.setLoss({ tauTop: lossBox.checked && tau > 0 ? tau : 0, power: Number(lossPowerSelect.value) });
314+ showParams();
315+}
316+for (const el of [lossBox, lossTauInput, lossPowerSelect]) el.addEventListener("change", applyLoss);
296317 perturbBox.addEventListener("change", reset);
297318
298319 fillModes();
src/driven-plate/modemap.tsmodified+8−7View file
@@ -11,7 +11,7 @@ import type { DrivenPlate } from "./drive";
1111 import { MAX_BINS, type History } from "./history";
1212
1313 const DECADES = 12;
14-const L = 92, R = 8, T = 30, B = 22; // margins (px)
14+const L = 92, R = 8, T = 44, B = 22; // margins (px)
1515
1616 // Inferno-like colour map, dark (low) to pale yellow (high), as a lookup table.
1717 const STOPS = [
@@ -119,16 +119,17 @@ export class ModeMap {
119119 ctx.fillStyle = "#7a7d82";
120120 ctx.textAlign = "center";
121121 ctx.fillText(`time (driven periods), 0 to ${tMax.toFixed(0)}`, L + pw / 2, Hc - 6);
122- const cbW = Math.min(160, pw / 3), cbX = L + pw - cbW - 66;
122+ // Colour bar on its own line under the hover text, right-aligned.
123+ const cbW = Math.min(160, pw / 3), cbX = L + pw - cbW - 66, cbY = 24;
123124 for (let i = 0; i < cbW; i++) {
124125 ctx.fillStyle = cssColor(i / (cbW - 1));
125- ctx.fillRect(cbX + i, 8, 1.5, 10);
126+ ctx.fillRect(cbX + i, cbY, 1.5, 10);
126127 }
127128 ctx.fillStyle = "#7a7d82";
128129 ctx.textAlign = "right";
129- ctx.fillText(`1e${top - DECADES}`, cbX - 4, 17);
130+ ctx.fillText(`1e${top - DECADES}`, cbX - 4, cbY + 9);
130131 ctx.textAlign = "left";
131- ctx.fillText(`1e${top} E_nl`, cbX + cbW + 4, 17);
132+ ctx.fillText(`1e${top} E_nl`, cbX + cbW + 4, cbY + 9);
132133
133134 // Hover: the mode and bin under the cursor.
134135 const h = this.hover;
@@ -138,12 +139,12 @@ export class ModeMap {
138139 const c = Math.min(bins.length - 1, Math.floor(((h.x - L) / pw) * bins.length));
139140 const r = Math.min(count - 1, Math.floor((T + ph - h.y) / rowH));
140141 const e = bins[c].modes[this.rows[r]];
141- ctx.fillText(`(${M.m[r]}, ${M.n[r]}) · ${(M.omega[r] / sim.omega).toFixed(2)} × driven · t ${bins[c].t.toFixed(0)} · ${e.toExponential(1)} E_nl`, L, 17);
142+ ctx.fillText(`(${M.m[r]}, ${M.n[r]}) · ${(M.omega[r] / sim.omega).toFixed(2)} × driven · t ${bins[c].t.toFixed(0)} · ${e.toExponential(1)} E_nl`, L, 13);
142143 ctx.strokeStyle = "rgba(255, 255, 255, 0.7)";
143144 ctx.strokeRect(L, yOf(r) - rowH / 2, pw, rowH);
144145 } else {
145146 ctx.fillStyle = "#7a7d82";
146- ctx.fillText("rows: modes by pitch (× driven); hover to read", L, 17);
147+ ctx.fillText("rows: modes by pitch (× driven); hover to read", L, 13);
147148 }
148149 }
149150 }
↑↓ move↵openescclose