weather.js
159 lines · in natura-src-1.0.0.zip as natura/weather.js · a module · its entry in the module reference · its test: weather-test.js
// weather.js -- the first real NATURA system (DESIGN.md step 1): dynamic weather that CHANGES THE
// WORLD IT FALLS ON. A grid atmosphere over a terrain: drifting pressure systems drive a circulating
// wind (geostrophic), which advects TEMPERATURE and HUMIDITY (Stam semi-Lagrangian backtrace -- the
// fluids-and-flow recipe); air forced UP over mountains cools and rains (orographic), leaving the lee
// dry (rain-shadow emerges, it is not scripted); rain ACCUMULATES into a moisture field other systems
// (erosion, vegetation) read -- the coupling that makes weather change the ground.
//
// A Natura system spec (init/step/fields/events), deterministic + node-testable. Reads the clock for
// season + day/night. Terrain comes in as an elevation sampler so weather stands alone and tests alone;
// wiring the real Terrain Engine as the coupled terrain system is a thin later step.
// UMD (node test / window.Weather).
(function (root, factory) {
if (typeof module !== 'undefined' && module.exports) module.exports = factory();
else root.Weather = factory();
})(typeof self !== 'undefined' ? self : this, function () {
'use strict';
const clamp = (v, a, b) => v < a ? a : v > b ? b : v;
// bilinear sample of a grid field at fractional cell (fi,fj), clamped to the edges
function sample(F, N, fi, fj) {
fi = clamp(fi, 0, N - 1.001); fj = clamp(fj, 0, N - 1.001);
const i0 = fi | 0, j0 = fj | 0, i1 = i0 + 1, j1 = j0 + 1, tx = fi - i0, ty = fj - j0;
const a = F[j0 * N + i0], b = F[j0 * N + i1], c = F[j1 * N + i0], d = F[j1 * N + i1];
return a * (1 - tx) * (1 - ty) + b * tx * (1 - ty) + c * (1 - tx) * ty + d * tx * ty;
}
// saturation humidity rises with temperature (warm air holds more water) -- a soft Clausius-ish curve
const satOf = (T) => 0.15 + 0.7 * Math.max(0, T) / (Math.max(0, T) + 10); // saturation well above typical humidity -> clear sky is the default
// make a Weather system spec. config: { N, seaLevel, lapse, elevation(x,y)->e, windK, evap, systems }
function makeWeather(config) {
config = config || {};
const N = config.N || 48;
const seaLevel = config.seaLevel != null ? config.seaLevel : 0;
const lapse = config.lapse != null ? config.lapse : 26; // deg cooling per unit elevation
const windK = config.windK != null ? config.windK : 0.9; // pressure-gradient -> wind gain
const evapK = config.evap != null ? config.evap : 0.9; // sea evaporation gain
const prevU = config.prevailU != null ? config.prevailU : N * 0.06; // a steady prevailing wind (the westerlies)
const prevV = config.prevailV != null ? config.prevailV : N * 0.012; // -- gives a persistent rain-shadow
const nSys = config.systems || 4; // moving pressure systems
const latGrad = config.latGradient != null ? config.latGradient : 0; // deg C cooler from equator (y=0) to pole (|y|=1)
const elevFn = config.elevation || ((x, y) => 0.2); // default: uniform low land
const gx = (i) => (i / (N - 1)) * 2 - 1; // cell -> world [-1,1]
// TORUS support (config.wrap): a boundless atmosphere. Advection, wind gradients, and the pressure
// systems wrap at the edges, so a storm leaving one side re-enters the other and the rain-shadow
// zones are seamless -- the "nearly infinite space" the world needs.
const wrap = !!config.wrap;
const wgi = (i) => wrap ? ((i % N) + N) % N : clamp(i, 0, N - 1); // wrap/clamp a grid index
const td2 = (a, b) => { let d = a - b; if (wrap) { if (d > 1) d -= 2; else if (d < -1) d += 2; } return d; }; // toroidal world delta (width 2)
function smpW(F, fi, fj) { // wrap-aware bilinear sample (else edge-clamped)
if (!wrap) return sample(F, N, fi, fj);
fi = ((fi % N) + N) % N; fj = ((fj % N) + N) % N;
const i0 = fi | 0, j0 = fj | 0, i1 = (i0 + 1) % N, j1 = (j0 + 1) % N, tx = fi - i0, ty = fj - j0;
const a = F[j0 * N + i0], b = F[j0 * N + i1], c = F[j1 * N + i0], d = F[j1 * N + i1];
return a * (1 - tx) * (1 - ty) + b * tx * (1 - ty) + c * (1 - tx) * ty + d * tx * ty;
}
return {
name: config.name || 'weather',
domain: 'physical',
scale: 'weather',
params: config,
init(rng) {
const elev = new Float32Array(N * N), land = new Uint8Array(N * N);
for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) {
const e = elevFn(gx(i), gx(j)); elev[j * N + i] = e; land[j * N + i] = e > seaLevel ? 1 : 0;
}
const sys = [];
for (let k = 0; k < nSys; k++) sys.push({
x: rng() * 2 - 1, y: rng() * 2 - 1,
vx: (rng() - 0.5) * 0.05, vy: (rng() - 0.5) * 0.05, // slow drift (world units / sec-ish, scaled by dt)
str: (rng() < 0.5 ? -1 : 1) * (0.6 + rng() * 0.8), // low (-) or high (+)
rad: 0.35 + rng() * 0.4, age: 0, life: 3 + rng() * 6
});
const z = () => new Float32Array(N * N);
const st = { N, elev, land, sys, rng,
temp: z(), hum: z(), rain: z(), cloud: z(), moist: z(), windU: z(), windV: z(), pres: z() };
// seed temperature/humidity to a plausible field
for (let k = 0; k < N * N; k++) { st.temp[k] = 15 - lapse * Math.max(0, elev[k]) - latGrad * Math.abs(gx((k / N) | 0)); st.hum[k] = land[k] ? 0.08 : 0.28; }
return st;
},
step(st, ctx, dt) {
const s = clamp(dt / 3600, 0, 1.2); // sim advances in ~hour-ish steps
const { N, elev, land } = st;
// 1. drift the pressure systems; retire the old, spawn fresh -> dynamic weather over days
for (const p of st.sys) {
p.x += p.vx * s * 6; p.y += p.vy * s * 6; p.age += s / 24;
if (wrap) { if (p.x > 1) p.x -= 2; else if (p.x < -1) p.x += 2; if (p.y > 1) p.y -= 2; else if (p.y < -1) p.y += 2; }
else { if (p.x < -1.3 || p.x > 1.3) p.vx *= -1; if (p.y < -1.3 || p.y > 1.3) p.vy *= -1; }
if (p.age > p.life) { p.x = st.rng() * 2 - 1; p.y = st.rng() * 2 - 1; p.str = (st.rng() < 0.5 ? -1 : 1) * (0.6 + st.rng() * 0.8); p.age = 0; p.life = 3 + st.rng() * 6; }
}
// 2. pressure field + 3. geostrophic wind (circulates lows/highs) -------------------------------
const P = st.pres;
for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) {
let pr = 0; const x = gx(i), y = gx(j);
for (const p of st.sys) { const dx = td2(x, p.x), dy = td2(y, p.y); pr += p.str * Math.exp(-(dx * dx + dy * dy) / (p.rad * p.rad)); }
P[j * N + i] = pr;
}
const seasonT = 15 + 10 * Math.sin((ctx.clock.yearPhase() - 0.15) * 2 * Math.PI); // warm summer, cold winter
const dayT = 4 * (ctx.clock.daylight() - 0.4); // warmer by day
const nT = new Float32Array(N * N), nH = new Float32Array(N * N);
for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) {
const k = j * N + i;
// wind = 90deg-rotated pressure gradient (geostrophic): circulates around systems
const il = wgi(i - 1), ir = wgi(i + 1), jd = wgi(j - 1), ju = wgi(j + 1);
const dPx = (P[j * N + ir] - P[j * N + il]) * 0.5, dPy = (P[ju * N + i] - P[jd * N + i]) * 0.5;
st.windU[k] = prevU - dPy * windK * N * 0.5; st.windV[k] = prevV + dPx * windK * N * 0.5; // prevailing + geostrophic (grid-cells/sec-ish)
}
// 4+5. ADVECT temp & humidity along the wind (semi-Lagrangian backtrace -- Stam) -----------------
for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) {
const k = j * N + i;
const fi = i - st.windU[k] * s, fj = j - st.windV[k] * s;
const baseT = seasonT + dayT - lapse * Math.max(0, elev[k]) - latGrad * Math.abs(gx(j));
nT[k] = smpW(st.temp, fi, fj) * 0.82 + baseT * 0.18; // advect + relax to equilibrium
let h = smpW(st.hum, fi, fj);
if (!land[k]) h += evapK * s * clamp(st.temp[k] / 25, 0.2, 1.4) * 0.07; // sea evaporates (warmer -> more)
// SINK: a base drizzle/loss PLUS subsidence drying under HIGH pressure (descending air is dry) --
// this is what keeps clear sky the default and confines cloud to the lows/fronts, not global haze.
nH[k] = h * (1 - (0.022 + Math.max(0, P[k]) * 0.05) * s);
}
st.temp = nT;
// 6. condensation -> RAIN: saturation excess + orographic uplift (wind blowing UPHILL) -----------
const rain = st.rain;
for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) {
const k = j * N + i;
const il = wgi(i - 1), ir = wgi(i + 1), jd = wgi(j - 1), ju = wgi(j + 1);
const dEx = (elev[j * N + ir] - elev[j * N + il]) * 0.5, dEy = (elev[ju * N + i] - elev[jd * N + i]) * 0.5;
const uphill = Math.max(0, st.windU[k] * dEx + st.windV[k] * dEy); // air forced up the slope
const sat = satOf(nT[k]);
let cond = Math.max(0, nH[k] - sat) * 0.85 + uphill * nH[k] * 0.08; // saturation + orographic
cond = Math.min(cond, nH[k]);
nH[k] -= cond;
rain[k] = cond;
st.moist[k] += cond * s; // ACCUMULATE -> couples to the ground
st.cloud[k] = clamp((nH[k] / (sat + 1e-3) - 0.6) / 0.4, 0, 1); // cloud only as humidity nears saturation (sharp onset)
}
st.hum = nH;
return st;
},
// what other systems read: rain (instantaneous), moisture (accumulated), wind, temp, cloud, pressure
fields(st) {
return { N: st.N, rain: st.rain, moisture: st.moist, temp: st.temp,
cloud: st.cloud, pressure: st.pres, windU: st.windU, windV: st.windV, land: st.land };
},
// a storm event when a strong low is well inside the domain
events(st) {
const out = [];
for (const p of st.sys) if (p.str < -1.0 && Math.abs(p.x) < 0.9 && Math.abs(p.y) < 0.9) out.push({ kind: 'storm', x: +p.x.toFixed(2), y: +p.y.toFixed(2), str: +p.str.toFixed(2) });
return out;
}
};
}
return { makeWeather, satOf };
});