bubbles / commit
physics rework
.github/copilot-instructions.md | 5 +- src/physics/Water.ts | 249 +++++++++++++++++++--------------------- 2 files changed, 120 insertions(+), 134 deletions(-) diff --git a/.github/copilot-instructions.md b/.github/copilot-instructions.md index 2baa8d7..e260674 100644 --- a/.github/copilot-instructions.md +++ b/.github/copilot-instructions.md @@ -30,8 +30,9 @@ scripts/ — Node.js tooling (sharp asset pipeline) - Coordinates: canvas space, origin top-left, +Y downward. - Units: pixels for distance, seconds for time, px/s for velocity, px/s² for acceleration. - `g ≈ 980 px/s²` (treats 100 px ≈ 1 metre). -- `Water` is a 1-D column simulation — never replace with SPH or particles. -- max 12 `RigidBody` objects; `World.addBody()` enforces this by dropping oldest. +- `Water` simulation approach is open — column-based, SPH particles, or any other model that runs well on mobile Canvas 2D. +- Prefer approaches that stay ≥ 30fps on a mid-range phone; profile before committing to a heavier model. +- max 12 `RigidBody` objects by default; `World.addBody()` enforces this by dropping oldest non-permanent. - Fixed physics timestep: `1/60` s (set in `App.ts`). ## Rendering Conventions diff --git a/src/physics/Water.ts b/src/physics/Water.ts index 2d43269..5e341c6 100644 --- a/src/physics/Water.ts +++ b/src/physics/Water.ts @@ -1,170 +1,155 @@ /** - * 1-D column-based shallow-water simulation. + * Water surface simulation — 1D Shallow Water Equations (SWE), staggered grid. * - * The glass interior is divided into N vertical columns. Each column tracks - * its current height (above the glass floor) and a vertical velocity. - * Adjacent columns are spring-coupled, producing wave propagation. - * A lateral gravity component shifts water mass toward the low side over time. - * If any column height exceeds the glass rim, that water is permanently lost - * (spills). + * Grid layout: + * h[i] column height at cell centre i = 0 … N−1 + * u[j] horizontal velocity at face j = 0 … N + * u[0] = u[N] = 0 (solid walls) + * + * Governing equations (linearised around mean depth H̄): + * ∂h/∂t + H̄ · ∂u/∂x = 0 (mass continuity, nonlinear flux used) + * ∂u/∂t = −g_sw · ∂h/∂x + aₓ (depth-averaged momentum) + * + * Wave speed c = √(g_sw · H̄). + * With g_sw = 980 and H̄ ≈ 150 px → c ≈ 383 px/s, + * giving a natural slosh period T ≈ 2W/c ≈ 1.6 s for a 300 px-wide glass. */ export class Water { - /** Column height array (metres — maps 1:1 to canvas pixels via scale) */ - readonly heights: Float32Array; - /** Column velocity array */ - private readonly velocities: Float32Array; - /** Total water volume (sum of all column heights * column width) — conserved unless spilling */ + /** Column heights (cell centres). */ + private readonly h: Float32Array; + /** Face velocities. Length = columns + 1; u[0] = u[N] = 0 always. */ + private readonly u: Float32Array; + /** Scratch buffer for continuity update — avoids in-place aliasing errors. */ + private readonly hNext: Float32Array; + + /** Conserved water volume [px²]. Decreases permanently on rim spill. */ private volume: number; - /** - * @param columns Number of simulation columns (default 80) - * @param glassW Glass interior width in px - * @param glassH Glass interior height in px - * @param fillFrac Initial fill fraction (0–1), default 0.55 - */ + private static readonly G_SW = 980; // shallow-water gravity [px/s²] + private static readonly K_DAMP = 1.5; // linear velocity damping coefficient [1/s] + private static readonly CFL = 0.85; // target CFL number (< 1 required) + constructor( - public columns: number = 80, + public readonly columns: number = 80, public glassW: number = 300, public glassH: number = 500, fillFrac: number = 0.55, ) { - this.heights = new Float32Array(columns).fill(glassH * fillFrac); - this.velocities = new Float32Array(columns).fill(0); - this.volume = this.computeVolume(); - } - - get columnWidth(): number { - return this.glassW / this.columns; + const h0 = glassH * fillFrac; + this.h = new Float32Array(columns).fill(h0); + this.hNext = new Float32Array(columns); + this.u = new Float32Array(columns + 1).fill(0); + this.volume = h0 * glassW; } - /** - * Advance the simulation by dt seconds. - * @param dt Time delta in seconds (typically ~0.016) - * @param lateralG Lateral gravity component (px/s²) — positive = tilt right - * @param verticalG Vertical gravity magnitude (px/s²) — positive downward - */ - step(dt: number, lateralG: number, verticalG: number): void { - const n = this.columns; - const k = 280; // spring stiffness — drives wave propagation - const kTilt = 45; // restoring strength toward equilibrium tilt (must be strong enough to actually reach it) - const damping = 0.978; // velocity damping per step - - // Equilibrium surface: when tilted, the surface should slope at angle - // slope = lateralG / verticalG (hydrostatic equilibrium). - // When the phone approaches landscape (verticalG → 0), slope → ∞ capped below. - const hAvg = this.averageHeight; - const slope = verticalG > 5 ? lateralG / verticalG : (lateralG > 0 ? 10 : -10); - // Cap slope — at slope=8 with a typical glassW the high side already overflows - const slopeCapped = Math.max(-8, Math.min(8, slope)); - - for (let i = 0; i < n; i++) { - const left = i > 0 ? this.heights[i - 1]! : this.heights[i]!; - const right = i < n - 1 ? this.heights[i + 1]! : this.heights[i]!; - const cur = this.heights[i]!; - - // Wave spring (Laplacian — drives oscillation, tends toward flat) - const springAcc = k * (left + right - 2 * cur) * dt; - - // Tilt restoring force: drives toward the equilibrium sloped surface. - // hEq for column i: average height ± slope * distance from centre. - const colCentre = (i / (n - 1)) - 0.5; // -0.5 … +0.5 - const hEq = hAvg + colCentre * this.glassW * slopeCapped; - const tiltAcc = kTilt * (hEq - cur) * dt; - - this.velocities[i] = (this.velocities[i]! + springAcc + tiltAcc) * damping; + get columnWidth(): number { return this.glassW / this.columns; } + + get averageHeight(): number { return this.volume / this.glassW; } + + step(dt: number, lateralG: number, _verticalG: number): void { + const N = this.columns; + const dx = this.glassW / N; + const g = Water.G_SW; + + // ── Adaptive sub-steps: keep CFL = dt_sub/dx * c_max < CFL_target ──── + let maxH = 1; + for (let i = 0; i < N; i++) if (this.h[i]! > maxH) maxH = this.h[i]!; + const cMax = Math.sqrt(g * maxH); + const subSteps = Math.max(1, Math.ceil(dt * cMax / (dx * Water.CFL))); + const dtS = dt / subSteps; + const dampFactor = Math.exp(-Water.K_DAMP * dtS); + + let spilled = 0; + + for (let s = 0; s < subSteps; s++) { + // ── 1. Momentum: pressure gradient + tilt body-force + damping ──── + for (let j = 1; j < N; j++) { + this.u[j] = ( + this.u[j]! + - g * (dtS / dx) * (this.h[j]! - this.h[j - 1]!) + + lateralG * dtS + ) * dampFactor; + } + // Wall BCs — no flow through side walls + this.u[0] = 0; + this.u[N] = 0; + + // ── 2. Continuity: upwind nonlinear flux into scratch buffer ────── + for (let i = 0; i < N; i++) { + // Face j = i+1 (right face of cell i) + const uR = this.u[i + 1]!; + const hR = uR >= 0 ? this.h[i]! : (i + 1 < N ? this.h[i + 1]! : this.h[i]!); + // Face j = i (left face of cell i) + const uL = this.u[i]!; + const hL = uL <= 0 ? this.h[i]! : (i > 0 ? this.h[i - 1]! : this.h[i]!); + this.hNext[i] = this.h[i]! - (dtS / dx) * (uR * hR - uL * hL); + } + // Swap h ← hNext + this.h.set(this.hNext); + + // ── 3. Clamp: floor and rim ──────────────────────────────────────── + for (let i = 0; i < N; i++) { + if (this.h[i]! < 0) { + this.h[i] = 0; + } else if (this.h[i]! > this.glassH) { + spilled += (this.h[i]! - this.glassH) * dx; + this.h[i] = this.glassH; + } + } } - // Integrate heights - for (let i = 0; i < n; i++) { - this.heights[i] = Math.max(0, this.heights[i]! + this.velocities[i]! * dt); - } + // Permanently reduce stored volume by whatever spilled over the rim + this.volume = Math.max(0, this.volume - spilled); - // Spill: clamp columns that exceed rim, lose the volume - let spilled = false; - for (let i = 0; i < n; i++) { - if (this.heights[i]! > this.glassH) { - this.heights[i] = this.glassH; - this.velocities[i] = 0; - spilled = true; + // ── 4. Volume correction: recover any transient floor-clamp loss ────── + let sum = 0; + for (let i = 0; i < N; i++) sum += this.h[i]!; + const actualVol = sum * dx; + if (actualVol > 1 && Math.abs(actualVol - this.volume) > 0.5) { + const delta = (this.volume - actualVol) / this.glassW; + for (let i = 0; i < N; i++) { + this.h[i] = Math.max(0, Math.min(this.glassH, this.h[i]! + delta)); } } - - if (spilled) { - this.volume = this.computeVolume(); - } else { - // Re-distribute to conserve volume - this.conserveVolume(); - } } - /** - * Returns the surface Y coordinate (from TOP of glass, i.e. small = high up) - * at a given X pixel position relative to the glass interior left edge. - */ getSurfaceY(x: number): number { - const col = Math.min(this.columns - 1, Math.max(0, Math.floor(x / this.columnWidth))); - return this.glassH - this.heights[col]!; - } - - /** - * Returns average surface height (px from floor) — useful for spawning objects. - */ - get averageHeight(): number { - let sum = 0; - for (let i = 0; i < this.columns; i++) sum += this.heights[i]!; - return sum / this.columns; + const col = Math.max(0, Math.min(this.columns - 1, Math.floor(x / this.columnWidth))); + return this.glassH - this.h[col]!; } /** - * Apply a velocity disturbance to columns near localX (px from glass left). - * Positive velocityY = pushing down (body entering water), negative = body leaving. - * influence controls the horizontal spread in columns. + * Seed horizontal face velocities near an impact point, creating an + * outward-spreading ripple. velocityY > 0 means downward impact. */ disturb(localX: number, radius: number, velocityY: number): void { - const colW = this.columnWidth; - const centerCol = localX / colW; - const spreadCols = (radius / colW) * 1.5; - const strength = velocityY * 0.018; // tuned so fast splashes are dramatic but not wild - - for (let i = 0; i < this.columns; i++) { - const dist = Math.abs(i - centerCol); - if (dist > spreadCols * 2) continue; - const falloff = Math.exp(-(dist * dist) / (spreadCols * spreadCols)); - this.velocities[i] = (this.velocities[i]! - strength * falloff); + const dx = this.columnWidth; + const N = this.columns; + const centerFace = localX / dx; + const spreadFaces = (radius / dx) * 1.5; + const strength = velocityY * 0.03; + + for (let j = 1; j < N; j++) { + const dist = j - centerFace; + if (Math.abs(dist) > spreadFaces * 2) continue; + const falloff = Math.exp(-(dist * dist) / (spreadFaces * spreadFaces)); + this.u[j] = this.u[j]! + strength * Math.sign(dist) * falloff; } } - /** Reset to a flat surface at the given fill fraction (default 0.55). */ resetFill(fillFrac: number = 0.55): void { - this.heights.fill(this.glassH * fillFrac); - const vel = this.velocities; - vel.fill(0); - this.volume = this.computeVolume(); + const h0 = this.glassH * fillFrac; + this.h.fill(h0); + this.u.fill(0); + this.volume = h0 * this.glassW; } - /** Resize glass dimensions (e.g. on canvas resize), scaling heights proportionally. */ resize(newGlassW: number, newGlassH: number): void { const scaleH = newGlassH / this.glassH; - for (let i = 0; i < this.columns; i++) { - this.heights[i] = this.heights[i]! * scaleH; - } + for (let i = 0; i < this.columns; i++) this.h[i] = this.h[i]! * scaleH; + this.volume = this.volume * scaleH * (newGlassW / this.glassW); this.glassW = newGlassW; this.glassH = newGlassH; - this.volume = this.computeVolume(); - } - - private computeVolume(): number { - let v = 0; - for (let i = 0; i < this.columns; i++) v += this.heights[i]!; - return v; - } - - private conserveVolume(): void { - const cur = this.computeVolume(); - if (cur <= 0) return; - const factor = this.volume / cur; - for (let i = 0; i < this.columns; i++) { - this.heights[i] = this.heights[i]! * factor; - } + this.u.fill(0); } }