Skip to main content

graph_explorer_layout/
force_sim.rs

1//! A persistent force simulation: node positions and velocities that survive
2//! across frames, so a graph change nudges the layout instead of recomputing it.
3//!
4//! Parameter names mirror `d3-force` deliberately — anyone who has tuned a
5//! force graph already knows what `alpha_decay` and `velocity_decay` do.
6//!
7//! Pure: no wasm, no GPU, and no time source. Ticks are driven by the caller,
8//! never self-scheduled, which is what makes the whole thing host-testable.
9
10use crate::QuadTree;
11use graph_explorer_core::{Graph, NodeId};
12use serde::Deserialize;
13use std::collections::HashMap;
14
15#[derive(Debug, Clone, Deserialize)]
16#[serde(deny_unknown_fields)]
17pub struct SimulationSpec {
18    /// When false the simulation never ANIMATES — but a Force-source
19    /// load/relayout still runs a synchronous settle so disabled hosts get a
20    /// finished static layout (see `settle_forced`).
21    #[serde(default = "yes")]
22    pub enabled: bool,
23    /// Per-tick multiplicative decay of `alpha`. The default reaches
24    /// `alpha_min` in ~300 ticks: `1 - 0.001^(1/300)`.
25    #[serde(default = "d_alpha_decay")]
26    pub alpha_decay: f32,
27    /// Below this the sim is considered settled and stops ticking, so a graph
28    /// at rest costs nothing.
29    #[serde(default = "d_alpha_min")]
30    pub alpha_min: f32,
31    /// Friction: velocity is multiplied by `1 - velocity_decay` each tick.
32    #[serde(default = "d_velocity_decay")]
33    pub velocity_decay: f32,
34    /// `alpha` restored on a graph change. Lower than 1 so a reveal nudges
35    /// rather than re-scrambling.
36    #[serde(default = "d_reheat_alpha")]
37    pub reheat_alpha: f32,
38    /// Barnes-Hut accuracy. 0 is exact and `O(n²)`; larger approximates more.
39    #[serde(default = "d_theta")]
40    pub theta: f32,
41    /// Ideal edge length, in world units.
42    #[serde(default = "d_link_distance")]
43    pub link_distance: f32,
44    /// Many-body strength. Negative repels.
45    #[serde(default = "d_charge")]
46    pub charge_strength: f32,
47    /// Pull toward the centroid. Without it, disconnected components drift
48    /// apart forever — the old (now-deleted) one-shot layout had no such force.
49    #[serde(default = "d_center")]
50    pub center_strength: f32,
51    /// What happens to a node when the user lets go of it.
52    #[serde(default)]
53    pub drag_release: DragRelease,
54}
55
56/// What becomes of a dragged node once the pointer is released.
57///
58/// This is interaction policy rather than physics, but it lives here because
59/// what it selects between is whether the node stays in the simulation's pin
60/// set — and because a host tuning drag feel is already tuning the rest of the
61/// simulation spec.
62#[derive(Debug, Clone, Copy, PartialEq, Eq, Default, Deserialize)]
63#[serde(rename_all = "kebab-case")]
64pub enum DragRelease {
65    /// Release the node. Physics pulls it back toward equilibrium with
66    /// everything else and the neighbourhood re-settles around it.
67    ///
68    /// The default: it is the only one of the three that cannot accumulate
69    /// state, so a host that never thinks about pinning never acquires a graph
70    /// full of immovable nodes.
71    #[default]
72    Free,
73    /// Hold the node where it was dropped, indefinitely. Lets a user arrange a
74    /// layout by hand. Requires the host to offer some way to un-stick a node,
75    /// or pins accumulate with no gesture to release them.
76    Pin,
77    /// Hold the node where it was dropped, until the next reheat frees it — so
78    /// the placement survives the settle it caused, and a later reveal or
79    /// relayout starts clean.
80    PinUntilReheat,
81}
82
83fn yes() -> bool { true }
84fn d_alpha_decay() -> f32 { 0.0228 }
85fn d_alpha_min() -> f32 { 0.001 }
86fn d_velocity_decay() -> f32 { 0.4 }
87fn d_reheat_alpha() -> f32 { 0.3 }
88fn d_theta() -> f32 { 0.9 }
89fn d_link_distance() -> f32 { 120.0 }
90fn d_charge() -> f32 { -800.0 }
91fn d_center() -> f32 { 0.02 }
92
93impl Default for SimulationSpec {
94    fn default() -> Self {
95        Self {
96            enabled: yes(),
97            alpha_decay: d_alpha_decay(),
98            alpha_min: d_alpha_min(),
99            velocity_decay: d_velocity_decay(),
100            reheat_alpha: d_reheat_alpha(),
101            theta: d_theta(),
102            link_distance: d_link_distance(),
103            charge_strength: d_charge(),
104            center_strength: d_center(),
105            drag_release: DragRelease::default(),
106        }
107    }
108}
109
110pub struct ForceSimulation {
111    spec: SimulationSpec,
112    index: HashMap<NodeId, usize>,
113    ids: Vec<NodeId>,
114    pos: Vec<[f32; 2]>,
115    vel: Vec<[f32; 2]>,
116    links: Vec<(usize, usize)>,
117    /// Nodes held in place, each at the point it is held at. Captured when the
118    /// pin is set, so a node returns to exactly that point rather than to
119    /// wherever the last integration happened to leave it.
120    ///
121    /// Keyed by id rather than index, because `sync` re-sorts the id list and
122    /// therefore reshuffles every index whenever a node arrives or leaves. An
123    /// index-keyed pin has to be re-resolved through the id map on every sync
124    /// or it silently holds the wrong node; keying by id removes that step
125    /// rather than repeating it per pin. The cost is one index lookup per
126    /// pinned node per tick, against a tick that already builds a Barnes-Hut
127    /// tree over every node.
128    ///
129    /// A set, not a slot: navigation holds the node the user is standing on,
130    /// and dragging holds the node under the cursor. Those are independent and
131    /// concurrent.
132    pinned: HashMap<NodeId, [f32; 2]>,
133    alpha: f32,
134    /// Monotonic count of nodes that have entered without a seed. Persisting it
135    /// on the struct — rather than restarting per `sync` — is what stops two
136    /// arrivals in *different* syncs from being handed the same spiral point
137    /// and stacking, which is the exact failure the spiral exists to prevent.
138    spawn_seq: usize,
139}
140
141/// Upper bound for `theta`. Beyond roughly 1.25 — and as low as 0.75 for a
142/// tight cluster with a distant outlier, which is precisely the shape an
143/// expanded aggregate makes — a cell that *contains* body `i` can pass the
144/// approximation test, and `i` then repels itself through its own cell's mass.
145/// Error stays within 21% at 1.0, so that is where the config boundary sits.
146/// Clamping here rather than in `walk` keeps the check off the hot path.
147const THETA_MAX: f32 = 1.0;
148
149/// Floor for `alpha_min`. `is_active` is `alpha >= alpha_min`, and `alpha`
150/// decays multiplicatively toward zero without ever reaching it — so at
151/// `alpha_min = 0` the predicate is permanently true and the simulation never
152/// settles. That is not a slow layout, it is three broken guarantees: the
153/// host's render loop runs a Barnes-Hut tick every frame at rest, the
154/// settle-detecting `!tick()` edge never fires (so a caller that rebuilds its
155/// scene there never does, and its target goes stale indefinitely), and
156/// `settle(max)` burns its entire budget and still reports active, which is
157/// what a reduced-motion caller uses to mean "snap, do not animate".
158///
159/// `1e-4` is the smallest value that still terminates promptly: from a full
160/// reheat at the default decay it is reached in ~400 ticks.
161const ALPHA_MIN_FLOOR: f32 = 1e-4;
162
163/// Ceiling for `alpha_min`, guarding the opposite failure. `alpha_min` above
164/// the alpha a reheat restores means `is_active` is false the moment it is
165/// set, so every reheat is silently a no-op and the simulation appears dead.
166/// `0.1` leaves headroom under the default `reheat_alpha` of 0.3. A host that
167/// also lowers `reheat_alpha` below its `alpha_min` can still arrange that
168/// no-op, but it fails safe — no motion, no CPU — rather than wedging.
169const ALPHA_MIN_CEIL: f32 = 0.1;
170
171/// Floor for `alpha_decay`. The schedule is `alpha -= alpha * alpha_decay`, so
172/// at zero `alpha` never decreases: the same never-settles failure as
173/// `alpha_min = 0`, reached from the other side.
174///
175/// The value is derived, not chosen. `settle(600)` is the reduced-motion
176/// contract — settle synchronously, then snap, with no motion afterwards — and
177/// it only holds if 600 ticks is actually enough. Worst case is a full reheat
178/// (`alpha = 1.0`, the ceiling `reheat_alpha` is clamped to) decaying to the
179/// smallest legal `alpha_min`, `ALPHA_MIN_FLOOR`, which takes
180/// `ln(1e-4) / ln(1 - alpha_decay)` ticks: 571 at `0.016`, against a budget of
181/// 600. Anything slower silently converts reduced motion into permanent
182/// motion. The price is that a host cannot ask for a settle longer than about
183/// 9.5 seconds at 60fps, which is a fair ceiling on how long a graph may wobble.
184const ALPHA_DECAY_MIN: f32 = 0.016;
185
186/// Ceiling for `alpha_decay`. At 1.0 `alpha` drops to zero in a single tick,
187/// which is already the fastest meaningful schedule; above it `alpha` goes
188/// negative, which is not a faster settle, just a nonsensical state.
189const ALPHA_DECAY_MAX: f32 = 1.0;
190
191/// Floor for `velocity_decay`, which enters the integrator as
192/// `friction = 1 - velocity_decay`. At zero there is no damping at all: a node
193/// keeps whatever velocity it accumulated and coasts in a straight line until
194/// `alpha` decays out, freezing it wherever it happened to be rather than at a
195/// force fixpoint. Below zero `friction > 1` and this explicit-Euler step
196/// *amplifies* velocity every tick — exponential divergence, NaN within the
197/// usual few hundred ticks. `0.01` guarantees energy actually leaves the
198/// system.
199const VELOCITY_DECAY_MIN: f32 = 0.01;
200
201/// Ceiling for `velocity_decay`. 1.0 is full damping — velocity is discarded
202/// each tick and the step degenerates to plain gradient descent, which is
203/// stable and occasionally what a host wants. Above it `friction` is negative
204/// and velocity flips sign every tick, so the layout oscillates instead of
205/// converging (and diverges outright past 2.0).
206const VELOCITY_DECAY_MAX: f32 = 1.0;
207
208/// Radius scale for an unseeded spawn, in world units.
209const SPAWN_RADIUS: f32 = 30.0;
210
211/// How many spiral rings a *late* arrival may use before the radius wraps.
212/// `SPAWN_RADIUS * sqrt(11.5) ≈ 102`, so an arrival lands within about one
213/// link distance of the centroid instead of at `30·√(global rank)` — which for
214/// a 500-node graph was 671 units, 5.6 link distances, and nowhere near the
215/// "near the centroid" the doc contract promises. The angle keeps advancing
216/// off the unwrapped sequence, so wrapping the radius never repeats a point.
217const ARRIVAL_RINGS: usize = 12;
218
219impl ForceSimulation {
220    pub fn new(spec: SimulationSpec) -> Self {
221        Self {
222            spec: Self::clamped(spec),
223            index: HashMap::new(),
224            ids: Vec::new(),
225            pos: Vec::new(),
226            vel: Vec::new(),
227            links: Vec::new(),
228            pinned: HashMap::new(),
229            alpha: 0.0,
230            spawn_seq: 0,
231        }
232    }
233
234    pub fn set_spec(&mut self, spec: SimulationSpec) { self.spec = Self::clamped(spec); }
235    pub fn spec(&self) -> &SimulationSpec { &self.spec }
236
237    /// Bring a host-supplied spec inside the range the forces can actually
238    /// integrate, and the alpha schedule inside the range that can actually
239    /// terminate. Done once, where the spec is adopted, so no per-tick branch
240    /// pays for it — and here rather than at the call site because this is the
241    /// config boundary: `new` and `set_spec` are the only two ways a spec gets
242    /// in, and `set_spec` is reachable from a host in a single call.
243    ///
244    /// The alpha-schedule bounds are not taste. `is_active` is the sole thing
245    /// standing between a settled graph and a render loop that ticks forever,
246    /// and it is one comparison against two of these values; see the constants
247    /// for what each unbounded end breaks.
248    fn clamped(mut spec: SimulationSpec) -> SimulationSpec {
249        spec.theta = spec.theta.clamp(0.0, THETA_MAX);
250        spec.alpha_min = spec.alpha_min.clamp(ALPHA_MIN_FLOOR, ALPHA_MIN_CEIL);
251        spec.alpha_decay = spec.alpha_decay.clamp(ALPHA_DECAY_MIN, ALPHA_DECAY_MAX);
252        spec.velocity_decay =
253            spec.velocity_decay.clamp(VELOCITY_DECAY_MIN, VELOCITY_DECAY_MAX);
254        // Nothing above 1 is reachable by `reheat_full`, and a negative alpha is
255        // just an immediately-settled one wearing a disguise. Bounding it keeps
256        // `alpha` inside `[0, 1]` for every code path, which is what makes the
257        // tick-count arithmetic in the constants above true.
258        spec.reheat_alpha = spec.reheat_alpha.clamp(0.0, 1.0);
259        spec
260    }
261
262    /// Adopt the graph's topology, preserving position and velocity for every
263    /// id already present. New ids take `seed` when given — a revealed ghost
264    /// enters exactly where it was drawn, so physics continues from the
265    /// preview rather than teleporting. Does not reheat; callers decide that.
266    pub fn sync(&mut self, graph: &Graph, seed: &HashMap<NodeId, [f32; 2]>) {
267        let mut ids: Vec<NodeId> = graph.nodes().map(|n| n.id.clone()).collect();
268        ids.sort(); // deterministic index assignment, so layouts reproduce
269
270        let mut pos = Vec::with_capacity(ids.len());
271        let mut vel = Vec::with_capacity(ids.len());
272        let mut index = HashMap::with_capacity(ids.len());
273        let centroid = self.centroid();
274        // A cold start is the one case where spiral rank *is* arrival count, so
275        // `30·√N` is the right spread and is kept exactly. Every later sync is
276        // an arrival into an existing layout, where that radius would be a
277        // function of where the id happens to sort rather than of anything
278        // physical, so those wrap instead.
279        let cold_start = self.index.is_empty();
280        let mut spawned_here = 0usize;
281
282        for (i, id) in ids.iter().enumerate() {
283            match self.index.get(id) {
284                Some(&old) => { pos.push(self.pos[old]); vel.push(self.vel[old]); }
285                None => {
286                    match seed.get(id) {
287                        Some(&at) => pos.push(at),
288                        None => {
289                            let ring = if cold_start {
290                                spawned_here
291                            } else {
292                                self.spawn_seq % ARRIVAL_RINGS
293                            };
294                            pos.push(Self::phyllotaxis(self.spawn_seq, ring, centroid));
295                            self.spawn_seq += 1;
296                            spawned_here += 1;
297                        }
298                    }
299                    vel.push([0.0, 0.0]);
300                }
301            }
302            index.insert(id.clone(), i);
303        }
304
305        self.links = graph
306            .edges()
307            .filter_map(|e| Some((*index.get(&e.source)?, *index.get(&e.target)?)))
308            .collect();
309        // `Graph::edges()` iterates a HashMap, so its order differs between two
310        // independently built copies of the same graph. Float addition is not
311        // associative, so an unsorted `links` makes the spring accumulation —
312        // and therefore the whole settled layout — vary run to run. Sorting
313        // costs one pass per topology change and restores the reproducibility
314        // the sorted `ids` above was already reaching for. Equal pairs are
315        // interchangeable (identical contributions), so a stable tie is fine.
316        self.links.sort_unstable();
317
318        self.ids = ids;
319        self.pos = pos;
320        self.vel = vel;
321        self.index = index;
322        // Pins need no re-resolution — they are keyed by id, and the ids did
323        // not change. They do need pruning: a pin is an instruction about a
324        // node that is present, so one whose node has left the graph is
325        // dropped rather than kept dormant. Otherwise a node that leaves and
326        // later returns would come back mysteriously frozen.
327        self.pinned.retain(|id, _| self.index.contains_key(id));
328    }
329
330    /// Deterministic spread for a node with no seed — a sunflower spiral, so
331    /// unseeded arrivals never stack on one point.
332    ///
333    /// `seq` drives the angle and `ring` the radius. On a cold start the two
334    /// are equal and this is the classic `30·√i` sunflower. For a later arrival
335    /// `ring` wraps (see [`ARRIVAL_RINGS`]) so the radius stays near the
336    /// centroid, while `seq` keeps climbing — the golden angle is an irrational
337    /// fraction of a turn, so a never-repeating angle is on its own enough to
338    /// guarantee two spawns never coincide.
339    fn phyllotaxis(seq: usize, ring: usize, center: [f32; 2]) -> [f32; 2] {
340        let r = SPAWN_RADIUS * (ring as f32 + 0.5).sqrt();
341        let a = seq as f32 * 2.399_963; // golden angle
342        [center[0] + r * a.cos(), center[1] + r * a.sin()]
343    }
344
345    fn centroid(&self) -> [f32; 2] {
346        if self.pos.is_empty() { return [0.0, 0.0]; }
347        let n = self.pos.len() as f32;
348        let s = self.pos.iter().fold([0.0f32, 0.0], |a, p| [a[0] + p[0], a[1] + p[1]]);
349        [s[0] / n, s[1] / n]
350    }
351
352    /// Re-energise after a graph change. Sets `alpha` unconditionally —
353    /// `enabled` is enforced solely by `is_active`/`tick`, so a disabled sim
354    /// still accumulates the alpha a later `settle_forced` needs to have
355    /// anything to act on (reduced-motion and physics-disabled hosts reheat
356    /// exactly like a live one, then force a settle instead of ticking).
357    pub fn reheat(&mut self) {
358        self.alpha = self.spec.reheat_alpha;
359    }
360
361    /// Re-energise from scratch (`alpha = 1`), for a cold start. See `reheat`
362    /// on why this does not gate on `enabled`.
363    pub fn reheat_full(&mut self) {
364        self.alpha = 1.0;
365    }
366
367    /// Hold `id` at its current position. Unknown ids are ignored — a pin
368    /// describes a node that is present now, and does not arm itself for one
369    /// that arrives later.
370    pub fn pin(&mut self, id: &str) {
371        if let Some(&i) = self.index.get(id) {
372            self.pinned.insert(id.to_string(), self.pos[i]);
373        }
374    }
375
376    /// Hold `id` at `at`, moving it there immediately. Unknown ids are ignored.
377    ///
378    /// This is the drag entry point: each pointer move re-places the hold, and
379    /// the surrounding graph follows through the ordinary force integration.
380    /// Velocity is cleared so the node does not carry momentum from wherever
381    /// the simulation was pushing it before the grab.
382    pub fn pin_at(&mut self, id: &str, at: [f32; 2]) {
383        if let Some(&i) = self.index.get(id) {
384            self.pos[i] = at;
385            self.vel[i] = [0.0, 0.0];
386            self.pinned.insert(id.to_string(), at);
387        }
388    }
389
390    /// Release `id`. Unknown or unpinned ids are ignored.
391    pub fn unpin(&mut self, id: &str) {
392        self.pinned.remove(id);
393    }
394
395    /// Release every pin.
396    pub fn unpin_all(&mut self) {
397        self.pinned.clear();
398    }
399
400    pub fn is_pinned(&self, id: &str) -> bool {
401        self.pinned.contains_key(id)
402    }
403
404    pub fn pinned_ids(&self) -> impl Iterator<Item = &NodeId> {
405        self.pinned.keys()
406    }
407
408    pub fn is_active(&self) -> bool { self.spec.enabled && self.alpha >= self.spec.alpha_min }
409
410    /// Whether `alpha` is still above the settle threshold — `is_active`
411    /// without the `enabled` term, so a disabled sim can be asked whether it
412    /// was frozen mid-settle (and so needs a `settle_forced` to finish).
413    pub fn is_unsettled(&self) -> bool { self.alpha >= self.spec.alpha_min }
414
415    pub fn position(&self, id: &str) -> Option<[f32; 2]> {
416        self.index.get(id).map(|&i| self.pos[i])
417    }
418
419    pub fn velocity(&self, id: &str) -> Option<[f32; 2]> {
420        self.index.get(id).map(|&i| self.vel[i])
421    }
422
423    pub fn positions(&self) -> impl Iterator<Item = (&NodeId, [f32; 2])> + '_ {
424        self.ids.iter().zip(self.pos.iter().copied())
425    }
426
427    /// Tick to convergence, bounded by `max` so a pathological graph cannot
428    /// hang the caller. Used by reduced-motion, which snaps rather than
429    /// animating.
430    pub fn settle(&mut self, max: usize) {
431        for _ in 0..max {
432            if !self.tick() { return; }
433        }
434    }
435
436    /// `settle`, but runs even when `spec.enabled` is false — for hosts that
437    /// need a finished layout without any animation (reduced motion, or
438    /// physics disabled while the Force layout source is active). The flag is
439    /// restored unconditionally before returning; `tick` runs no host code,
440    /// so a scoped set/restore is panic-safe in practice.
441    pub fn settle_forced(&mut self, max: usize) {
442        let was = self.spec.enabled;
443        self.spec.enabled = true;
444        self.settle(max);
445        self.spec.enabled = was;
446    }
447
448    /// Advance one step. Returns false once settled (`alpha < alpha_min`), so
449    /// a caller can stop scheduling frames when the graph is at rest.
450    pub fn tick(&mut self) -> bool {
451        if !self.is_active() {
452            return false;
453        }
454        let n = self.pos.len();
455        if n == 0 {
456            self.alpha = 0.0;
457            return false;
458        }
459
460        let mut force = vec![[0.0f32, 0.0]; n];
461
462        // Many-body repulsion, Barnes-Hut approximated. Rebuilt each tick
463        // because every position moved.
464        let tree = QuadTree::build(&self.pos);
465        for (i, f) in force.iter_mut().enumerate() {
466            let r = tree.repulsion(i, &self.pos, self.spec.theta, self.spec.charge_strength);
467            f[0] += r[0];
468            f[1] += r[1];
469        }
470
471        // Link springs: pull each edge toward `link_distance`.
472        //
473        // Both the stiffness and the split are normalised by degree, exactly as
474        // d3-force does, and it is load-bearing rather than cosmetic. A flat
475        // 0.5/0.5 split gives a degree-k node a summed stiffness of k/2, and
476        // this integrator is explicit Euler: once `friction * k/2 > 2` the node
477        // overshoots further every tick and the layout diverges. Measured with
478        // a flat split: stable to degree 18, exponential blow-up by 20, NaN at
479        // 26 — and the demo's seed graph has a node with 26 neighbours.
480        // `1/min(deg)` bounds each node's total stiffness at ~1 instead, which
481        // is stable for any degree.
482        let mut deg = vec![0u32; n];
483        for &(a, b) in &self.links {
484            deg[a] += 1;
485            deg[b] += 1;
486        }
487        for &(a, b) in &self.links {
488            let dx = self.pos[b][0] - self.pos[a][0];
489            let dy = self.pos[b][1] - self.pos[a][1];
490            let d = (dx * dx + dy * dy).sqrt().max(1e-3);
491            let strength = 1.0 / deg[a].min(deg[b]).max(1) as f32;
492            // Positive when too far apart, so both ends move toward each other.
493            let pull = (d - self.spec.link_distance) / d * strength;
494            // Share of the correction taken by `b`. The busier end moves less,
495            // so a hub stays put and its leaves swing around it.
496            let (da, db) = (deg[a] as f32, deg[b] as f32);
497            let share_b = da / (da + db);
498            force[a][0] += dx * pull * (1.0 - share_b);
499            force[a][1] += dy * pull * (1.0 - share_b);
500            force[b][0] -= dx * pull * share_b;
501            force[b][1] -= dy * pull * share_b;
502        }
503
504        // Weak centering, so disconnected components stay in frame.
505        let c = self.centroid();
506        for (i, f) in force.iter_mut().enumerate() {
507            f[0] += (c[0] - self.pos[i][0]) * self.spec.center_strength;
508            f[1] += (c[1] - self.pos[i][1]) * self.spec.center_strength;
509        }
510
511        // Integrate. Force is scaled by alpha so motion damps out with the
512        // schedule, and velocity is damped by friction — d3's Verlet-ish step.
513        let friction = 1.0 - self.spec.velocity_decay;
514        let alpha = self.alpha;
515        for ((p, v), f) in self.pos.iter_mut().zip(self.vel.iter_mut()).zip(force.iter()) {
516            v[0] = (v[0] + f[0] * alpha) * friction;
517            v[1] = (v[1] + f[1] * alpha) * friction;
518            p[0] += v[0];
519            p[1] += v[1];
520        }
521
522        // Pinned nodes exert force but never move. Restoring the position
523        // (rather than only zeroing velocity) is required: a pin that drifted
524        // a little each tick would ratchet its own hold point away.
525        //
526        // Destructured so `pinned` can be iterated while `pos`/`vel` are
527        // written; borrowing `self.pinned` whole would conflict.
528        let Self { pinned, index, pos, vel, .. } = self;
529        for (id, at) in pinned.iter() {
530            if let Some(&i) = index.get(id) {
531                pos[i] = *at;
532                vel[i] = [0.0, 0.0];
533            }
534        }
535
536        self.alpha -= self.alpha * self.spec.alpha_decay;
537        self.is_active()
538    }
539}
540
541#[cfg(test)]
542mod tests {
543    use super::*;
544    use graph_explorer_core::{Edge, Graph, Node};
545    use std::collections::HashMap;
546
547    fn graph_of(nodes: &[&str], edges: &[(&str, &str, &str)]) -> Graph {
548        let mut g = Graph::default();
549        for n in nodes { g.add_node(Node::new(*n)); }
550        for (id, a, b) in edges { g.add_edge(Edge::new(*id, *a, *b)); }
551        g
552    }
553
554    fn seeded(pairs: &[(&str, [f32; 2])]) -> HashMap<String, [f32; 2]> {
555        pairs.iter().map(|(k, v)| ((*k).to_string(), *v)).collect()
556    }
557
558    #[test]
559    fn is_unsettled_tracks_alpha_not_enabled() {
560        let mut sim = ForceSimulation::new(SimulationSpec::default());
561        sim.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &seeded(&[]));
562        sim.reheat_full();
563        assert!(sim.is_unsettled());
564        sim.settle(10_000);
565        assert!(!sim.is_unsettled());
566    }
567
568    #[test]
569    fn defaults_match_the_documented_schedule() {
570        let s = SimulationSpec::default();
571        assert!(s.enabled);
572        assert_eq!(s.alpha_min, 0.001);
573        assert_eq!(s.reheat_alpha, 0.3);
574        assert_eq!(s.theta, 0.9);
575        // ~300 ticks from alpha 1 to alpha_min, which is d3's schedule.
576        assert!((s.alpha_decay - 0.0228).abs() < 1e-4);
577    }
578
579    #[test]
580    fn spec_rejects_unknown_fields() {
581        // `deny_unknown_fields` is what turns a host's typo into an error
582        // instead of a silently ignored setting.
583        let r: Result<SimulationSpec, _> = serde_json::from_str(r#"{"charge": -30}"#);
584        assert!(r.is_err(), "unknown field must be rejected");
585    }
586
587    #[test]
588    fn alpha_decays_to_min_in_about_three_hundred_ticks_then_stops() {
589        let mut s = ForceSimulation::new(SimulationSpec::default());
590        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &HashMap::new());
591        s.reheat_full();
592        let mut ticks = 0;
593        while s.tick() { ticks += 1; assert!(ticks < 1000, "runaway"); }
594        assert!((250..=350).contains(&ticks), "settled in {ticks} ticks, expected ~300");
595        assert!(!s.is_active(), "a settled sim must report inactive");
596        assert!(!s.tick(), "ticking a settled sim stays settled");
597    }
598
599    #[test]
600    fn reheat_reenergises_a_settled_sim() {
601        let mut s = ForceSimulation::new(SimulationSpec::default());
602        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &HashMap::new());
603        while s.tick() {}
604        assert!(!s.is_active());
605        s.reheat();
606        assert!(s.is_active(), "reheat must restart the schedule");
607    }
608
609    #[test]
610    fn sync_preserves_position_and_velocity_for_known_nodes() {
611        // Warm start is the whole reason a one-shot layout could not be reused:
612        // without it, every reveal teleports the graph.
613        let mut s = ForceSimulation::new(SimulationSpec::default());
614        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &seeded(&[("a", [100.0, 0.0]), ("b", [-100.0, 0.0])]));
615        s.reheat();
616        s.tick();
617        let a_before = s.position("a").unwrap();
618        let v_before = s.velocity("a").unwrap();
619        assert!(v_before[0] != 0.0 || v_before[1] != 0.0, "expected motion to compare against");
620
621        // A third node arrives — a reveal.
622        s.sync(&graph_of(&["a", "b", "c"], &[("e", "a", "b")]), &seeded(&[("c", [0.0, 50.0])]));
623        assert_eq!(s.position("a").unwrap(), a_before, "existing position preserved");
624        assert_eq!(s.velocity("a").unwrap(), v_before, "existing velocity preserved");
625        assert_eq!(s.position("c").unwrap(), [0.0, 50.0], "new node enters at its seed");
626    }
627
628    #[test]
629    fn sync_drops_nodes_that_left_the_graph() {
630        let mut s = ForceSimulation::new(SimulationSpec::default());
631        s.sync(&graph_of(&["a", "b"], &[]), &HashMap::new());
632        assert!(s.position("b").is_some());
633        s.sync(&graph_of(&["a"], &[]), &HashMap::new());
634        assert!(s.position("b").is_none(), "a removed node must leave the sim");
635    }
636
637    #[test]
638    fn an_unseeded_new_node_still_gets_a_finite_distinct_position() {
639        let mut s = ForceSimulation::new(SimulationSpec::default());
640        s.sync(&graph_of(&["a", "b", "c"], &[]), &HashMap::new());
641        let ps: Vec<[f32; 2]> = ["a", "b", "c"].iter().map(|i| s.position(i).unwrap()).collect();
642        for p in &ps { assert!(p[0].is_finite() && p[1].is_finite()); }
643        assert!(ps[0] != ps[1] && ps[1] != ps[2], "unseeded nodes must not stack");
644    }
645
646    #[test]
647    fn a_pinned_node_holds_still_while_its_neighbours_move() {
648        let mut s = ForceSimulation::new(SimulationSpec::default());
649        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &seeded(&[("a", [0.0, 0.0]), ("b", [400.0, 0.0])]));
650        s.pin("a");
651        s.reheat_full();
652        for _ in 0..30 { s.tick(); }
653        assert_eq!(s.position("a").unwrap(), [0.0, 0.0], "pinned node must not drift");
654        assert!(s.position("b").unwrap()[0] < 400.0, "its neighbour should have been pulled in");
655    }
656
657    #[test]
658    fn a_pinned_node_still_exerts_force() {
659        // Pinning fixes position, not influence — otherwise the graph would
660        // collapse through the node the user is standing on.
661        //
662        // `b` starts *off* both axes so this is a test about force rather than
663        // about motion. Asserting only that `b` moved, or that one coordinate
664        // grew, is satisfied by a model that slides the whole graph a constant
665        // step each tick with no repulsion at all. Repulsion from a fixed `a`
666        // has a signature that drift does not: it is purely radial, so `b` must
667        // travel straight out along the a->b ray and stay on it.
668        let mut s = ForceSimulation::new(SimulationSpec::default());
669        s.sync(&graph_of(&["a", "b"], &[]), &seeded(&[("a", [0.0, 0.0]), ("b", [3.0, 4.0])]));
670        s.pin("a");
671        s.reheat_full();
672        for _ in 0..20 { s.tick(); }
673
674        let (a, b) = (s.position("a").unwrap(), s.position("b").unwrap());
675        assert_eq!(a, [0.0, 0.0], "the pinned node must not have moved");
676
677        let d = ((b[0] - a[0]).powi(2) + (b[1] - a[1]).powi(2)).sqrt();
678        assert!(d > 5.0, "b sits {d} from a, no further than it started; nothing repelled it");
679
680        // `a` is at the origin, so the ray is just `b` normalised: (0.6, 0.8).
681        let unit = [(b[0] - a[0]) / d, (b[1] - a[1]) / d];
682        assert!(
683            (unit[0] - 0.6).abs() < 0.01 && (unit[1] - 0.8).abs() < 0.01,
684            "b reached {b:?}, off the a->b ray (direction {unit:?}); that is drift, not repulsion"
685        );
686    }
687
688    /// `sync` re-resolves the pin through the node id, because sorting the ids
689    /// reshuffles every index as soon as one is added or removed. That is the
690    /// subtlest logic in the warm start and it fails *silently* — a wrong index
691    /// holds the wrong node rather than panicking — so it needs its own test.
692    #[test]
693    fn a_pin_survives_the_index_shift_from_a_new_node() {
694        let mut s = ForceSimulation::new(SimulationSpec::default());
695        s.sync(&graph_of(&["m", "z"], &[]), &seeded(&[("m", [10.0, 0.0]), ("z", [200.0, 0.0])]));
696        s.pin("m"); // m is index 0 of [m, z]
697
698        // "a" sorts before "m", so m slides from index 0 to index 1.
699        s.sync(&graph_of(&["a", "m", "z"], &[]), &seeded(&[("a", [0.0, 300.0])]));
700        s.reheat_full();
701        for _ in 0..60 { s.tick(); }
702
703        assert_eq!(
704            s.position("m").unwrap(), [10.0, 0.0],
705            "the pin must still hold m, at m's own hold point, after the reindex"
706        );
707        assert!(s.position("z").unwrap() != [200.0, 0.0], "z was not pinned and should have moved");
708        assert!(s.position("a").unwrap() != [0.0, 300.0], "a was not pinned and should have moved");
709    }
710
711    /// A pinned node that leaves the graph must drop the pin, not hand it to
712    /// whichever survivor inherits its index — that would freeze an innocent
713    /// node the user never pinned.
714    #[test]
715    fn a_pinned_node_leaving_the_graph_drops_the_pin() {
716        let mut s = ForceSimulation::new(SimulationSpec::default());
717        s.sync(
718            &graph_of(&["m", "q", "z"], &[]),
719            &seeded(&[("m", [0.0, 0.0]), ("q", [40.0, 0.0]), ("z", [0.0, 40.0])]),
720        );
721        s.pin("m"); // index 0 of [m, q, z]
722
723        // m leaves; "q" inherits index 0. A pin tracked by index would hold q.
724        s.sync(&graph_of(&["q", "z"], &[]), &HashMap::new());
725        s.reheat_full();
726        for _ in 0..60 { s.tick(); }
727
728        assert!(s.position("m").is_none(), "m left the graph");
729        assert!(
730            s.position("q").unwrap() != [40.0, 0.0],
731            "q inherited m's index and is frozen; the pin transferred instead of dropping"
732        );
733        assert!(s.position("z").unwrap() != [0.0, 40.0], "z should be free to move");
734    }
735
736    #[test]
737    fn unpinning_releases_the_node() {
738        let mut s = ForceSimulation::new(SimulationSpec::default());
739        s.sync(&graph_of(&["a", "b"], &[]), &seeded(&[("a", [0.0, 0.0]), ("b", [5.0, 0.0])]));
740        s.pin("a");
741        s.reheat_full();
742        for _ in 0..10 { s.tick(); }
743        s.unpin("a");
744        s.reheat_full();
745        for _ in 0..20 { s.tick(); }
746        assert!(s.position("a").unwrap() != [0.0, 0.0], "released node should move");
747    }
748
749    #[test]
750    fn pinning_an_unknown_id_is_a_no_op_not_a_panic() {
751        let mut s = ForceSimulation::new(SimulationSpec::default());
752        s.sync(&graph_of(&["a"], &[]), &HashMap::new());
753        s.pin("ghost-that-is-not-here");
754        s.reheat_full();
755        s.tick();
756        assert!(s.position("a").unwrap()[0].is_finite());
757    }
758
759    /// The property the old single `fixed: Option<usize>` slot could not
760    /// express at all. Navigation holds the node the user is standing on and a
761    /// drag holds the node under the cursor; both at once is the ordinary case,
762    /// not an edge case.
763    #[test]
764    fn two_nodes_can_be_pinned_at_once() {
765        let mut s = ForceSimulation::new(SimulationSpec::default());
766        s.sync(
767            &graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c")]),
768            &seeded(&[("a", [0.0, 0.0]), ("b", [400.0, 0.0]), ("c", [0.0, 400.0])]),
769        );
770        s.pin("a");
771        s.pin("c");
772        s.reheat_full();
773        for _ in 0..60 { s.tick(); }
774
775        assert_eq!(s.position("a").unwrap(), [0.0, 0.0], "a must hold");
776        assert_eq!(s.position("c").unwrap(), [0.0, 400.0], "c must hold");
777        assert!(s.position("b").unwrap() != [400.0, 0.0], "b was free and should have moved");
778    }
779
780    #[test]
781    fn unpinning_one_leaves_the_others_held() {
782        let mut s = ForceSimulation::new(SimulationSpec::default());
783        s.sync(
784            &graph_of(&["a", "b", "c"], &[("e", "a", "b")]),
785            &seeded(&[("a", [0.0, 0.0]), ("b", [400.0, 0.0]), ("c", [0.0, 400.0])]),
786        );
787        s.pin("a");
788        s.pin("c");
789        s.unpin("a");
790        s.reheat_full();
791        for _ in 0..60 { s.tick(); }
792
793        assert!(s.position("a").unwrap() != [0.0, 0.0], "a was released and should move");
794        assert_eq!(s.position("c").unwrap(), [0.0, 400.0], "c must still hold");
795        assert!(!s.is_pinned("a"));
796        assert!(s.is_pinned("c"));
797    }
798
799    /// `pin_at` is the drag entry point: the hold point moves with the cursor,
800    /// and the node tracks it rather than snapping back to where it was grabbed.
801    #[test]
802    fn pin_at_moves_the_hold_point() {
803        let mut s = ForceSimulation::new(SimulationSpec::default());
804        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &seeded(&[("a", [0.0, 0.0]), ("b", [120.0, 0.0])]));
805        s.pin("a");
806        s.reheat_full();
807        for _ in 0..10 { s.tick(); }
808        assert_eq!(s.position("a").unwrap(), [0.0, 0.0]);
809
810        // Three "pointer moves", each followed by ticks, as a real drag would.
811        for step in 1..=3 {
812            let at = [step as f32 * 50.0, step as f32 * 10.0];
813            s.pin_at("a", at);
814            for _ in 0..5 { s.tick(); }
815            assert_eq!(s.position("a").unwrap(), at, "a must sit exactly where it is held");
816        }
817        // ...and the neighbour followed rather than staying put.
818        assert!(s.position("b").unwrap() != [120.0, 0.0], "b should have been dragged along");
819    }
820
821    /// A dropped pin must not lie in wait. A node that leaves the graph and
822    /// later comes back is a new arrival, not a resumption — coming back
823    /// mysteriously frozen would be a silent, hard-to-attribute bug.
824    #[test]
825    fn a_pin_does_not_resurrect_when_its_node_returns() {
826        let mut s = ForceSimulation::new(SimulationSpec::default());
827        s.sync(&graph_of(&["m", "z"], &[]), &seeded(&[("m", [10.0, 0.0]), ("z", [200.0, 0.0])]));
828        s.pin("m");
829        assert!(s.is_pinned("m"));
830
831        s.sync(&graph_of(&["z"], &[]), &HashMap::new()); // m leaves; pin is pruned
832        assert!(!s.is_pinned("m"), "the pin should have been dropped with the node");
833
834        s.sync(&graph_of(&["m", "z"], &[]), &seeded(&[("m", [10.0, 0.0])])); // m returns
835        assert!(!s.is_pinned("m"), "a returning node must not come back pinned");
836        s.reheat_full();
837        for _ in 0..60 { s.tick(); }
838        assert!(s.position("m").unwrap() != [10.0, 0.0], "m should be free to move");
839    }
840
841    #[test]
842    fn unpin_all_releases_everything() {
843        let mut s = ForceSimulation::new(SimulationSpec::default());
844        s.sync(&graph_of(&["a", "b"], &[]), &seeded(&[("a", [0.0, 0.0]), ("b", [50.0, 0.0])]));
845        s.pin("a");
846        s.pin("b");
847        s.unpin_all();
848        assert_eq!(s.pinned_ids().count(), 0);
849        s.reheat_full();
850        for _ in 0..40 { s.tick(); }
851        assert!(s.position("a").unwrap() != [0.0, 0.0]);
852        assert!(s.position("b").unwrap() != [50.0, 0.0]);
853    }
854
855    #[test]
856    fn disabled_spec_never_activates() {
857        let spec = SimulationSpec { enabled: false, ..SimulationSpec::default() };
858        let mut s = ForceSimulation::new(spec);
859        s.sync(&graph_of(&["a", "b"], &[]), &HashMap::new());
860        s.reheat();
861        assert!(!s.is_active(), "a disabled sim must never run");
862        assert!(!s.tick());
863    }
864
865    #[test]
866    fn settle_converges_and_respects_its_bound() {
867        let mut s = ForceSimulation::new(SimulationSpec::default());
868        s.sync(&graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c")]), &HashMap::new());
869        s.reheat_full();
870        s.settle(1000);
871        assert!(!s.is_active(), "settle should reach convergence");
872
873        s.reheat_full();
874        s.settle(5);
875        assert!(s.is_active(), "a bounded settle must stop early rather than spin");
876    }
877
878    #[test]
879    fn settle_forced_runs_while_disabled_and_restores_the_flag() {
880        let spec = SimulationSpec { enabled: false, ..SimulationSpec::default() };
881        let mut s = ForceSimulation::new(spec);
882        s.sync(&graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c")]), &HashMap::new());
883        s.reheat_full();
884        let before: Vec<[f32; 2]> = s.positions().map(|(_, p)| p).collect();
885        s.settle_forced(600);
886        assert!(!s.spec().enabled, "enabled flag must be restored");
887        assert!(!s.is_active(), "disabled sim reports inactive after the forced settle");
888        let after: Vec<[f32; 2]> = s.positions().map(|(_, p)| p).collect();
889        assert_ne!(before, after, "forced settle must actually move nodes");
890    }
891
892    #[test]
893    fn settle_forced_on_an_enabled_sim_behaves_like_settle() {
894        let mut s = ForceSimulation::new(SimulationSpec::default());
895        s.sync(&graph_of(&["a", "b"], &[("e1", "a", "b")]), &HashMap::new());
896        s.reheat_full();
897        s.settle_forced(1000);
898        assert!(s.spec().enabled);
899        assert!(!s.is_active(), "converged");
900    }
901
902    #[test]
903    fn a_triangle_settles_to_roughly_equilateral() {
904        let mut s = ForceSimulation::new(SimulationSpec::default());
905        let g = graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c"), ("e3", "c", "a")]);
906        s.sync(&g, &seeded(&[("a", [0.0, 0.0]), ("b", [10.0, 0.0]), ("c", [0.0, 10.0])]));
907        s.reheat_full();
908        s.settle(2000);
909
910        let d = |x: &str, y: &str| {
911            let (p, q) = (s.position(x).unwrap(), s.position(y).unwrap());
912            ((p[0] - q[0]).powi(2) + (p[1] - q[1]).powi(2)).sqrt()
913        };
914        let (ab, bc, ca) = (d("a", "b"), d("b", "c"), d("c", "a"));
915        let mean = (ab + bc + ca) / 3.0;
916        for (name, side) in [("ab", ab), ("bc", bc), ("ca", ca)] {
917            assert!((side - mean).abs() / mean < 0.15, "{name}={side} vs mean {mean}: not equilateral");
918        }
919        assert!(mean > 20.0, "sides collapsed to {mean}; repulsion is not holding them apart");
920    }
921
922    #[test]
923    fn linked_nodes_are_pulled_toward_link_distance() {
924        let mut s = ForceSimulation::new(SimulationSpec::default());
925        // Start far apart; the spring should reel them in.
926        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &seeded(&[("a", [0.0, 0.0]), ("b", [900.0, 0.0])]));
927        s.reheat_full();
928        s.settle(2000);
929        let (p, q) = (s.position("a").unwrap(), s.position("b").unwrap());
930        let d = ((p[0] - q[0]).powi(2) + (p[1] - q[1]).powi(2)).sqrt();
931        assert!(d < 400.0, "linked nodes stayed {d} apart; the spring is not pulling");
932    }
933
934    #[test]
935    fn a_disconnected_component_does_not_drift_to_infinity() {
936        // Pure repulsion pushes unlinked components apart forever. The
937        // centering force is what bounds them, and nothing else tests it.
938        let mut s = ForceSimulation::new(SimulationSpec::default());
939        s.sync(&graph_of(&["a", "b", "c", "d"], &[("e", "a", "b")]), &HashMap::new());
940        s.reheat_full();
941        s.settle(3000);
942        for id in ["a", "b", "c", "d"] {
943            let p = s.position(id).unwrap();
944            assert!(p[0].is_finite() && p[1].is_finite(), "{id} went non-finite");
945            assert!(p[0].abs() < 5000.0 && p[1].abs() < 5000.0, "{id} drifted to {p:?}");
946        }
947    }
948
949    /// `a_disconnected_component_does_not_drift_to_infinity` does not actually
950    /// test the centering force: with `center_strength = 0` its worst
951    /// coordinate is 495, well inside the 5000 bound, so it passes either way.
952    /// One alpha schedule is simply too short for the drift to show. Repeated
953    /// reheats are what separate the two: unbounded growth without centering
954    /// (measured 504 -> 972 -> 1269 -> 1504 over ten reheats, still climbing)
955    /// versus a fixed point with it (291 -> 269 -> 265 -> 265).
956    #[test]
957    fn centering_bounds_spread_across_repeated_reheats() {
958        let mut s = ForceSimulation::new(SimulationSpec::default());
959        s.sync(&graph_of(&["a", "b", "c", "d"], &[("e", "a", "b")]), &HashMap::new());
960        for _ in 0..10 {
961            s.reheat_full();
962            s.settle(3000);
963        }
964        let c = s.centroid();
965        let radius = s
966            .positions()
967            .map(|(_, p)| ((p[0] - c[0]).powi(2) + (p[1] - c[1]).powi(2)).sqrt())
968            .fold(0.0f32, f32::max);
969        assert!(
970            radius < 600.0,
971            "spread reached {radius} after ten reheats; the centering force is not bounding it"
972        );
973    }
974
975    #[test]
976    fn an_isolated_single_node_does_not_move() {
977        let mut s = ForceSimulation::new(SimulationSpec::default());
978        s.sync(&graph_of(&["only"], &[]), &seeded(&[("only", [42.0, -7.0])]));
979        s.reheat_full();
980        s.settle(500);
981        let p = s.position("only").unwrap();
982        assert!((p[0] - 42.0).abs() < 0.01 && (p[1] + 7.0).abs() < 0.01,
983            "a lone node has nothing to push against and should sit still, got {p:?}");
984    }
985
986    #[test]
987    fn identical_inputs_produce_identical_layouts() {
988        // The old one-shot layout was seeded for reproducibility; the sim must not
989        // regress that, or tests and screenshots become non-deterministic.
990        let run = || {
991            let mut s = ForceSimulation::new(SimulationSpec::default());
992            s.sync(&graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c")]), &HashMap::new());
993            s.reheat_full();
994            s.settle(500);
995            ["a", "b", "c"].iter().map(|i| s.position(i).unwrap()).collect::<Vec<_>>()
996        };
997        assert_eq!(run(), run());
998    }
999
1000    /// A flat 0.5/0.5 spring split gives a degree-k node a summed stiffness of
1001    /// k/2, which explicit Euler cannot integrate: measured, it was stable to
1002    /// degree 18, diverging by 20, and NaN at 26. The demo's seed graph has a
1003    /// node with 26 neighbours, so this is the shipping configuration, not a
1004    /// pathological one.
1005    #[test]
1006    fn a_high_degree_hub_does_not_explode() {
1007        let ids: Vec<String> = (0..=26).map(|i| format!("n{i:03}")).collect();
1008        let refs: Vec<&str> = ids.iter().map(|s| s.as_str()).collect();
1009        let edges: Vec<(String, String, String)> = (1..=26)
1010            .map(|i| (format!("e{i}"), "n000".to_string(), format!("n{i:03}")))
1011            .collect();
1012        let eref: Vec<(&str, &str, &str)> =
1013            edges.iter().map(|(i, a, b)| (i.as_str(), a.as_str(), b.as_str())).collect();
1014
1015        let mut s = ForceSimulation::new(SimulationSpec::default());
1016        s.sync(&graph_of(&refs, &eref), &HashMap::new());
1017        s.reheat_full();
1018        s.settle(3000);
1019
1020        let hub = s.position("n000").unwrap();
1021        assert!(hub[0].is_finite() && hub[1].is_finite(), "hub went non-finite: {hub:?}");
1022        for id in &refs {
1023            let p = s.position(id).unwrap();
1024            assert!(p[0].is_finite() && p[1].is_finite(), "{id} went non-finite: {p:?}");
1025            let r = ((p[0] - hub[0]).powi(2) + (p[1] - hub[1]).powi(2)).sqrt();
1026            assert!(r < 2000.0, "{id} flew to {r} from the hub; the spring is diverging");
1027        }
1028    }
1029
1030    /// Degree normalisation must not simply freeze every spring: a star's
1031    /// spokes should still be reeled in toward `link_distance`.
1032    #[test]
1033    fn a_hub_still_reels_its_spokes_to_roughly_link_distance() {
1034        let ids: Vec<String> = (0..=12).map(|i| format!("n{i:03}")).collect();
1035        let refs: Vec<&str> = ids.iter().map(|s| s.as_str()).collect();
1036        let edges: Vec<(String, String, String)> = (1..=12)
1037            .map(|i| (format!("e{i}"), "n000".to_string(), format!("n{i:03}")))
1038            .collect();
1039        let eref: Vec<(&str, &str, &str)> =
1040            edges.iter().map(|(i, a, b)| (i.as_str(), a.as_str(), b.as_str())).collect();
1041
1042        let mut s = ForceSimulation::new(SimulationSpec::default());
1043        s.sync(&graph_of(&refs, &eref), &HashMap::new());
1044        s.reheat_full();
1045        s.settle(3000);
1046
1047        let hub = s.position("n000").unwrap();
1048        for i in 1..=12 {
1049            let p = s.position(&format!("n{i:03}")).unwrap();
1050            let r = ((p[0] - hub[0]).powi(2) + (p[1] - hub[1]).powi(2)).sqrt();
1051            assert!((30.0..600.0).contains(&r), "spoke n{i:03} settled at {r} from the hub");
1052        }
1053    }
1054
1055    /// `Graph::edges()` iterates a HashMap, so two independently built copies
1056    /// of one graph can enumerate edges in different orders. Without a sort in
1057    /// `sync` the non-associative f32 spring sum diverges the layouts — this
1058    /// reproduced 26 times in 60 runs before `links.sort_unstable()`.
1059    #[test]
1060    fn link_order_is_independent_of_edge_hash_order() {
1061        let build = || {
1062            let mut s = ForceSimulation::new(SimulationSpec::default());
1063            s.sync(
1064                &graph_of(
1065                    &["a", "b", "c", "d", "e", "f"],
1066                    &[("e1", "a", "b"), ("e2", "b", "c"), ("e3", "c", "d"),
1067                      ("e4", "d", "e"), ("e5", "e", "f"), ("e6", "f", "a")],
1068                ),
1069                &HashMap::new(),
1070            );
1071            s.links.clone()
1072        };
1073        for _ in 0..50 {
1074            assert_eq!(build(), build(), "links order must not depend on hash iteration");
1075        }
1076    }
1077
1078    /// Build a wide ring of seeded nodes, so the centroid is well defined and
1079    /// far from any one member.
1080    fn big_seeded_ring(n: usize) -> (Vec<String>, HashMap<String, [f32; 2]>) {
1081        let ids: Vec<String> = (0..n).map(|k| format!("n{k:04}")).collect();
1082        let seed = ids
1083            .iter()
1084            .enumerate()
1085            .map(|(k, id)| {
1086                let a = k as f32 / n as f32 * std::f32::consts::TAU;
1087                (id.clone(), [400.0 * a.cos(), 400.0 * a.sin()])
1088            })
1089            .collect();
1090        (ids, seed)
1091    }
1092
1093    /// The spiral radius used to be `30·√(global rank)`, so a single unseeded
1094    /// arrival into a 500-node graph landed 671 units out — 5.6 link distances,
1095    /// nowhere near the "near the centroid" the doc contract promises.
1096    #[test]
1097    fn a_late_unseeded_arrival_lands_near_the_centroid() {
1098        let (ids, seed) = big_seeded_ring(500);
1099        let refs: Vec<&str> = ids.iter().map(|s| s.as_str()).collect();
1100        let mut s = ForceSimulation::new(SimulationSpec::default());
1101        s.sync(&graph_of(&refs, &[]), &seed);
1102        let c = s.centroid();
1103
1104        let mut with_new = refs.clone();
1105        with_new.push("zzz-arrives-last");
1106        s.sync(&graph_of(&with_new, &[]), &HashMap::new());
1107
1108        let p = s.position("zzz-arrives-last").unwrap();
1109        let r = ((p[0] - c[0]).powi(2) + (p[1] - c[1]).powi(2)).sqrt();
1110        let link = s.spec().link_distance;
1111        assert!(
1112            r < 2.0 * link,
1113            "arrival landed {r} from the centroid ({:.1} link distances), not near it",
1114            r / link
1115        );
1116    }
1117
1118    /// Radius came from the id's position in the *sorted* id list, so the same
1119    /// node entering the same graph spawned at r≈21 or r≈671 purely according
1120    /// to how its id happened to sort. Nothing physical justifies that.
1121    #[test]
1122    fn an_arrival_position_does_not_depend_on_where_its_id_sorts() {
1123        let radius_for = |new_id: &str| {
1124            let (ids, seed) = big_seeded_ring(500);
1125            let refs: Vec<&str> = ids.iter().map(|s| s.as_str()).collect();
1126            let mut s = ForceSimulation::new(SimulationSpec::default());
1127            s.sync(&graph_of(&refs, &[]), &seed);
1128            let c = s.centroid();
1129            let mut with_new = refs.clone();
1130            with_new.push(new_id);
1131            s.sync(&graph_of(&with_new, &[]), &HashMap::new());
1132            let p = s.position(new_id).unwrap();
1133            ((p[0] - c[0]).powi(2) + (p[1] - c[1]).powi(2)).sqrt()
1134        };
1135        // "aaa-new" sorts before every "nNNNN"; "zzz-new" sorts after them all.
1136        let (first, last) = (radius_for("aaa-new"), radius_for("zzz-new"));
1137        assert!(
1138            (first - last).abs() < 1e-3,
1139            "spawn radius depends on id ordering: {first} when sorting first, {last} when last"
1140        );
1141    }
1142
1143    /// The spawn counter lives on the struct, not on the `sync` call. A
1144    /// per-sync counter would hand every solitary arrival `seq = 0` and stack
1145    /// them all on one point — the precise failure the spiral exists to avoid.
1146    #[test]
1147    fn arrivals_in_separate_syncs_do_not_stack() {
1148        let mut s = ForceSimulation::new(SimulationSpec::default());
1149        s.sync(&graph_of(&["base"], &[]), &seeded(&[("base", [0.0, 0.0])]));
1150
1151        let mut live = vec!["base".to_string()];
1152        let mut arrivals = Vec::new();
1153        for k in 0..40 {
1154            let id = format!("late{k:03}");
1155            live.push(id.clone());
1156            let refs: Vec<&str> = live.iter().map(|x| x.as_str()).collect();
1157            s.sync(&graph_of(&refs, &[]), &HashMap::new()); // one arrival per sync
1158            arrivals.push(s.position(&id).unwrap());
1159        }
1160
1161        for i in 0..arrivals.len() {
1162            for j in (i + 1)..arrivals.len() {
1163                let (p, q) = (arrivals[i], arrivals[j]);
1164                let d = ((p[0] - q[0]).powi(2) + (p[1] - q[1]).powi(2)).sqrt();
1165                assert!(d > 1.0, "arrivals {i} and {j} spawned {d} apart, effectively stacked");
1166            }
1167        }
1168        // And they must all still be near the centroid, not spiralling outward.
1169        let link = s.spec().link_distance;
1170        for (k, p) in arrivals.iter().enumerate() {
1171            let r = (p[0] * p[0] + p[1] * p[1]).sqrt();
1172            assert!(r < 2.0 * link, "arrival {k} spawned {r} out; the radius is not bounded");
1173        }
1174    }
1175
1176    /// A cold start is the one case where spiral rank really is arrival count,
1177    /// and its `30·√N` spread is correct today. Bounding the *arrival* radius
1178    /// must not quietly shrink it.
1179    #[test]
1180    fn a_cold_start_keeps_its_full_spiral_spread() {
1181        let ids: Vec<String> = (0..500).map(|k| format!("n{k:04}")).collect();
1182        let refs: Vec<&str> = ids.iter().map(|s| s.as_str()).collect();
1183        let mut s = ForceSimulation::new(SimulationSpec::default());
1184        s.sync(&graph_of(&refs, &[]), &HashMap::new());
1185
1186        let c = s.centroid();
1187        let max = s
1188            .positions()
1189            .map(|(_, p)| ((p[0] - c[0]).powi(2) + (p[1] - c[1]).powi(2)).sqrt())
1190            .fold(0.0f32, f32::max);
1191        let expected = 30.0 * 499.5f32.sqrt(); // the outermost seat, ~670
1192        assert!(
1193            (max - expected).abs() < 10.0,
1194            "cold-start spread is {max}, expected ~{expected}; the sunflower was truncated"
1195        );
1196    }
1197
1198    /// Above theta ≈ 1.25 — and ≈ 0.75 for a tight cluster with a distant
1199    /// outlier, which is exactly the shape an expanded aggregate makes — a cell
1200    /// that contains body `i` can pass the approximation test, and `i` starts
1201    /// repelling itself through its own cell's mass. The spec boundary is where
1202    /// that gets excluded, so no per-tick branch has to.
1203    #[test]
1204    fn theta_is_clamped_when_the_spec_is_adopted() {
1205        let s = ForceSimulation::new(SimulationSpec { theta: 4.0, ..SimulationSpec::default() });
1206        assert_eq!(s.spec().theta, 1.0, "new() must clamp an out-of-range theta");
1207
1208        let mut t = ForceSimulation::new(SimulationSpec::default());
1209        t.set_spec(SimulationSpec { theta: -3.0, ..SimulationSpec::default() });
1210        assert_eq!(t.spec().theta, 0.0, "set_spec must clamp a negative theta to exact O(n^2)");
1211
1212        let u = ForceSimulation::new(SimulationSpec { theta: 0.9, ..SimulationSpec::default() });
1213        assert_eq!(u.spec().theta, 0.9, "an in-range theta must be left exactly alone");
1214    }
1215
1216    #[test]
1217    fn the_alpha_schedule_is_clamped_when_the_spec_is_adopted() {
1218        // Both entry points, since `set_spec` is the one a host reaches in a
1219        // single call and is therefore the one that has to be safe.
1220        let s = ForceSimulation::new(SimulationSpec {
1221            alpha_min: 0.0,
1222            alpha_decay: 0.0,
1223            velocity_decay: -1.0,
1224            reheat_alpha: 7.0,
1225            ..SimulationSpec::default()
1226        });
1227        assert_eq!(s.spec().alpha_min, ALPHA_MIN_FLOOR);
1228        assert_eq!(s.spec().alpha_decay, ALPHA_DECAY_MIN);
1229        assert_eq!(s.spec().velocity_decay, VELOCITY_DECAY_MIN);
1230        assert_eq!(s.spec().reheat_alpha, 1.0);
1231
1232        let mut t = ForceSimulation::new(SimulationSpec::default());
1233        t.set_spec(SimulationSpec {
1234            alpha_min: 5.0,
1235            alpha_decay: 9.0,
1236            velocity_decay: 3.0,
1237            reheat_alpha: -2.0,
1238            ..SimulationSpec::default()
1239        });
1240        assert_eq!(t.spec().alpha_min, ALPHA_MIN_CEIL);
1241        assert_eq!(t.spec().alpha_decay, ALPHA_DECAY_MAX);
1242        assert_eq!(t.spec().velocity_decay, VELOCITY_DECAY_MAX);
1243        assert_eq!(t.spec().reheat_alpha, 0.0);
1244
1245        // The defaults sit inside every bound and must survive untouched, or
1246        // the clamps would be silently retuning the shipped schedule.
1247        let d = ForceSimulation::new(SimulationSpec::default());
1248        let want = SimulationSpec::default();
1249        assert_eq!(d.spec().alpha_min, want.alpha_min);
1250        assert_eq!(d.spec().alpha_decay, want.alpha_decay);
1251        assert_eq!(d.spec().velocity_decay, want.velocity_decay);
1252        assert_eq!(d.spec().reheat_alpha, want.reheat_alpha);
1253    }
1254
1255    #[test]
1256    fn a_zero_alpha_min_cannot_wedge_the_sim_active_forever() {
1257        // `is_active` is `alpha >= alpha_min` and `alpha` decays multiplicatively
1258        // toward zero without reaching it, so an unclamped `alpha_min: 0` is
1259        // permanently true: the host's render loop would tick Barnes-Hut every
1260        // frame at rest, and the settle edge that callers rebuild their scene on
1261        // would never fire.
1262        let mut s = ForceSimulation::new(SimulationSpec {
1263            alpha_min: 0.0,
1264            ..SimulationSpec::default()
1265        });
1266        s.sync(&graph_of(&["a", "b"], &[("e", "a", "b")]), &HashMap::new());
1267        s.reheat_full();
1268        s.settle(5000);
1269        assert!(!s.is_active(), "a zero alpha_min must not leave the sim active forever");
1270        assert!(!s.tick(), "and the next tick must still report settled");
1271    }
1272
1273    #[test]
1274    fn every_legal_schedule_settles_within_the_reduced_motion_budget() {
1275        // The guarantee `ALPHA_DECAY_MIN` is derived from. A reduced-motion
1276        // caller settles synchronously with `settle(600)` and then snaps; if 600
1277        // ticks is not enough the sim is still active afterwards and animates
1278        // every frame — "no visible motion" inverted into permanent motion, plus
1279        // 600 wasted Barnes-Hut ticks per reheat. Each spec below pushes a
1280        // different schedule value to the end that used to break it.
1281        let hostile = [
1282            ("zero alpha_min", SimulationSpec { alpha_min: 0.0, ..Default::default() }),
1283            ("zero alpha_decay", SimulationSpec { alpha_decay: 0.0, ..Default::default() }),
1284            ("negative alpha_decay", SimulationSpec { alpha_decay: -1.0, ..Default::default() }),
1285            ("zero velocity_decay", SimulationSpec { velocity_decay: 0.0, ..Default::default() }),
1286            ("oversized reheat_alpha", SimulationSpec { reheat_alpha: 9.0, ..Default::default() }),
1287            (
1288                "all of them at once",
1289                SimulationSpec {
1290                    alpha_min: 0.0,
1291                    alpha_decay: 0.0,
1292                    velocity_decay: -3.0,
1293                    reheat_alpha: 9.0,
1294                    ..Default::default()
1295                },
1296            ),
1297        ];
1298        for (what, spec) in hostile {
1299            let mut s = ForceSimulation::new(spec);
1300            s.sync(
1301                &graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c")]),
1302                &HashMap::new(),
1303            );
1304            // `reheat_full`, not `reheat`: alpha 1.0 is the worst case the
1305            // budget has to cover.
1306            s.reheat_full();
1307            s.settle(600);
1308            assert!(!s.is_active(), "{what}: still active after settle(600)");
1309            assert!(!s.tick(), "{what}: the next tick must still report settled");
1310        }
1311    }
1312
1313    #[test]
1314    fn a_negative_velocity_decay_cannot_diverge_the_layout() {
1315        // `friction = 1 - velocity_decay`, so a negative decay makes friction
1316        // exceed 1 and this explicit-Euler step amplifies velocity every tick
1317        // instead of damping it — unclamped, the positions reach NaN.
1318        let mut s = ForceSimulation::new(SimulationSpec::default());
1319        s.set_spec(SimulationSpec { velocity_decay: -5.0, ..SimulationSpec::default() });
1320        s.sync(
1321            &graph_of(&["a", "b", "c"], &[("e1", "a", "b"), ("e2", "b", "c")]),
1322            &seeded(&[("a", [0.0, 0.0]), ("b", [10.0, 0.0]), ("c", [0.0, 10.0])]),
1323        );
1324        s.reheat_full();
1325        s.settle(2000);
1326        for id in ["a", "b", "c"] {
1327            let p = s.position(id).expect("node present");
1328            assert!(p[0].is_finite() && p[1].is_finite(), "{id} diverged to {p:?}");
1329        }
1330    }
1331
1332    #[test]
1333    fn the_alpha_decay_floor_leaves_headroom_under_the_settle_budget() {
1334        // Pins the derivation in `ALPHA_DECAY_MIN`'s doc, so that lowering the
1335        // floor (or raising `ALPHA_MIN_FLOOR`) without redoing the arithmetic
1336        // fails here rather than silently in reduced motion.
1337        let ticks = (ALPHA_MIN_FLOOR.ln() / (1.0 - ALPHA_DECAY_MIN).ln()).ceil();
1338        assert!(ticks <= 600.0, "worst-case settle is {ticks} ticks, over the 600 budget");
1339    }
1340
1341    #[test]
1342    fn coincident_start_positions_do_not_produce_nan() {
1343        let mut s = ForceSimulation::new(SimulationSpec::default());
1344        s.sync(&graph_of(&["a", "b", "c"], &[]),
1345               &seeded(&[("a", [0.0, 0.0]), ("b", [0.0, 0.0]), ("c", [0.0, 0.0])]));
1346        s.reheat_full();
1347        s.settle(200);
1348        for id in ["a", "b", "c"] {
1349            let p = s.position(id).unwrap();
1350            assert!(p[0].is_finite() && p[1].is_finite(), "{id} = {p:?}");
1351        }
1352    }
1353}