feat(70): the propagation tier — RWR, heat diffusion, resistance, SimRank, and the triad census

Round 69 asked which algorithms could take MCL's dense-matmul lead;
this round builds the five families the scientific-use assessment
ranked highest, aimed at network biology's propagation methods.  Seven
public methods, each with a CPU reference (the spec), executor
routing, closed-form Node specs, and — for the five dense forms — GPU
kernels, parity specs with verified controls, and benchmark rows.

The tier's structural decision: seed forms vs dense forms.
randomWalkWithRestart({seeds}) and heatDiffusion({seeds}) — the
everyday bio calls — are O(E)-per-step sparse walks, so they are
CPU-only and an explicit 'gpu' rejects, pointing at the dense sibling
(the weighted-betweenness no-path precedent).  The dense forms are
the GPU tier:

- randomWalkWithRestartProximity: S = c(I − (1−c)W)⁻¹ by Neumann
  iteration, one matmul per step; CPU = one sparse solve per column.
- heatKernel: exp(−tL) by scaling-and-squaring (Taylor chain + s
  squarings); CPU applies the same scaled operator 2^s times per
  column — the same power, sparsely.  Pinned against the pair and
  triangle matrix exponentials in closed form.
- effectiveResistance (+ commuteTime): B = L + J-blocks (the 1/n_c
  shift cancels out of every resistance), CPU f64 Gauss–Jordan, GPU
  Newton–Schulz — O(n³) both sides, so the GPU wins at every density.
  The NS converge compare is *relative* (NS_COMPARE): an absolute
  bound sits under f32 noise wherever the inverse's entries are large
  (measured: 96/96 iterations ran at n=1024; relative freezing took
  it 642 → 200 ms).
- simRank: S′ = C·Q·S·Qᵀ, two matmuls per step, diagonal pinned;
  pinned against the 4-cycle fixed point x = C(1+x)/2.
- motifCensus: the 16 Holland–Leinhardt classes ('030T' = the
  feed-forward loop) as closed forms over seven trace primitives
  shared verbatim by both executors — CPU wedge walks, GPU four
  matmuls + Hadamard folds.  The load-bearing spec is a brute-force
  differential: an independent classifier over every triple of six
  random digraphs, exact equality per class.

Density gates as in round 69 (sparse CPU walks own sparse graphs);
resistance gates on size alone.  Conventions: simple graph (parallel
edges collapse / sum weights, loops drop); RWR sinks absorb
(documented leak); heat and resistance demand positive weights
(TypeError).  Six new guards, all pinned; throw gate at zero unrun;
JSDoc gates 100%; 2245 Node specs green.

Parity: five new live specs, each proven able to fail by degrading
its kernel (C skewed, restart diagonal skewed, a Taylor term dropped,
2I skewed to 2.01I, a fold mask swapped) — all five controls failed,
all 20 specs restored green on apple metal-3.

Measured on the M2 (one-off; the archive run remains the RX 580's):
heatKernel 932× at n=1024 — the largest ratio any family has measured
— rwrProximity 119×, simRank 45×, resistance 9×, census 12× at 2048.
Bench rows sized so a cell stays in MCL's cost class (the dense CPU
references are minutes at n=2048; sizes stop at 1024, heat pins
time = 0.02).

Claude-Session: https://claude.ai/code/session_01T6ZDogLZVCQ7QRM1aviH4h

Author: maxkfranz

CHANGELOG.md

@@ -117,19 +117,34 @@ that compile and then behave differently.
   (O(E) per iteration — orders of magnitude on sparse graphs) and the
   hierarchical merge engine went flat-typed, leaving the GPU no edge
   to win there.
-- **Three new algorithm families, designed matmul-first for the GPU
-  tier** (round 69), all on the same async `executor` contract and
-  with no v3 counterpart: **`eles.triangleCount()`** (per-node
-  triangle counts, local clustering coefficients, total triangles and
-  transitivity — A²∘A on the GPU), **`eles.neighborhoodSimilarity()`**
-  (pairwise Jaccard / cosine / overlap coefficients over neighbor
-  sets — A·Aᵀ on the GPU) and **`eles.katzCentrality()`** (attenuated
-  walk counting; like PageRank its sparse CPU iteration owns `'auto'`
-  and the GPU path serves an explicit `'gpu'`).  The whole-collection
-  `closenessCentralityNormalized` joined the async tier in the same
-  round: its GPU path rides the blocked Floyd–Warshall kernels and
-  folds each distance row on the device, reading back n floats
-  instead of the n² matrix.
+- **Eight new algorithm families, designed matmul-first for the GPU
+  tier** (rounds 69–70), all on the same async `executor` contract
+  and with no v3 counterpart.  Round 69: **`eles.triangleCount()`**
+  (per-node triangle counts, local clustering coefficients, total
+  triangles and transitivity — A²∘A on the GPU),
+  **`eles.neighborhoodSimilarity()`** (pairwise Jaccard / cosine /
+  overlap coefficients over neighbor sets — A·Aᵀ on the GPU) and
+  **`eles.katzCentrality()`** (attenuated walk counting; like
+  PageRank its sparse CPU iteration owns `'auto'` and the GPU path
+  serves an explicit `'gpu'`).  Round 70, aimed at network-biology
+  workloads: **`eles.randomWalkWithRestart()`** (seed-set network
+  propagation — the disease-gene-prioritization primitive) and
+  **`eles.randomWalkWithRestartProximity()`** (the all-pairs
+  proximity matrix, a Neumann matmul iteration on the GPU),
+  **`eles.heatDiffusion()`** / **`eles.heatKernel()`** (HotNet-style
+  heat propagation; exp(−t·L) by scaling-and-squaring on the GPU),
+  **`eles.effectiveResistance()`** (resistance distance and commute
+  time off the Laplacian pseudo-inverse — f64 elimination on the CPU,
+  Newton–Schulz matmul iteration on the GPU; O(n³) on both sides, so
+  the GPU wins at every density), **`eles.simRank()`** (the Jeh–Widom
+  recursive similarity, two matmuls per iteration) and
+  **`eles.motifCensus()`** (the sixteen-class Holland–Leinhardt triad
+  census — '030T' is the feed-forward loop — computed from seven
+  trace primitives and pinned by a brute-force classifier spec).
+  The whole-collection `closenessCentralityNormalized` joined the
+  async tier in round 69: its GPU path rides the blocked
+  Floyd–Warshall kernels and folds each distance row on the device,
+  reading back n floats instead of the n² matrix.
 
 ### Changed
 

EXECUTIVE_SUMMARY.md

@@ -38,9 +38,9 @@ The v4 rewrite: a columnar model and a WebGPU renderer, per
 
 | | |
 |---|---|
-| Automated tests | 2,223 unit · 427 module · 24 soak · 386 browser (some skip for want of a WebGPU adapter) |
-| Documented API | 366 members over 48 sections, gated at 100% |
-| Visual regression | 46 goldens compared **exactly** — zero differing pixels · 45 live v3-vs-v4 pixel-parity scenes, 7 of them close-ups at zoom 3–4 · 11 numeric routing-parity scenes · 15 CPU-vs-GPU algorithm-parity scenes |
+| Automated tests | 2,245 unit · 427 module · 24 soak · 396 browser (some skip for want of a WebGPU adapter) |
+| Documented API | 373 members over 48 sections, gated at 100% |
+| Visual regression | 46 goldens compared **exactly** — zero differing pixels · 45 live v3-vs-v4 pixel-parity scenes, 7 of them close-ups at zoom 3–4 · 11 numeric routing-parity scenes · 20 CPU-vs-GPU algorithm-parity scenes |
 | Benchmarks | 25 suites, 4 published profiles · **all 366 v3-comparative pairs read v4-faster** (geometric mean 13.7×, minimum 1.03×) · GPU algorithm executors 13× geo-mean over their CPU reference |
 | Style parity | v4 accepts 157 of v3's 291 style property names; the rest dropped by decision |
 | Bundle | 691 KiB minified / 185 KiB gzipped — ~1.5× v3 (410 / 126 KiB); the WGSL shaders, which v3 has no equivalent of, are minified at build time |
@@ -166,6 +166,20 @@ The v4 rewrite: a columnar model and a WebGPU renderer, per
   - Buys three algorithms v3 never had, on kernels the suite already trusts;
     each of the five new parity specs was proven able to fail by degrading its
     kernel; crossover numbers await the benchmark machine.
+- **12 Aug** — the propagation tier: network biology's algorithms
+  - Five more families, chosen by scientific usefulness: random walk with
+    restart (seed propagation — the disease-gene-prioritization primitive —
+    plus the all-pairs proximity matrix), heat-kernel diffusion (HotNet-style
+    exp(−tL), seed and all-pairs forms), effective resistance / commute time
+    (Laplacian pseudo-inverse: f64 elimination on the CPU, Newton–Schulz
+    matmuls on the GPU — O(n³) both sides, so the GPU wins at every density),
+    SimRank, and the sixteen-class triad census ('030T' is the biology
+    literature's feed-forward loop), whose closed forms are pinned by a
+    brute-force classify-every-triple spec.  Seed forms are CPU-only by
+    design — O(E) walks with nothing for a kernel to win.
+  - Buys the network-biology propagation toolbox on the existing kernel
+    machinery; five more parity specs, each proven able to fail; measured on
+    an M2: RWR proximity 119×, SimRank 45× at n=1024.
 
 ---
 
@@ -183,9 +197,11 @@ The v4 rewrite: a columnar model and a WebGPU renderer, per
   `executor: 'cpu' | 'gpu' | 'auto'` — including, since 12 Aug, the
   whole-collection `closenessCentralityNormalized`; the single-root form
   stays synchronous.
-- **Three algorithm families v3 never had**: `triangleCount`,
-  `neighborhoodSimilarity` and `katzCentrality`, on the same executor
-  contract.
+- **Eight algorithm families v3 never had**, on the same executor
+  contract: `triangleCount`, `neighborhoodSimilarity`, `katzCentrality`,
+  `randomWalkWithRestart` (+ its all-pairs proximity form),
+  `heatDiffusion`/`heatKernel`, `effectiveResistance`, `simRank` and
+  `motifCensus`.
 - **`cy.collection()` throws if passed an argument**; **`cy.$()` and
   `cy.byId()`** restored as aliases.
 

MIGRATING.md

@@ -424,14 +424,21 @@ slot-native CPU walks.  Only the whole-collection
 `closenessCentralityNormalized` moved to the async tier — it is the O(n³)
 all-pairs computation, exactly what the tier exists for.
 
-**v4 also adds three algorithm families v3 never had**, on the same async
+**v4 also adds eight algorithm families v3 never had**, on the same async
 `executor` contract: `triangleCount()` (per-node triangles, local
 clustering coefficients, transitivity), `neighborhoodSimilarity()`
-(pairwise Jaccard / cosine / overlap coefficients over neighbor sets) and
-`katzCentrality()` (attenuated walk counting).  Nothing to migrate — they
-are new surface — but note they read the collection as a simple graph:
-parallel edges collapse, loops are excluded, and `triangleCount` ignores
-direction outright.
+(pairwise Jaccard / cosine / overlap coefficients over neighbor sets),
+`katzCentrality()` (attenuated walk counting), `randomWalkWithRestart()`
+and `randomWalkWithRestartProximity()` (network propagation from a seed
+set, and the all-pairs proximity matrix), `heatDiffusion()` and
+`heatKernel()` (heat-kernel propagation, exp(−t·L)),
+`effectiveResistance()` (the graph as a resistor network, with
+`commuteTime`), `simRank()` (recursive neighborhood similarity) and
+`motifCensus()` (the sixteen-class triad census; '030T' is the
+feed-forward loop).  Nothing to migrate — they are new surface — but
+note they read the collection as a simple graph: parallel edges collapse
+(summing weights where weights are read), loops are excluded, and the
+triangle/heat/resistance families ignore direction outright.
 
 **State is a condition, not a selector.** See "Styling element state" above:
 `:selected`, `:active`, `:locked` and the rest are `when` conditions on a

PLAN.md

@@ -19854,3 +19854,158 @@ closeness should also gain a Brandes-style batched-BFS unweighted path
 (the FW route was chosen because the CPU reference is FW-shaped, so
 parity is structural — a BFS path would be faster still on sparse
 graphs but needs its own reference).
+
+## Round 70 — the propagation tier: network biology's algorithms (maintainer goal, 2026-08-12)
+
+Round 69 asked "which algorithms could take MCL's lead"; the follow-up
+question was which of the candidates matter for *scientific* use, and
+the answer reordered the queue: network biology has converged on
+propagation methods, so the round builds the five families the
+maintainer picked from that assessment — random walk with restart,
+heat-kernel diffusion, effective resistance / commute time, SimRank,
+and the triad census.  Seven public methods, each with a CPU reference
+(the spec), executor routing, hand-computed Node specs, and — for the
+five dense forms — GPU kernels, live parity specs with verified
+controls, and benchmark rows.
+
+### 70.1 — the seed/dense split, decided once for the tier
+
+RWR and heat diffusion each have two natural forms.  The *seed* form —
+propagate from a seed collection, answer a score per node — is the
+everyday bio call, and it is an O(E)-per-step sparse walk: the pageRank
+verdict applies, no kernel can win, so `randomWalkWithRestart` and
+`heatDiffusion` are CPU-only and an explicit `executor: 'gpu'` rejects
+with a message pointing at the dense sibling (the weighted-betweenness
+no-path precedent).  The *all-pairs* form — the full proximity/kernel
+matrix — is iterated dense products, the MCL shape, and that is where
+the GPU tier lives: `randomWalkWithRestartProximity` (Neumann
+iteration S′ = (1−c)W·S + cI, one matmul per step),
+`heatKernel` (exp(−tL) by scaling-and-squaring: the scaled Taylor
+series as a matmul chain, then s squarings), `simRank`
+(S′ = C·Q·S·Qᵀ, two matmuls per step) and the census's trace products.
+The CPU references for the all-pairs forms are deliberately *sparse*
+(one seed-vector solve per column; per-column Taylor applications;
+per-row/column neighbor averages), so the density gates are honest:
+the m ≥ n²/16–n²/32 gates of round 69, same reasoning, and 'auto'
+stays on the CPU for sparse graphs however large.
+
+`effectiveResistance` is the exception and the headline: the Laplacian
+pseudo-inverse has no sparse shortcut, so the CPU reference is dense
+f64 Gauss–Jordan at O(n³) and the GPU runs Newton–Schulz —
+X ← X(2I − BX), nothing but matmuls, quadratically convergent for the
+positive-definite B = L + J-blocks this family builds (the per-
+component 1/n_c shift makes B invertible while cancelling out of every
+resistance difference).  Like MCL it wins at every density.  f32
+bounds the achievable accuracy on ill-conditioned systems, documented;
+the parity bound is relative 5e-3 against the f64 elimination.
+
+### 70.2 — semantics worth recording
+
+- RWR: W column-normalized by out-weight; the undirected default walks
+  both ways; a directed sink *absorbs* (its column leaks, scores can
+  sum under 1) rather than redistributing — documented, simpler on
+  both executors, irrelevant on undirected bio graphs.  `seeds` is
+  required for the seed form and uniform over its nodes; the fixed
+  point is c(I − (1−c)W)⁻¹p₀, pinned in specs by the two-node closed
+  form c/(1−(1−c)²).
+- Heat: the combinatorial weighted Laplacian, undirected, positive
+  weights enforced (a negative conductance is not a heat problem —
+  TypeError).  Both executors share the same approximation constants
+  (‖tL/2^s‖∞ ≤ ½, ten Taylor terms); the CPU applies the scaled
+  operator 2^s times per column, the GPU squares s times — the same
+  power.  Specs pin the pair and triangle matrix exponentials in
+  closed form and heat conservation through the scaling path at t=10.
+- SimRank: Jeh–Widom with in-neighborhoods under `directed: true`,
+  all neighbors otherwise (the library's undirected default; bio
+  graphs are undirected).  Diagonal pinned to 1 per iteration; empty
+  neighborhoods answer 0.  The 4-cycle fixed point x = C(1+x)/2 →
+  x = 2/3 at C = 0.8 pins the maths in specs.
+- Census: sixteen closed forms over seven trace primitives
+  (S₁ = ΣC²∘C … S₇ = ΣCCᵀ∘M over the asymmetric and mutual masks),
+  the dyad totals and six degree-pair sums, shared verbatim by both
+  executors (`censusFromPrimitives`) — the executors can only
+  disagree if a matmul disagrees with a wedge walk.  **The formulas
+  themselves are the risk**, so the load-bearing spec is a
+  brute-force differential: an independent classifier written from
+  the class definitions, run over every triple of six random digraphs
+  sweeping sparse to dense, exact equality demanded per class (plus
+  Σ = C(n,3)).  It passed on the first complete run of the closed
+  forms, and it is the spec that would catch a sign or orientation
+  error nothing else can see.  `directed: false` files every edge as
+  mutual, so the undirected census (empty / one-edge / path /
+  triangle) is the same code path reading 003/102/201/300.
+
+### 70.3 — verification
+
+Node tier: `test/algorithms-propagation.mjs` (24 specs) and
+`test/algorithms-motifs.mjs` (15 assertions per seed across six
+seeds), all closed-form or brute-force; the executor sweep extended by
+seven entries; throw gate at zero unrun (six new guards: seeds,
+restartProbability, time, two positive-weight conductance guards,
+dampingFactor); JSDoc gates 100%; types and full `test:js` green.
+
+Parity tier: five new live specs — simRank (1e-4 plus an exact-1
+diagonal, both directions), rwrProximity (1e-4 plus column
+conservation ≥ 0.999), heatKernel (1e-4, symmetry, row conservation,
+at t = 2 so several squarings run), effectiveResistance (relative
+5e-3, the unit-resistor identity, and exact Infinity agreement across
+components on a fixture that has both), motifCensus (all sixteen
+counts exactly equal on a 64-node random digraph, with a
+populated-classes discrimination check).  **Every spec was run once
+with its kernel deliberately degraded and failed**: simRank with C
+skewed 1%, rwrProximity with the restart diagonal skewed, heatKernel
+with the k=2 Taylor term dropped, resistance with the Newton–Schulz
+2I skewed to 2.01I, the census with S₆ folded against the wrong mask.
+All five controls failed; all 20 specs green restored.
+
+### 70.4 — benchmark rows, and a sizing lesson
+
+Five families joined `algorithms-gpu-bench.mjs`: resistance on the
+plain fixture (both sides O(n³) — the MCL-class row), the four
+iterated-product families on the dense fixture.  The first draft
+priced the dense families at n = 2048 and had to be walked back: a
+*bench cell* pays REPS×(cpu+gpu) calls, and the dense CPU references
+are MCL-cost already at n = 1024 (rwrProximity is one sparse solve per
+column; heatKernel's scaling exponent grows with t·degree, so its
+bench row pins time = 0.02).  Sizes stop at 1024 for those three
+(census stays to 2048 — its CPU walk is O(Σ deg²), far cheaper), and
+the iteration knobs are pinned in the rows per the round-33.2 rule.
+The suite edit moves the harness fingerprint again; the round-69 note
+about `EQUIVALENT_HARNESSES` applies to the new hash the same way,
+and only one entry is needed once a run under the final hash is
+published.
+
+**Open**: the crossover sweep on the benchmark machine (density-gate
+constants re-tuned from it, and the resistance family's parity bound
+revisited on the RX 580's f32); a `normalized` Laplacian option for
+the heat family; motif significance tooling (the census exists so it
+can be run per randomized network — the ensemble driver itself is
+app-level and stays out of scope).
+
+### 70.5 — measured on the M2 (Metal, one-off; the archive run stays the RX 580's)
+
+| family | n=256 | n=512 | n=1024 |
+|---|---:|---:|---:|
+| effectiveResistance | 3.3× | 7.8× | 9.0× |
+| simRank (10 iters, dense) | 14× | 33× | 45× |
+| rwrProximity (dense) | 21× | 68× | **119×** |
+| heatKernel (dense) | 82× | 291× | **932×** |
+| motifCensus (dense, to n=2048) | — | 3.9× | 8.7× / 12.3× @2048 |
+
+heatKernel's 932× is the largest ratio any family has measured — its
+CPU reference pays 2^s operator applications per column while the GPU
+pays s squarings total, so the scaling exponent multiplies the CPU
+side only.  Two findings from the measurement worth their notes:
+
+- **A converge tolerance below f32's noise floor buys nothing and
+  costs everything.**  Newton–Schulz at an absolute 1e-5 ran 96 of 96
+  encoded iterations at n=1024 (642 ms): the inverse's entries grow as
+  1/λ₂, so on any weakly-connected graph the iterate's float noise
+  exceeds an absolute bound forever and the no-diff converge never
+  fires.  The compare is *relative* now (`NS_COMPARE`,
+  |Δ| > tol·max(1, |x|)) — 642 → 200 ms at n=1024, 9× over the f64
+  elimination, and the parity spec still passes at 5e-3 relative.
+- **A bench cell's budget is REPS × the slow side.**  The first row
+  draft priced the dense families at n = 2048 and a single
+  rwrProximity CPU call there runs minutes; sizes stopped at 1024 and
+  the heat row pinned time = 0.02 (its CPU cost scales with 2^s).

benchmark/algorithms-gpu-bench.mjs

@@ -276,6 +276,59 @@ const FAMILIES = [
       cy.elements().katzCentrality({ executor })
         .then((r) => r.katz(cy.nodes()[0]))`,
   },
+  // round 70: effective resistance is O(n³) on BOTH executors (a dense
+  // inverse has no sparse shortcut), so like MCL it is priced on the
+  // plain fixture and wins at every density
+  {
+    key: 'effectiveResistance',
+    sizes: [256, 512, 1024],
+    kind: 'graph',
+    op: `(cy, executor) =>
+      cy.elements().effectiveResistance({ executor })
+        .then((r) => r.resistance(cy.nodes()[0], cy.nodes()[1]))`,
+  },
+  // round 70: the iterated-product families, priced dense (their 'auto'
+  // gates route to the GPU only there — sparse CPU walks own the rest).
+  // Iteration knobs are pinned (the round-33.2 rule), and sizes stop at
+  // 1024 because the dense CPU references are MCL-cost there already —
+  // simRank pays 2·n·m per iteration and rwrProximity one sparse solve
+  // per column; heatKernel's time is small so its scaling exponent (and
+  // with it the CPU's 2^s operator applications) stays bounded
+  {
+    key: 'simRank',
+    sizes: [256, 512, 1024],
+    kind: 'graph-dense',
+    op: `(cy, executor) =>
+      cy.elements().simRank({ executor, maxIterations: 10 })
+        .then((r) => r.similarity(cy.nodes()[0], cy.nodes()[1]))`,
+  },
+  {
+    key: 'rwrProximity',
+    sizes: [256, 512, 1024],
+    kind: 'graph-dense',
+    op: `(cy, executor) =>
+      cy.elements().randomWalkWithRestartProximity({
+        executor,
+        restartProbability: 0.3,
+        tolerance: 0.00001,
+      }).then((r) => r.proximity(cy.nodes()[0], cy.nodes()[1]))`,
+  },
+  {
+    key: 'heatKernel',
+    sizes: [256, 512, 1024],
+    kind: 'graph-dense',
+    op: `(cy, executor) =>
+      cy.elements().heatKernel({ executor, time: 0.02 })
+        .then((r) => r.heat(cy.nodes()[0], cy.nodes()[1]))`,
+  },
+  {
+    key: 'motifCensus',
+    sizes: [512, 1024, 2048],
+    kind: 'graph-dense',
+    op: `(cy, executor) =>
+      cy.elements().motifCensus({ executor })
+        .then((r) => r.counts['030T'])`,
+  },
 ];
 
 const jobs = [];