reaction-diffusion.js
128 lines · in natura-src-1.0.0.zip as natura/reaction-diffusion.js · a module · its entry in the module reference · its test: reaction-diffusion-test.js
// reaction-diffusion.js -- Natura REACTION-DIFFUSION: a two-chemical field that self-organises into Turing patterns.
// Render-free + deterministic (seeded rng, no Date/Math.random). UNIVERSAL pattern-former: a leopard's spots, a
// pufferfish's maze, coral, fingerprints, travelling waves -- ONE two-species PDE tuned differently, exactly Turing's
// 1952 idea that markings are not drawn by any plan but EMERGE. A proper Natura SYSTEM (init/step/fields) -- it fits
// the framework the way weather/vegetation/ecosystem do (a stepped field on a grid). The renderer maps V to colour;
// this only evolves the chemistry. Captured from the "It Paints Itself" piece so every future piece shares one solver.
// UMD (node / window.ReactionDiffusion).
(function (root, factory) {
if (typeof module !== 'undefined' && module.exports) module.exports = factory();
else root.ReactionDiffusion = factory();
})(typeof self !== 'undefined' ? self : this, function () {
'use strict';
// The PEARSON regimes -- the same Gray-Scott law at different (feed,kill), each a different skin. "More kinds" is
// literal here: one solver, a whole atlas of patterns. (feed = U replenished, kill = V removed at feed+kill.)
const PRESETS = {
spots: { feed: 0.0367, kill: 0.0649 }, // mitosis -- dividing dots
coral: { feed: 0.0545, kill: 0.0620 }, // growing coral fingers
maze: { feed: 0.0290, kill: 0.0565 }, // fingerprint labyrinth
waves: { feed: 0.0140, kill: 0.0450 }, // travelling waves
worms: { feed: 0.0580, kill: 0.0630 }, // wandering worms
holes: { feed: 0.0390, kill: 0.0580 }, // inverse spots
};
function makeReactionDiffusion(config) {
config = config || {};
const N = config.N || 96; // grid side (torus)
const Du = config.Du != null ? config.Du : 0.16, Dv = config.Dv != null ? config.Dv : 0.08; // measured: 0.16/0.08 keeps ALL six regimes alive unattended (the (feed,kill) are Pearson-atlas values, tuned for this pair). The pieces run a hotter 1.0/0.5 -- gorgeous under continuous re-seeding but the low-feed 'waves' regime dies unattended -- so a piece that wants that exact look passes Du/Dv through config.
const p = PRESETS[config.preset] || null;
const feed0 = config.feed != null ? config.feed : (p ? p.feed : 0.0367);
const kill0 = config.kill != null ? config.kill : (p ? p.kill : 0.0649);
const seedBlobs = config.seedBlobs != null ? config.seedBlobs : 12;
const stepsPerDay = config.stepsPerDay != null ? config.stepsPerDay : 900; // Gray-Scott sub-steps advanced per sim-DAY
return {
name: config.name || 'rd',
domain: 'physical',
scale: 'micro',
params: { N: N, Du: Du, Dv: Dv, feed: feed0, kill: kill0 },
init: function (rng, params, clock) {
const n2 = N * N;
const U = new Float32Array(n2), V = new Float32Array(n2);
for (let i = 0; i < n2; i++) { U[i] = 1; V[i] = 0; }
// seed a scatter of V blobs -- the perturbation the instability grows from (deterministic via rng)
for (let b = 0; b < seedBlobs; b++) {
const cx = (rng() * N) | 0, cy = (rng() * N) | 0, r = 3 + rng() * 4, ri = Math.ceil(r); // IPI's blob size; ri is the INTEGER loop bound (float dy would index a string key and write nothing)
for (let dy = -ri; dy <= ri; dy++) for (let dx = -ri; dx <= ri; dx++) {
if (dx * dx + dy * dy > r * r) continue;
const x = ((cx + dx) % N + N) % N, y = ((cy + dy) % N + N) % N, k = y * N + x;
V[k] = 0.5; U[k] = 0.25;
}
}
return { N: N, U: U, V: V, Un: new Float32Array(n2), Vn: new Float32Array(n2), feed: feed0, kill: kill0, Du: Du, Dv: Dv };
},
// step: advance the Gray-Scott reaction-diffusion. dt is in SECONDS (the Natura convention -- World.step passes
// seconds, UNITS.DAY = 86400); convert to sim-days and run stepsPerDay sub-steps per day. Clamp the day-fraction
// like vegetation.js so a huge dt can't ask for 77M sub-steps. Double-buffered; a 9-point Laplacian (orthogonal
// 0.2, diagonal 0.05, centre -1) on the torus. feed/kill may be steered live.
step: function (st, ctx, dt) {
const n = st.N, dd = Math.max(0, Math.min(2, (dt || 0) / 86400)); // sim-days this step
const sub = Math.round(stepsPerDay * dd);
if (sub <= 0) return st;
const f = st.feed, k = st.kill, Du = st.Du, Dv = st.Dv;
let U = st.U, V = st.V, Un = st.Un, Vn = st.Vn;
for (let s = 0; s < sub; s++) {
for (let y = 0; y < n; y++) {
const yn = ((y - 1 + n) % n) * n, yp = ((y + 1) % n) * n, y0 = y * n;
for (let x = 0; x < n; x++) {
const xl = (x - 1 + n) % n, xr = (x + 1) % n, i = y0 + x;
const lu = 0.2 * (U[y0 + xl] + U[y0 + xr] + U[yn + x] + U[yp + x]) + 0.05 * (U[yn + xl] + U[yn + xr] + U[yp + xl] + U[yp + xr]) - U[i];
const lv = 0.2 * (V[y0 + xl] + V[y0 + xr] + V[yn + x] + V[yp + x]) + 0.05 * (V[yn + xl] + V[yn + xr] + V[yp + xl] + V[yp + xr]) - V[i];
const uvv = U[i] * V[i] * V[i];
let un = U[i] + (Du * lu - uvv + f * (1 - U[i]));
let vn = V[i] + (Dv * lv + uvv - (f + k) * V[i]);
Un[i] = un < 0 ? 0 : un > 1 ? 1 : un;
Vn[i] = vn < 0 ? 0 : vn > 1 ? 1 : vn;
}
}
const tU = U; U = Un; Un = tU; const tV = V; V = Vn; Vn = tV; // swap buffers
}
st.U = U; st.V = V; st.Un = Un; st.Vn = Vn;
return st;
},
fields: function (st) {
const n = st.N, V = st.V;
return {
N: n, V: V, U: st.U,
// sample the pattern (V concentration) at normalized (x,y) in [0,1) with wrap + bilinear read
at: function (x, y) {
const fx = ((x % 1) + 1) % 1 * n, fy = ((y % 1) + 1) % 1 * n;
const x0 = fx | 0, y0 = fy | 0, x1 = (x0 + 1) % n, y1 = (y0 + 1) % n, tx = fx - x0, ty = fy - y0;
const a = V[y0 * n + x0], b = V[y0 * n + x1], c = V[y1 * n + x0], d = V[y1 * n + x1];
return (a * (1 - tx) + b * tx) * (1 - ty) + (c * (1 - tx) + d * tx) * ty;
},
// the LOCAL SLOPE of V (central difference on the torus) -- what a relief render dots with a light
// vector to shade the pattern into raised ridges (It Paints Itself's signature look). {gx,gy}.
gradAt: function (x, y) {
const xi = (((x % 1) + 1) % 1 * n) | 0, yi = (((y % 1) + 1) % 1 * n) | 0;
const xl = (xi - 1 + n) % n, xr = (xi + 1) % n, yl = (yi - 1 + n) % n, yr = (yi + 1) % n;
return { gx: (V[yi * n + xr] - V[yi * n + xl]) * 0.5, gy: (V[yr * n + xi] - V[yl * n + xi]) * 0.5 };
},
};
},
events: function () { return []; },
// let a piece steer the chemistry live (move through the Pearson plane) without rebuilding
setFK: function (st, feed, kill) { st.feed = feed; st.kill = kill; },
// PAINT V into the field at normalized (nx,ny) -- the interactive drag-to-seed of It Paints Itself, and the
// way a piece injects a fresh perturbation mid-run (a fingertip, a falling drop). radius in cells, v in [0,1].
seedAt: function (st, nx, ny, radius, v) {
const n = st.N, cx = (((nx % 1) + 1) % 1 * n) | 0, cy = (((ny % 1) + 1) % 1 * n) | 0;
const r = radius != null ? radius : 4, ri = Math.ceil(r), vv = v != null ? v : 0.5;
for (let dy = -ri; dy <= ri; dy++) for (let dx = -ri; dx <= ri; dx++) {
if (dx * dx + dy * dy > r * r) continue;
const x = ((cx + dx) % n + n) % n, y = ((cy + dy) % n + n) % n, k = y * n + x;
st.V[k] = vv; st.U[k] = 1 - vv;
}
},
};
}
return { makeReactionDiffusion: makeReactionDiffusion, PRESETS: PRESETS };
});