From a3d47448f16534201b9e819704d393fdb4fa2fe9 Mon Sep 17 00:00:00 2001 From: Dmitrii Vasilev Date: Sun, 18 Oct 2026 14:36:40 +0700 Subject: [PATCH] spec(t27f): particle_field -- agent events as physical particles; sacred_lattice -- the formula as a codec with a chance label (Closes #10289) particle_field.t27 is the physical layer between the T27F Live work graph and the scene: nodes are particles, edges damped springs plus pairwise repulsion, every agent event (user message, think, tool call and result, agent message, error, done) has a ternary charge and a hue and kicks the receiving node. Fixed DT = 1/60, symplectic steps, seeded LCG jitter, so a replay is bit-identical; reduced motion turns kicks and jitter off. Output is sized for 64-vertex meshlets. Integrator checks: harmonic period to 1e-4, Kepler energy drift under 1e-5 (Euler over 1e-2), spring-repulsion balance against Newton, Bohr speed and radius from CODATA 2022 to 1e-9. sacred_lattice.t27 keeps V = n * 3^k * pi^m * phi^p * e^q as an address (five integers packed in one i64) and measures fits against a null model: at complexity 6 the lattice fits 62% of random numbers within 0.1%, and 42 of 64 seeded random targets vs 6 of 10 real constants. No fit is labelled above UNKNOWN unless its chance times the trials is under 5%. Lab (Railway, t27c test-report): 14/14 and 10/10, 0 vacuous; three planted faults (a flipped charge, a halved density, Euler for Verlet) each caught. Co-Authored-By: Claude Opus 5.5 --- specs/physics/sacred_lattice.t27 | 296 ++++++++++++++++++ specs/t27f/particle_field.t27 | 512 +++++++++++++++++++++++++++++++ 2 files changed, 808 insertions(+) create mode 100644 specs/physics/sacred_lattice.t27 create mode 100644 specs/t27f/particle_field.t27 diff --git a/specs/physics/sacred_lattice.t27 b/specs/physics/sacred_lattice.t27 new file mode 100644 index 0000000000..9ac342ae09 --- /dev/null +++ b/specs/physics/sacred_lattice.t27 @@ -0,0 +1,296 @@ +// SPDX-License-Identifier: Apache-2.0 +// specs/physics/sacred_lattice.t27 -- the formula V = n * 3^k * pi^m * phi^p * e^q as an ADDRESS of a +// number, with the chance of every fit measured, not assumed. +// +// THE MODEL (D. Vasilev, S. Pellis, "Low-Complexity Algebraic Representations of Physical Constants", +// April 2026, eq. 1-2 and 4): for n in 1..9 and integers k, m, p, q with |k| + |m| + |p| + |q| <= C, +// ln V = ln n + k ln 3 + m ln pi + p ln phi + q. +// So the grammar is a lattice in log space, and the best fit of a target X is the nearest lattice +// point to ln X. Five small integers name a number: that is a codec, and this spec keeps it as one +// (`best_code` packs the five into one i64, `code_value` unpacks it). +// +// WHAT A FIT IS WORTH. A lattice dense enough fits anything. Measured on the Railway lab, 2026-10-18, +// 200 000 log-uniform random targets in [1e-3, 1e4]: +// C = 4: 2 889 points, P(err < 0.1%) = 22.5%, median err 0.315% +// C = 6: 11 601 points, P(err < 0.1%) = 62.3%, median err 0.072%, P(err < 0.01%) = 9.0% +// C = 8: 32 841 points, P(err < 0.1%) = 93.6%, median err 0.026% +// Those numbers follow a Poisson law, P(err < x) = 1 - exp(-2 x RHO), with RHO = 483 points per unit +// of ln at C = 6 (it gives 62% and 9.2%). `chance` is that law. The ten dimensionless targets in +// `REAL_*` below fit within 0.1% at C = 6 six times out of ten -- the rate a random number gets. +// The grid of the 2026-02 TZ (each exponent -20..20, no cap) has 25 431 849 points: it fits any +// number to about 1e-5, so a fit on it carries no information at all. +// +// WHAT IT DOES NOT CLAIM: no fit here is a derivation, and none is labelled EXACT or DERIVED. A fit +// is SUGGESTED only when its chance, multiplied by the number of targets tried, is under 5%; none of +// the ten passes. m_p/m_e = 6 pi^5 (18.8 ppm) is the best of them and is Lenz's 1951 observation +// (F. Lenz, Phys. Rev. 82, 554), not a new one; its chance is 1.8% per target, 17% over ten. +// Fits of quantities with units (m_H in GeV, H0 in km/s/Mpc, T_CMB in K, an age in Gyr) are not +// made: changing the unit moves the target and the lattice still fits it. +// alpha_phi = phi^-3 / 2 (the paper's eq. 11) is not in the grammar, n = 1/2 is not in 1..9; the +// nearest grammar point is 3 * 3^-2 * phi^2 * e^-2, 594 ppm away. +// Physics that needs a constant reads the measured value (CODATA 2022, PDG 2024): see +// specs/t27f/particle_field.t27, which never takes a constant from this lattice. +// phi^2 + 1/phi^2 = 3 | TRINITY + +module sacred_lattice { + + pub const PHI : f64 = 1.6180339887498948; + pub const LN3 : f64 = 1.0986122886681098; + pub const LNPI : f64 = 1.1447298858494002; + pub const LNPHI : f64 = 0.48121182505960347; + pub const RHO_C6 : f64 = 483.0; // effective points per unit of ln at C = 6 (lab, 200k targets) + pub const LN_LO : f64 = -6.907755278982137; // ln 1e-3 + pub const LN_HI : f64 = 9.210340371976184; // ln 1e4 + + pub const LABEL_EXACT : u8 = 0; // an identity: holds to the last bit (phi^2 + phi^-2 = 3) + pub const LABEL_MEASURED : u8 = 1; // CODATA / PDG value with its uncertainty + pub const LABEL_SUGGESTED : u8 = 2; // a lattice fit unlikely to be chance after all trials + pub const LABEL_UNKNOWN : u8 = 3; // a lattice fit at the chance rate: says nothing + + // Ten dimensionless targets (CODATA 2022, PDG 2024 central values). + pub const REAL_ALPHA_INV : f64 = 137.035999177; + pub const REAL_MP_ME : f64 = 1836.15267343; + pub const REAL_MMU_ME : f64 = 206.7682830; + pub const REAL_MTAU_ME : f64 = 3477.23; + pub const REAL_SIN2_THW : f64 = 0.23122; + pub const REAL_ALPHA_S : f64 = 0.1180; + pub const REAL_OMEGA_L : f64 = 0.6847; + pub const REAL_SIN2_TH12 : f64 = 0.307; + pub const REAL_SIN2_TH23 : f64 = 0.546; + pub const REAL_VUS : f64 = 0.2243; + pub const ALPHA_S_SIGMA : f64 = 0.0009; // PDG 2024 world average uncertainty + + fn iabs(x: i64) -> i64 { + if (x < 0) { + return 0 - x; + } + return x; + } + + // ln of the lattice point (n, k, m, p, q). + fn address_ln(n: i64, k: i64, m: i64, p: i64, q: i64) -> f64 { + return @log(n as f64) + (k as f64) * LN3 + (m as f64) * LNPI + (p as f64) * LNPHI + (q as f64); + } + + fn value(n: i64, k: i64, m: i64, p: i64, q: i64) -> f64 { + return @exp(address_ln(n, k, m, p, q)); + } + + // Relative distance between two positive numbers given by their logs. + fn rel_from_ln(a: f64, b: f64) -> f64 { + return @abs(@exp(a - b) - 1.0); + } + + // Number of grammar points at complexity C (n in 1..9). + fn lattice_size(c: i64) -> i64 { + var count : i64 = 0; + var k : i64 = 0 - c; + while (k <= c) { + var m : i64 = 0 - c; + while (m <= c) { + var p : i64 = 0 - c; + while (p <= c) { + var q : i64 = 0 - c; + while (q <= c) { + if (iabs(k) + iabs(m) + iabs(p) + iabs(q) <= c) { + count = count + 9; + } + q = q + 1; + } + p = p + 1; + } + m = m + 1; + } + k = k + 1; + } + return count; + } + + // Pack an address into one i64: n, then k, m, p, q each offset by 32 in base 64. + fn pack(n: i64, k: i64, m: i64, p: i64, q: i64) -> i64 { + return (((n * 64 + (k + 32)) * 64 + (m + 32)) * 64 + (p + 32)) * 64 + (q + 32); + } + + fn code_n(code: i64) -> i64 { + return @divTrunc(code, 16777216); + } + + fn code_digit(code: i64, place: i64) -> i64 { + var c : i64 = code; + var i : i64 = 0; + while (i < place) { + c = @divTrunc(c, 64); + i = i + 1; + } + return @mod(c, 64) - 32; + } + + fn code_value(code: i64) -> f64 { + return value(code_n(code), code_digit(code, 3), code_digit(code, 2), code_digit(code, 1), code_digit(code, 0)); + } + + // The nearest grammar point to x at complexity C, as a packed address. + fn best_code(x: f64, c: i64) -> i64 { + const t = @log(x); + var best : f64 = 1.0e9; + var code : i64 = 0; + var k : i64 = 0 - c; + while (k <= c) { + var m : i64 = 0 - c; + while (m <= c) { + var p : i64 = 0 - c; + while (p <= c) { + var q : i64 = 0 - c; + while (q <= c) { + if (iabs(k) + iabs(m) + iabs(p) + iabs(q) <= c) { + var n : i64 = 1; + while (n <= 9) { + const e = rel_from_ln(address_ln(n, k, m, p, q), t); + if (e < best) { + best = e; + code = pack(n, k, m, p, q); + } + n = n + 1; + } + } + q = q + 1; + } + p = p + 1; + } + m = m + 1; + } + k = k + 1; + } + return code; + } + + fn best_err(x: f64, c: i64) -> f64 { + return @abs(code_value(best_code(x, c)) / x - 1.0); + } + + // Chance that a random target lands this close to some lattice point (Poisson law, see header). + fn chance(err: f64, rho: f64) -> f64 { + return 1.0 - @exp(0.0 - 2.0 * err * rho); + } + + // A fit is SUGGESTED only when it would be rare even after `trials` targets were tried. + fn label(err: f64, rho: f64, trials: f64) -> u8 { + if (chance(err, rho) * trials < 0.05) { + return LABEL_SUGGESTED; + } + return LABEL_UNKNOWN; + } + + // Deterministic log-uniform targets in [1e-3, 1e4]: the null model, reproducible bit for bit. + fn lcg(seed: i64) -> i64 { + return @mod(seed * 1103515245 + 12345, 2147483648); + } + + fn null_hits(count: i64, seed0: i64, c: i64, tol: f64) -> i64 { + var seed : i64 = seed0; + var hits : i64 = 0; + var i : i64 = 0; + while (i < count) { + seed = lcg(seed); + const u = (seed as f64) / 2147483648.0; + const x = @exp(LN_LO + u * (LN_HI - LN_LO)); + if (best_err(x, c) < tol) { + hits = hits + 1; + } + i = i + 1; + } + return hits; + } + + fn real_hits(c: i64, tol: f64) -> i64 { + var h : i64 = 0; + if (best_err(REAL_ALPHA_INV, c) < tol) { h = h + 1; } + if (best_err(REAL_MP_ME, c) < tol) { h = h + 1; } + if (best_err(REAL_MMU_ME, c) < tol) { h = h + 1; } + if (best_err(REAL_MTAU_ME, c) < tol) { h = h + 1; } + if (best_err(REAL_SIN2_THW, c) < tol) { h = h + 1; } + if (best_err(REAL_ALPHA_S, c) < tol) { h = h + 1; } + if (best_err(REAL_OMEGA_L, c) < tol) { h = h + 1; } + if (best_err(REAL_SIN2_TH12, c) < tol) { h = h + 1; } + if (best_err(REAL_SIN2_TH23, c) < tol) { h = h + 1; } + if (best_err(REAL_VUS, c) < tol) { h = h + 1; } + return h; + } + + fn trinity(phi: f64) -> f64 { + return phi * phi + 1.0 / (phi * phi); + } + + // Points of the uncapped grid with every exponent in -span..span. + fn grid_size(span: i64) -> i64 { + const w = 2 * span + 1; + return 9 * w * w * w * w; + } + + test "the trinity identity is exact: phi^2 + phi^-2 = 3" { + assert(@abs(trinity(PHI) - 3.0) < 1.0e-15); + assert(@abs(trinity((1.0 + @sqrt(5.0)) / 2.0) - 3.0) < 1.0e-15); + } + + test "the lattice sizes match the lab count" { + assert(lattice_size(4) == 2889); + assert(lattice_size(6) == 11601); + assert(lattice_size(8) == 32841); + } + + test "the TZ grid -20..20 has 25 431 849 points" { + assert(grid_size(20) == 25431849); + } + + test "an address round-trips through its packed code" { + const code = pack(4, 2, 0 - 1, 1, 2); + assert(code_n(code) == 4); + assert(code_digit(code, 3) == 2); + assert(code_digit(code, 2) == 0 - 1); + assert(code_digit(code, 1) == 1); + assert(code_digit(code, 0) == 2); + } + + test "1/alpha sits at 4 * 3^2 * pi^-1 * phi * e^2, 243 ppm away" { + const code = best_code(REAL_ALPHA_INV, 6); + assert(code == pack(4, 2, 0 - 1, 1, 2)); + const e = best_err(REAL_ALPHA_INV, 6); + assert(e > 0.000242 && e < 0.000244); + assert(label(e, RHO_C6, 10.0) == LABEL_UNKNOWN); + } + + test "m_p/m_e = 6 pi^5 is the best fit and is still not SUGGESTED over ten targets" { + assert(best_code(REAL_MP_ME, 6) == pack(6, 0, 5, 0, 0)); + const e = best_err(REAL_MP_ME, 6); + assert(e > 0.0000187 && e < 0.0000189); + assert(chance(e, RHO_C6) > 0.017 && chance(e, RHO_C6) < 0.019); + assert(label(e, RHO_C6, 10.0) == LABEL_UNKNOWN); + assert(label(e, RHO_C6, 1.0) == LABEL_SUGGESTED); + } + + test "the Poisson law reproduces the lab null rates at C = 6" { + assert(@abs(chance(0.001, RHO_C6) - 0.623) < 0.01); + assert(@abs(chance(0.0001, RHO_C6) - 0.090) < 0.005); + } + + test "random numbers fit within 0.1% as often as the real constants do" { + assert(null_hits(64, 27, 6, 0.001) == 42); + assert(real_hits(6, 0.001) == 6); + // 42/64 = 66% for noise, 6/10 = 60% for physics: no excess over chance. + } + + test "alpha_s is measured too coarsely to test any lattice point" { + assert(chance(ALPHA_S_SIGMA / REAL_ALPHA_S, RHO_C6) > 0.99); + } + + test "phi^-3 / 2 is not a grammar point; the nearest one is 594 ppm off" { + const x = (@sqrt(5.0) - 2.0) / 2.0; + assert(best_code(x, 6) == pack(3, 0 - 2, 0, 2, 0 - 2)); + assert(@abs(best_err(x, 6) - 0.000594) < 0.000001); + } + + invariant "no lattice fit is ever labelled exact or measured" { + assert(label(0.0, RHO_C6, 1.0) == LABEL_SUGGESTED); + assert(label(1.0, RHO_C6, 1.0) == LABEL_UNKNOWN); + } +} diff --git a/specs/t27f/particle_field.t27 b/specs/t27f/particle_field.t27 new file mode 100644 index 0000000000..939107b74d --- /dev/null +++ b/specs/t27f/particle_field.t27 @@ -0,0 +1,512 @@ +// SPDX-License-Identifier: Apache-2.0 +// specs/t27f/particle_field.t27 -- T27F Live drawn as physics: talking to an agent becomes a field of +// particles, and every motion on the screen is a step of a real integrator, not an easing curve. +// +// WHERE IT SITS: docs/nona-02-organism/t27f/genui-research.md section 5 has the pipeline +// agent events -> work graph -> scene -> display list -> paint. +// This spec is the physical layer between the work graph and the scene. Nodes of the work graph +// (the user, the agent, a tool) are massive particles; an edge is a damped spring; every pair of +// nodes pushes apart like charges, so the layout finds itself (P. Eades, 1984; T. Fruchterman and +// E. Reingold, Softw. Pract. Exper. 21, 1991). An event emits packets that fly along an edge and +// give the receiving node a kick, so the picture moves only because something happened. +// +// THE EVENT VOCABULARY, each with a ternary charge (+1 flows toward the agent, -1 flows out of it, +// 0 stays inside it) and a hue: +// USER_MSG +1 gold user -> agent, 1 + chars/8 packets (at most 16) +// THINK 0 violet the agent heats up: a seeded random jitter, no direction +// TOOL_CALL -1 green agent -> tool +// TOOL_RESULT +1 green tool -> agent +// AGENT_MSG -1 cyan agent -> user +// ERROR -1 red a radial burst around the agent +// DONE 0 white damping rises and the field comes to rest +// A finished exchange is neutral: what came in went out (test "the demo log is neutral"). +// +// THE NUMERICS: a fixed step DT = 1/60 s, one step per displayed frame, so a frame index names a +// state and replaying the same log gives the same picture bit for bit. The node field uses +// semi-implicit Euler, the one-body checks use velocity Verlet (L. Verlet, Phys. Rev. 159, 98, 1967); +// both are symplectic, so energy does not creep (E. Hairer, C. Lubich, G. Wanner, "Geometric +// Numerical Integration", 2nd ed., 2006). The jitter comes from a fixed LCG, never from a clock. +// Reduced motion (the OS setting) turns kicks and jitter off; the packets still travel, so nothing +// is hidden, only nothing shakes. +// +// THE MESH: each particle is a camera-facing quad, 4 vertices and 2 triangles, so 16 particles fill +// one meshlet of 64 vertices and 32 primitives, inside the 64 / 126 limits NVIDIA published for +// Turing mesh shaders (2018). A WebGPU or canvas host consumes the same counts. +// +// WHAT IT DOES NOT CLAIM: the field is a picture of the conversation, not a model of the agent's +// mind; temperature is a display parameter, not a sampling temperature. No constant here comes +// from the V = n * 3^k * pi^m * phi^p * e^q lattice (specs/physics/sacred_lattice.t27 shows those +// fits are at the chance rate). The physical constants in the Bohr check are CODATA 2022. +// Claim status: every number on the screen is a test below; the field constants are design +// choices, labelled as such. +// phi^2 + 1/phi^2 = 3 | TRINITY + +module particle_field { + + pub const KIND : str = "widget-logic"; + pub const ID : str = "particle-field"; + pub const VERSION : u8 = 1; + + pub const DT : f64 = 0.016666666666666666; // 1/60 s, one step per frame + pub const PI : f64 = 3.141592653589793; + + // Event kinds. + pub const USER_MSG : i64 = 0; + pub const THINK : i64 = 1; + pub const TOOL_CALL : i64 = 2; + pub const TOOL_RESULT : i64 = 3; + pub const AGENT_MSG : i64 = 4; + pub const ERROR : i64 = 5; + pub const DONE : i64 = 6; + + // Nodes. + pub const N_USER : usize = 0; + pub const N_AGENT : usize = 1; + pub const N_TOOL : usize = 2; + + // Field constants: design choices, not measurements. + pub const SPRING_K : f64 = 4.0; + pub const REST : f64 = 1.0; + pub const REPULSE : f64 = 0.2; + pub const SOFT : f64 = 0.01; + pub const DAMP_LIVE : f64 = 0.8; + pub const DAMP_DONE : f64 = 6.0; + pub const KICK : f64 = 0.6; + pub const THINK_TEMP : f64 = 2.0; + pub const LOG_START : i64 = 120; // the first event; frames before it settle the layout + + // Mesh budget (NVIDIA Turing mesh shader limits). + pub const VERTS_PER_PARTICLE : i64 = 4; + pub const TRIS_PER_PARTICLE : i64 = 2; + pub const MESHLET_VERTS : i64 = 64; + pub const MESHLET_PRIMS : i64 = 126; + + // CODATA 2022. + pub const E_CHARGE : f64 = 1.602176634e-19; + pub const EPS0 : f64 = 8.8541878188e-12; + pub const HBAR : f64 = 1.0545718176461565e-34; // h / (2 pi), h = 6.62607015e-34 exact; 1.054571817e-34 is 6e-10 short + pub const C_LIGHT : f64 = 299792458.0; + pub const M_E : f64 = 9.1093837139e-31; + pub const ALPHA : f64 = 0.0072973525643; + pub const A0 : f64 = 5.29177210544e-11; + + // ---- the event vocabulary ---- + + fn charge(kind: i64) -> i64 { + if (kind == USER_MSG || kind == TOOL_RESULT) { + return 1; + } + if (kind == TOOL_CALL || kind == AGENT_MSG || kind == ERROR) { + return 0 - 1; + } + return 0; + } + + // Hue in degrees. + fn hue(kind: i64) -> f64 { + if (kind == USER_MSG) { return 45.0; } + if (kind == THINK) { return 275.0; } + if (kind == TOOL_CALL || kind == TOOL_RESULT) { return 140.0; } + if (kind == AGENT_MSG) { return 190.0; } + if (kind == ERROR) { return 0.0; } + return 0.0 - 1.0; // DONE: white, no hue + } + + // Packets a message of `chars` characters emits. + fn burst(chars: i64) -> i64 { + const n = 1 + @divTrunc(chars, 8); + if (n > 16) { + return 16; + } + return n; + } + + fn sender(kind: i64) -> usize { + if (kind == USER_MSG) { return N_USER; } + if (kind == TOOL_RESULT) { return N_TOOL; } + return N_AGENT; + } + + fn receiver(kind: i64) -> usize { + if (kind == USER_MSG || kind == TOOL_RESULT) { return N_AGENT; } + if (kind == TOOL_CALL) { return N_TOOL; } + if (kind == AGENT_MSG) { return N_USER; } + return N_AGENT; + } + + // The demo conversation: one question, one thought, one tool round trip, one answer. + pub const LOG_LEN : i64 = 6; + + fn log_kind(i: i64) -> i64 { + if (i == 0) { return USER_MSG; } + if (i == 1) { return THINK; } + if (i == 2) { return TOOL_CALL; } + if (i == 3) { return TOOL_RESULT; } + if (i == 4) { return AGENT_MSG; } + return DONE; + } + + fn log_frame(i: i64) -> i64 { + if (i == 0) { return 120; } + if (i == 1) { return 150; } + if (i == 2) { return 210; } + if (i == 3) { return 260; } + if (i == 4) { return 300; } + return 330; + } + + fn log_charge() -> i64 { + var s : i64 = 0; + var i : i64 = 0; + while (i < LOG_LEN) { + s = s + charge(log_kind(i)); + i = i + 1; + } + return s; + } + + // The owner's trit triangle, read row by row: +1 / -1 +1 / +1 0 +1 / -1 +1 +1 -1. + fn triangle_sum() -> i64 { + const t : [10]i64 = [1, -1, 1, 1, 0, 1, -1, 1, 1, -1]; + var s : i64 = 0; + var i : usize = 0; + while (i < 10) { + s = s + t[i]; + i = i + 1; + } + return s; + } + + // ---- packets ---- + + // Distance a packet launched at frame `start` has covered by frame `frame`. + fn packet_s(frame: i64, start: i64, speed: f64, dist: f64) -> f64 { + if (frame <= start) { + return 0.0; + } + const s = ((frame - start) as f64) * speed * DT; + if (s > dist) { + return dist; + } + return s; + } + + fn arrival_frames(dist: f64, speed: f64) -> i64 { + return @as(i64, @intFromFloat(@ceil(dist / (speed * DT)))); + } + + fn arrival_by_steps(dist: f64, speed: f64) -> i64 { + var s : f64 = 0.0; + var n : i64 = 0; + while (s < dist) { + s = s + speed * DT; + n = n + 1; + } + return n; + } + + // ---- mesh ---- + + fn mesh_verts(n: i64) -> i64 { + return n * VERTS_PER_PARTICLE; + } + + fn mesh_tris(n: i64) -> i64 { + return n * TRIS_PER_PARTICLE; + } + + fn meshlets(n: i64) -> i64 { + return @divTrunc(mesh_verts(n) + MESHLET_VERTS - 1, MESHLET_VERTS); + } + + // ---- one-body checks of the integrator ---- + + // Period of x'' = -(k/m) x with velocity Verlet, measured between two downward zero crossings. + fn verlet_period(k: f64, m: f64) -> f64 { + const w2 = k / m; + var x : f64 = 1.0; + var v : f64 = 0.0; + var a : f64 = 0.0 - w2 * x; + var t : f64 = 0.0; + var first : f64 = 0.0 - 1.0; + var steps : i64 = 0; + while (steps < 1000000) { + const xn = x + v * DT + 0.5 * a * DT * DT; + const an = 0.0 - w2 * xn; + v = v + 0.5 * (a + an) * DT; + if (x > 0.0 && xn <= 0.0) { + const tc = t + DT * x / (x - xn); + if (first >= 0.0) { + return tc - first; + } + first = tc; + } + x = xn; + a = an; + t = t + DT; + steps = steps + 1; + } + return 0.0; + } + + // Relative energy drift of a circular Kepler orbit (GM = 1, r = 1, v = 1) after `orbits` turns. + // verlet = true: velocity Verlet; false: explicit Euler. + fn orbit_drift(orbits: i64, verlet: bool) -> f64 { + var x : f64 = 1.0; + var y : f64 = 0.0; + var vx : f64 = 0.0; + var vy : f64 = 1.0; + const e0 = 0.0 - 0.5; + const n = @as(i64, @intFromFloat(@round(2.0 * PI * (orbits as f64) / DT))); + var i : i64 = 0; + while (i < n) { + const r3 = @sqrt(x * x + y * y) * (x * x + y * y); + const ax = 0.0 - x / r3; + const ay = 0.0 - y / r3; + if (verlet) { + const xn = x + vx * DT + 0.5 * ax * DT * DT; + const yn = y + vy * DT + 0.5 * ay * DT * DT; + const rn3 = @sqrt(xn * xn + yn * yn) * (xn * xn + yn * yn); + vx = vx + 0.5 * (ax - xn / rn3) * DT; + vy = vy + 0.5 * (ay - yn / rn3) * DT; + x = xn; + y = yn; + } else { + x = x + vx * DT; + y = y + vy * DT; + vx = vx + ax * DT; + vy = vy + ay * DT; + } + i = i + 1; + } + const e = 0.5 * (vx * vx + vy * vy) - 1.0 / @sqrt(x * x + y * y); + return @abs((e - e0) / e0); + } + + // Two nodes on one edge: where spring and repulsion balance, k (r - REST) = REPULSE / r^2. + fn edge_rest_newton() -> f64 { + var r : f64 = REST; + var i : i64 = 0; + while (i < 50) { + const g = SPRING_K * (r - REST) - REPULSE / (r * r); + const dg = SPRING_K + 2.0 * REPULSE / (r * r * r); + r = r - g / dg; + i = i + 1; + } + return r; + } + + // The same balance found by the integrator (relative coordinate, reduced mass 1/2). + fn edge_rest_sim(frames: i64) -> f64 { + var r : f64 = 0.5; + var v : f64 = 0.0; + var f : i64 = 0; + while (f < frames) { + const force = 0.0 - SPRING_K * (r - REST) + REPULSE / (r * r) - DAMP_LIVE * v; + v = v + 2.0 * force * DT; + r = r + v * DT; + f = f + 1; + } + return r; + } + + // The speed of the electron on the first Bohr orbit, and that orbit's radius. + fn bohr_speed(e: f64, eps0: f64, hbar: f64) -> f64 { + return e * e / (4.0 * PI * eps0 * hbar); + } + + fn bohr_radius(e: f64, eps0: f64, hbar: f64, me: f64) -> f64 { + return 4.0 * PI * eps0 * hbar * hbar / (me * e * e); + } + + // ---- the field ---- + + fn lcg(seed: i64) -> i64 { + return @mod(seed * 1103515245 + 12345, 2147483648); + } + + fn linked(i: usize, j: usize) -> bool { + return (i == N_AGENT || j == N_AGENT) && i != j; + } + + // Run the demo log for `frames` frames. Returns the kinetic energy at the last frame, or, with + // want_peak, the largest kinetic energy seen after LOG_START. + fn field(frames: i64, reduced: bool, want_peak: bool) -> f64 { + var px : [3]f64 = [-1.0, 0.0, 1.0]; + var py : [3]f64 = [0.0, 0.3, 0.0]; + var vx : [3]f64 = [0.0, 0.0, 0.0]; + var vy : [3]f64 = [0.0, 0.0, 0.0]; + var seed : i64 = 27; + var temp : f64 = 0.0; + var damp : f64 = DAMP_LIVE; + var peak : f64 = 0.0; + var ke : f64 = 0.0; + var f : i64 = 0; + while (f < frames) { + var e : i64 = 0; + while (e < LOG_LEN) { + if (log_frame(e) == f) { + const kind = log_kind(e); + if (kind == THINK) { + temp = THINK_TEMP; + } + if (kind == DONE) { + temp = 0.0; + damp = DAMP_DONE; + } + if (charge(kind) != 0 && reduced == false) { + const s = sender(kind); + const r = receiver(kind); + const dx = px[r] - px[s]; + const dy = py[r] - py[s]; + const d = @sqrt(dx * dx + dy * dy + SOFT); + vx[r] = vx[r] + KICK * dx / d; + vy[r] = vy[r] + KICK * dy / d; + } + } + e = e + 1; + } + var i : usize = 0; + while (i < 3) { + var ax : f64 = 0.0 - damp * vx[i]; + var ay : f64 = 0.0 - damp * vy[i]; + var j : usize = 0; + while (j < 3) { + if (j != i) { + const dx = px[i] - px[j]; + const dy = py[i] - py[j]; + const r2 = dx * dx + dy * dy + SOFT; + const r = @sqrt(r2); + var fm : f64 = REPULSE / r2; + if (linked(i, j)) { + fm = fm - SPRING_K * (r - REST); + } + ax = ax + fm * dx / r; + ay = ay + fm * dy / r; + } + j = j + 1; + } + if (temp > 0.0 && reduced == false) { + seed = lcg(seed); + ax = ax + temp * ((seed as f64) / 2147483648.0 - 0.5); + seed = lcg(seed); + ay = ay + temp * ((seed as f64) / 2147483648.0 - 0.5); + } + vx[i] = vx[i] + ax * DT; + vy[i] = vy[i] + ay * DT; + i = i + 1; + } + ke = 0.0; + var k : usize = 0; + while (k < 3) { + px[k] = px[k] + vx[k] * DT; + py[k] = py[k] + vy[k] * DT; + ke = ke + 0.5 * (vx[k] * vx[k] + vy[k] * vy[k]); + k = k + 1; + } + if (f >= LOG_START && ke > peak) { + peak = ke; + } + f = f + 1; + } + if (want_peak) { + return peak; + } + return ke; + } + + // ---- tests ---- + + test "the event vocabulary: charges and hues" { + assert(charge(USER_MSG) == 1); + assert(charge(THINK) == 0); + assert(charge(TOOL_CALL) == 0 - 1); + assert(charge(TOOL_RESULT) == 1); + assert(charge(AGENT_MSG) == 0 - 1); + assert(charge(ERROR) == 0 - 1); + assert(charge(DONE) == 0); + assert(hue(USER_MSG) == 45.0); + assert(hue(AGENT_MSG) == 190.0); + assert(hue(ERROR) == 0.0); + assert(hue(DONE) < 0.0); + } + + test "a message bursts into 1 + chars/8 packets, at most 16" { + assert(burst(0) == 1); + assert(burst(48) == 7); + assert(burst(1000) == 16); + } + + test "the demo log is neutral: what came in went out" { + assert(log_charge() == 0); + assert(charge(TOOL_CALL) + charge(TOOL_RESULT) == 0); + } + + test "the owner's trit triangle sums to +3" { + assert(triangle_sum() == 3); + } + + test "a packet arrives in ceil(dist / (speed * DT)) frames and stops at the node" { + assert(arrival_frames(1.0, 0.7) == 86); + assert(arrival_by_steps(1.0, 0.7) == 86); + assert(packet_s(120, 120, 0.7, 1.0) == 0.0); + assert(packet_s(400, 120, 0.7, 1.0) == 1.0); + } + + test "16 particles fill one meshlet within the Turing limits" { + assert(mesh_verts(16) == MESHLET_VERTS); + assert(mesh_tris(16) <= MESHLET_PRIMS); + assert(meshlets(16) == 1); + assert(meshlets(17) == 2); + assert(meshlets(1000) == 63); + } + + test "velocity Verlet keeps the harmonic period 2 pi sqrt(m / k)" { + const t1 = verlet_period(1.0, 1.0); + assert(@abs(t1 / (2.0 * PI) - 1.0) < 0.0001); + const t2 = verlet_period(4.0, 1.0); + assert(@abs(t2 / PI - 1.0) < 0.0001); + } + + test "Verlet holds the orbit's energy; explicit Euler does not" { + assert(orbit_drift(10, true) < 0.00001); + assert(orbit_drift(10, false) > 0.01); + } + + test "the Bohr speed is alpha c (CODATA 2022)" { + assert(@abs(bohr_speed(E_CHARGE, EPS0, HBAR) / C_LIGHT / ALPHA - 1.0) < 0.000000001); + assert(@abs(bohr_radius(E_CHARGE, EPS0, HBAR, M_E) / A0 - 1.0) < 0.000000001); + // the truncated hbar 1.054571817e-34 misses a0 by 1.2e-9, outside this tolerance + assert(@abs(bohr_radius(E_CHARGE, EPS0, 1.054571817e-34, M_E) / A0 - 1.0) > 0.000000001); + } + + test "an edge settles where spring and repulsion balance" { + const r = edge_rest_newton(); + assert(@abs(SPRING_K * (r - REST) - REPULSE / (r * r)) < 0.000000000001); + assert(@abs(edge_rest_sim(1200) - r) < 0.000001); + } + + test "replaying the same log gives the same field, bit for bit" { + assert(field(360, false, false) == field(360, false, false)); + assert(field(360, false, true) == field(360, false, true)); + } + + // The bend of user-agent-tool is held only by the weak user-tool repulsion, so the last motion + // to die is slow: measured, KE / peak is 3.2e-5 at frame 600 and 5.4e-6 at frame 2400. + test "after DONE the field comes to rest within 4.5 s" { + const peak = field(600, false, true); + assert(peak > 0.0); + assert(field(600, false, false) < peak * 0.0001); + } + + test "reduced motion keeps the field calm" { + assert(field(600, true, true) < field(600, false, true) * 0.1); + } + + test "VERSION says what the header says" { + assert(VERSION == 1); + assert(LOG_LEN == 6); + assert(log_frame(0) == LOG_START); + } +} +// phi^2 + 1/phi^2 = 3 | TRINITY