From 8e49c867604d20abbc744fcfe5b0880dda1397fa Mon Sep 17 00:00:00 2001 From: bojanstef <5675392+bojanstef@users.noreply.github.com> Date: Thu, 18 Jun 2026 22:15:03 -0400 Subject: [PATCH 1/3] fix(physics): stop forceAtlas2Based layout from rotating forever (#2193) The FA2 repulsion scaled the pairwise force by the receiving node's degree only, making it non-reciprocal (violating Newton's third law) and injecting a net torque each tick that damping caps but never removes -> the converged layout spins forever. Remove the degree factor to restore reciprocity. Adds a deterministic regression test measuring net angular velocity about the COM. --- .../physics/FA2BasedRepulsionSolver.js | 4 +- test/physics/forceatlas2-rotation.test.ts | 169 ++++++++++++++++++ 2 files changed, 170 insertions(+), 3 deletions(-) create mode 100644 test/physics/forceatlas2-rotation.test.ts diff --git a/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js b/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js index 16fe1dbbfe..93f3e3ae2f 100644 --- a/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js +++ b/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js @@ -38,14 +38,12 @@ class ForceAtlas2BasedRepulsionSolver extends BarnesHutSolver { ); } - const degree = node.edges.length + 1; // the dividing by the distance cubed instead of squared allows us to get the fx and fy components without sines and cosines // it is shorthand for gravityforce with distance squared and fx = dx/distance * gravityForce const gravityForce = (this.options.gravitationalConstant * parentBranch.mass * - node.options.mass * - degree) / + node.options.mass) / Math.pow(distance, 2); const fx = dx * gravityForce; const fy = dy * gravityForce; diff --git a/test/physics/forceatlas2-rotation.test.ts b/test/physics/forceatlas2-rotation.test.ts new file mode 100644 index 0000000000..e88870a615 --- /dev/null +++ b/test/physics/forceatlas2-rotation.test.ts @@ -0,0 +1,169 @@ +import { expect } from "chai"; +import Network from "../../lib/network/Network.js"; +import { canvasMockify } from "../canvas-mock.js"; + +/** + * Regression test for issue #2193: with the `forceAtlas2Based` solver an + * asymmetric graph never settles. After the layout has visually converged the + * whole graph keeps rotating about its center of mass forever. + * + * Root cause: `FA2BasedRepulsionSolver._calculateForces` scaled the pairwise + * repulsion by the *receiving* node's degree only, so the repulsion between two + * nodes of different degree was non-reciprocal (it violated Newton's third law). + * Each tick that injected a small net torque into the system; velocity damping + * caps but never removes it, so the converged layout spins at a terminal + * angular velocity. The bug only manifests on a degree-heterogeneous, + * asymmetric graph -- on a symmetric one the per-pair torques cancel, which is + * why it stayed undiagnosed. + * + * This test builds an asymmetric, strongly degree-heterogeneous graph (a few + * hubs in a cycle, each carrying very different numbers of leaves), converges + * it with a fixed timestep, and then measures the net angular velocity of the + * whole configuration about its center of mass over a window of steps, where + * omega = L / I, L = sum(rx*vy - ry*vx), I = sum(rx^2 + ry^2), and + * r = pos - centerOfMass. + * + * On the buggy code this is a *sustained* spin (|omega| ~= 1.4e-4 -- the window + * mean and window max coincide, i.e. it does not decay). With the fix the + * configuration is genuinely settled and |omega| drops by ~4 orders of + * magnitude (~1.1e-8). The threshold below sits well clear of both. + */ +describe("ForceAtlas2 rotation (issue #2193)", function (): void { + beforeEach(function () { + this.clearJSDOM = canvasMockify("
"); + this.container = document.getElementById("mynetwork"); + }); + + afterEach(function () { + this.clearJSDOM(); + + delete this.clearJSDOM; + delete this.container; + }); + + it("a converged asymmetric layout does not keep rotating", function (): void { + this.timeout(120000); + + // --- Build a deterministic, asymmetric, degree-heterogeneous graph. ------ + // A symmetric graph would NOT exhibit the bug (its per-pair torques cancel), + // so the heterogeneity of the hub degrees below is essential. Initial + // positions come from a seeded PRNG so the run is fully reproducible. + const rng = mulberry32(0xc0ffee); + const nodes: { id: number; x: number; y: number }[] = []; + const edges: { from: number; to: number }[] = []; + + const place = (id: number): void => { + nodes.push({ id, x: (rng() - 0.5) * 600, y: (rng() - 0.5) * 600 }); + }; + + const hubs = [0, 1, 2, 3, 4]; + hubs.forEach((h): void => place(h)); + // Connect the hubs in a cycle. + for (let i = 0; i < hubs.length; i++) { + edges.push({ from: hubs[i], to: hubs[(i + 1) % hubs.length] }); + } + // Give each hub a very different number of leaves -> heterogeneous degrees. + let nextId = hubs.length; + const leavesPerHub = [40, 20, 10, 4, 1]; + hubs.forEach((hub, idx): void => { + for (let k = 0; k < leavesPerHub[idx]; k++) { + const id = nextId++; + place(id); + edges.push({ from: hub, to: id }); + } + }); + + const network = new Network( + this.container, + { nodes, edges }, + { + physics: { + solver: "forceAtlas2Based", + // Constant timestep -> deterministic, comparable angular-velocity + // measurements. Drive the engine ourselves below instead of relying + // on the asynchronous, render-loop-based stabilization. + adaptiveTimestep: false, + stabilization: { enabled: false }, + }, + }, + ); + + const engine = network.physics; + engine.updatePhysicsData(); + + // Converge the layout. + for (let i = 0; i < 6000; i++) { + engine.physicsTick(); + } + + // Measure the net rotation over a window of further steps. Averaging cancels + // the small transient settling motion and isolates *sustained* rotation; + // the per-step maximum is also tracked to confirm the value is steady (a + // terminal spin) rather than a decaying transient. + let omegaSum = 0; + let omegaAbsMax = 0; + const window = 200; + for (let i = 0; i < window; i++) { + engine.physicsTick(); + const omega = angularVelocityAboutCenterOfMass(network); + omegaSum += omega; + omegaAbsMax = Math.max(omegaAbsMax, Math.abs(omega)); + } + const meanOmega = Math.abs(omegaSum / window); + + // Buggy code: meanOmega ~= 1.39e-4 and omegaAbsMax ~= 1.39e-4 (sustained). + // Fixed code: meanOmega ~= 1.14e-8 and omegaAbsMax ~= 1.14e-8 (settled). + // The threshold sits cleanly between the two by orders of magnitude. + expect(meanOmega).to.be.below(1e-5); + expect(omegaAbsMax).to.be.below(1e-5); + }); +}); + +/** + * Net angular velocity of the whole configuration about its center of mass, + * omega = L / I, using live node positions and physics velocities. + * @param network - the vis Network to inspect + * @returns the (signed) angular velocity about the center of mass + */ +function angularVelocityAboutCenterOfMass(network: any): number { + const nodes = network.body.nodes; + const velocities = network.physics.physicsBody.velocities; + const ids = network.physics.physicsBody.physicsNodeIndices; + + let cx = 0; + let cy = 0; + for (const id of ids) { + cx += nodes[id].x; + cy += nodes[id].y; + } + cx /= ids.length; + cy /= ids.length; + + let angularMomentum = 0; // L = sum( rx*vy - ry*vx ) + let momentOfInertia = 0; // I = sum( rx^2 + ry^2 ) + for (const id of ids) { + const rx = nodes[id].x - cx; + const ry = nodes[id].y - cy; + const v = velocities[id]; + angularMomentum += rx * v.y - ry * v.x; + momentOfInertia += rx * rx + ry * ry; + } + return angularMomentum / momentOfInertia; +} + +/** + * Small seeded PRNG (mulberry32) so initial node positions are deterministic + * and the test is fully reproducible. + * @param seed - 32-bit seed + * @returns a function returning the next pseudo-random number in [0, 1) + */ +function mulberry32(seed: number): () => number { + let a = seed >>> 0; + return function (): number { + a |= 0; + a = (a + 0x6d2b79f5) | 0; + let t = Math.imul(a ^ (a >>> 15), 1 | a); + t = (t + Math.imul(t ^ (t >>> 7), 61 | t)) ^ t; + return ((t ^ (t >>> 14)) >>> 0) / 4294967296; + }; +} From a7796e0ad47c0bb8c8c445ab365b9d3dcedbf0ee Mon Sep 17 00:00:00 2001 From: bojanstef <5675392+bojanstef@users.noreply.github.com> Date: Fri, 19 Jun 2026 11:37:48 -0400 Subject: [PATCH 2/3] test(physics): regenerate forceAtlas2Based snapshots for the rotation fix The fix changes forceAtlas2Based layouts (removing the spurious degree factor), so the 6 FA2 pinning snapshots whose graphs have edges change. Values taken from the PR CI run for environment consistency. The 2-disconnected FA2 snapshot is unchanged (degree-1 nodes are unaffected), as are all non-FA2 solver snapshots. --- __snapshots__/test/snapshots.test.ts.js | 196 ++++++++++++------------ 1 file changed, 98 insertions(+), 98 deletions(-) diff --git a/__snapshots__/test/snapshots.test.ts.js b/__snapshots__/test/snapshots.test.ts.js index 2498f17cd0..39633b36e7 100644 --- a/__snapshots__/test/snapshots.test.ts.js +++ b/__snapshots__/test/snapshots.test.ts.js @@ -193,12 +193,12 @@ exports['Physics snapshots Barnes Hut defaults binary tree 1'] = { exports['Physics snapshots Force Atlas 2 defaults 2 connected nodes 1'] = { "1": { - "x": -11, - "y": 58 + "x": -3, + "y": 47 }, "2": { - "x": -8, - "y": -68 + "x": 0, + "y": -53 } } @@ -215,172 +215,172 @@ exports['Physics snapshots Force Atlas 2 defaults 2 disconnected nodes 1'] = { exports['Physics snapshots Force Atlas 2 defaults 2 doubly connected nodes 1'] = { "1": { - "x": -18, - "y": 64 + "x": -13, + "y": 31 }, "2": { - "x": 3, - "y": -65 + "x": 0, + "y": -33 } } exports['Physics snapshots Force Atlas 2 defaults 6 node circle 1'] = { "1": { - "x": -23, - "y": -151 + "x": -9, + "y": -90 }, "2": { - "x": -133, - "y": -48 + "x": -77, + "y": -29 }, "3": { - "x": -111, - "y": 102 + "x": -61, + "y": 60 }, "4": { - "x": 34, - "y": 147 + "x": 21, + "y": 87 }, "5": { - "x": 145, - "y": 45 + "x": 86, + "y": 25 }, "6": { - "x": 120, - "y": -104 + "x": 72, + "y": -61 } } exports['Physics snapshots Force Atlas 2 defaults 6 node complete graph 1'] = { "1": { - "x": -160, - "y": -106 + "x": -60, + "y": -41 }, "2": { - "x": -171, - "y": 94 + "x": -64, + "y": 37 }, "3": { - "x": 175, - "y": 84 + "x": 65, + "y": 37 }, "4": { - "x": 11, - "y": 197 + "x": -1, + "y": 75 }, "5": { - "x": 39, - "y": -151 + "x": 2, + "y": -70 }, "6": { - "x": 110, - "y": -109 + "x": 63, + "y": -32 } } exports['Physics snapshots Force Atlas 2 defaults 6 node star 1'] = { "1": { - "x": 12, - "y": -2 + "x": 11, + "y": 1 }, "2": { - "x": -146, - "y": -12 + "x": -118, + "y": -11 }, "3": { - "x": -40, - "y": -150 + "x": -34, + "y": -117 }, "4": { - "x": 154, - "y": -73 + "x": 126, + "y": -57 }, "5": { - "x": 121, - "y": 116 + "x": 98, + "y": 96 }, "6": { - "x": -48, - "y": 146 + "x": -36, + "y": 121 } } exports['Physics snapshots Force Atlas 2 defaults binary tree 1'] = { "1": { - "x": 52, - "y": 48 - }, - "2": { - "x": -97, - "y": 87 - }, - "3": { - "x": 206, - "y": 34 - }, - "4": { - "x": -113, - "y": -86 - }, - "5": { - "x": -135, - "y": 245 - }, - "6": { - "x": 211, - "y": -132 - }, - "7": { - "x": 199, - "y": 204 - }, - "8": { - "x": -278, - "y": -65 - }, - "9": { - "x": -26, - "y": -228 + "x": 29, + "y": 51 }, "10": { - "x": -279, - "y": 227 + "x": -199, + "y": 151 }, "11": { - "x": -49, - "y": 358 + "x": -35, + "y": 245 }, "12": { - "x": 356, - "y": -104 + "x": 265, + "y": -65 }, "13": { - "x": 213, - "y": -277 + "x": 133, + "y": -186 }, "14": { - "x": 340, - "y": 156 + "x": 252, + "y": 95 }, "15": { - "x": 134, - "y": 333 + "x": 105, + "y": 227 }, "16": { - "x": -304, - "y": -209 + "x": -213, + "y": -147 }, "17": { - "x": -360, - "y": 55 + "x": -248, + "y": 38 }, "18": { - "x": 63, - "y": -344 + "x": 58, + "y": -224 }, "19": { - "x": -138, - "y": -327 + "x": -119, + "y": -230 + }, + "2": { + "x": -48, + "y": 56 + }, + "3": { + "x": 111, + "y": 16 + }, + "4": { + "x": -68, + "y": -55 + }, + "5": { + "x": -99, + "y": 168 + }, + "6": { + "x": 148, + "y": -92 + }, + "7": { + "x": 150, + "y": 130 + }, + "8": { + "x": -190, + "y": -44 + }, + "9": { + "x": -22, + "y": -170 } } From 91b1da40ad2f059b4746a983b2ccd236019950a7 Mon Sep 17 00:00:00 2001 From: Bojan Stefanovic <5675392+bojanstef@users.noreply.github.com> Date: Sat, 20 Jun 2026 12:52:41 -0400 Subject: [PATCH 3/3] fix(physics): apply forceAtlas2 degree weighting symmetrically (#2193) The earlier commits removed degree-based repulsion entirely to stop the perpetual rotation, but degree weighting is core to ForceAtlas2 (it pushes high-degree hubs to the periphery). The actual cause is that the degree factor was applied to only the receiving node, making the pairwise force non-reciprocal (violating Newton's third law) and injecting a net torque each tick that damping caps but never removes. Restore degree weighting in the reciprocal form from the FA2 paper (repulsion proportional to (deg(a)+1)(deg(b)+1)) by folding (degree+1) into the Barnes-Hut effective mass used both to build the tree and to evaluate the force. Hubs are still pushed to the periphery; the force is reciprocal again, so a converged layout no longer rotates. Rewrite the regression test to measure the actual best-fit rotation of node positions over a long, never-stabilizing run; the previous test measured a frozen velocity residual, which does not reflect the visible rotation. Validated across hub, barbell, scale-free and random graphs: net rotation drops 7x-5800x and the layout settles in every case. Regenerate the six affected forceAtlas2Based snapshots accordingly. --- __snapshots__/test/snapshots.test.ts.js | 194 +++++++++--------- .../physics/FA2BasedRepulsionSolver.js | 60 +++++- test/physics/forceatlas2-rotation.test.ts | 136 +++++++----- 3 files changed, 236 insertions(+), 154 deletions(-) diff --git a/__snapshots__/test/snapshots.test.ts.js b/__snapshots__/test/snapshots.test.ts.js index 39633b36e7..0ab84c0b2f 100644 --- a/__snapshots__/test/snapshots.test.ts.js +++ b/__snapshots__/test/snapshots.test.ts.js @@ -193,12 +193,12 @@ exports['Physics snapshots Barnes Hut defaults binary tree 1'] = { exports['Physics snapshots Force Atlas 2 defaults 2 connected nodes 1'] = { "1": { - "x": -3, - "y": 47 + "x": -14, + "y": 64 }, "2": { - "x": 0, - "y": -53 + "x": -10, + "y": -74 } } @@ -215,172 +215,172 @@ exports['Physics snapshots Force Atlas 2 defaults 2 disconnected nodes 1'] = { exports['Physics snapshots Force Atlas 2 defaults 2 doubly connected nodes 1'] = { "1": { - "x": -13, - "y": 31 + "x": -16, + "y": 71 }, "2": { "x": 0, - "y": -33 + "y": -73 } } exports['Physics snapshots Force Atlas 2 defaults 6 node circle 1'] = { "1": { - "x": -9, - "y": -90 + "x": -31, + "y": -191 }, "2": { - "x": -77, - "y": -29 + "x": -173, + "y": -61 }, "3": { - "x": -61, - "y": 60 + "x": -141, + "y": 129 }, "4": { - "x": 21, - "y": 87 + "x": 41, + "y": 187 }, "5": { - "x": 86, - "y": 25 + "x": 184, + "y": 58 }, "6": { - "x": 72, - "y": -61 + "x": 151, + "y": -131 } } exports['Physics snapshots Force Atlas 2 defaults 6 node complete graph 1'] = { "1": { - "x": -60, - "y": -41 + "x": -211, + "y": -144 }, "2": { - "x": -64, - "y": 37 + "x": -234, + "y": 117 }, "3": { - "x": 65, - "y": 37 + "x": 220, + "y": 137 }, "4": { - "x": -1, - "y": 75 + "x": -8, + "y": 264 }, "5": { - "x": 2, - "y": -70 + "x": 22, + "y": -247 }, "6": { - "x": 63, - "y": -32 + "x": 220, + "y": -116 } } exports['Physics snapshots Force Atlas 2 defaults 6 node star 1'] = { "1": { - "x": 11, - "y": 1 + "x": 9, + "y": -1 }, "2": { - "x": -118, - "y": -11 + "x": -196, + "y": -16 }, "3": { - "x": -34, - "y": -117 + "x": -51, + "y": -197 }, "4": { - "x": 126, - "y": -57 + "x": 193, + "y": -97 }, "5": { - "x": 98, - "y": 96 + "x": 155, + "y": 145 }, "6": { - "x": -36, - "y": 121 + "x": -68, + "y": 192 } } exports['Physics snapshots Force Atlas 2 defaults binary tree 1'] = { "1": { - "x": 29, - "y": 51 + "x": 73, + "y": 104 + }, + "2": { + "x": -129, + "y": 123 + }, + "3": { + "x": 269, + "y": 41 + }, + "4": { + "x": -166, + "y": -124 + }, + "5": { + "x": -202, + "y": 345 + }, + "6": { + "x": 327, + "y": -187 + }, + "7": { + "x": 312, + "y": 276 + }, + "8": { + "x": -396, + "y": -106 + }, + "9": { + "x": -61, + "y": -342 }, "10": { - "x": -199, - "y": 151 + "x": -382, + "y": 329 }, "11": { - "x": -35, - "y": 245 + "x": -86, + "y": 488 }, "12": { - "x": 265, - "y": -65 + "x": 501, + "y": -161 }, "13": { - "x": 133, - "y": -186 + "x": 296, + "y": -368 }, "14": { - "x": 252, - "y": 95 + "x": 484, + "y": 221 }, "15": { - "x": 105, - "y": 227 + "x": 215, + "y": 439 }, "16": { - "x": -213, - "y": -147 + "x": -437, + "y": -282 }, "17": { - "x": -248, - "y": 38 + "x": -499, + "y": 49 }, "18": { - "x": 58, - "y": -224 + "x": 81, + "y": -461 }, "19": { - "x": -119, - "y": -230 - }, - "2": { - "x": -48, - "y": 56 - }, - "3": { - "x": 111, - "y": 16 - }, - "4": { - "x": -68, - "y": -55 - }, - "5": { - "x": -99, - "y": 168 - }, - "6": { - "x": 148, - "y": -92 - }, - "7": { - "x": 150, - "y": 130 - }, - "8": { - "x": -190, - "y": -44 - }, - "9": { - "x": -22, - "y": -170 + "x": -182, + "y": -474 } } diff --git a/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js b/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js index 93f3e3ae2f..9aef81d36a 100644 --- a/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js +++ b/lib/network/modules/components/physics/FA2BasedRepulsionSolver.js @@ -16,6 +16,57 @@ class ForceAtlas2BasedRepulsionSolver extends BarnesHutSolver { this._rng = Alea("FORCE ATLAS 2 BASED REPULSION SOLVER"); } + /** + * The repulsion "charge" of a node: its mass weighted by its degree + * (`edges.length + 1`). ForceAtlas2 scales repulsion by node degree so that + * highly interconnected hubs are pushed to the periphery. + * + * The degree must enter the force symmetrically. ForceAtlas2's repulsion is + * `k * (deg(a)+1) * (deg(b)+1) / distance` (Jacomy et al. 2014), i.e. it + * depends on the product of both degrees, so the force on `a` from `b` equals + * the force on `b` from `a` (Newton's third law). Applying the degree to only + * the receiving node — as was done previously — makes the pairwise force + * non-reciprocal and injects a net torque every tick. Velocity damping caps + * but never removes it, so a converged layout spins about its center of mass + * forever (issue #2193). Folding the degree into the mass used both to build + * the Barnes-Hut tree and to evaluate the force keeps it reciprocal. + * @param {Node} node + * @returns {number} the degree-weighted mass of the node + * @private + */ + _getMass(node) { + return node.options.mass * (node.edges.length + 1); + } + + /** + * Update the mass of a branch using the degree-weighted mass (see _getMass), + * so the accumulated branch mass and its center of mass are consistent with + * the degree-weighted repulsion applied in _calculateForces. + * @param {object} parentBranch + * @param {Node} node + * @private + */ + _updateBranchMass(parentBranch, node) { + const mass = this._getMass(node); + const centerOfMass = parentBranch.centerOfMass; + const totalMass = parentBranch.mass + mass; + const totalMassInv = 1 / totalMass; + + centerOfMass.x = centerOfMass.x * parentBranch.mass + node.x * mass; + centerOfMass.x *= totalMassInv; + + centerOfMass.y = centerOfMass.y * parentBranch.mass + node.y * mass; + centerOfMass.y *= totalMassInv; + + parentBranch.mass = totalMass; + const biggestSize = Math.max( + Math.max(node.height, node.radius), + node.width, + ); + parentBranch.maxWidth = + parentBranch.maxWidth < biggestSize ? biggestSize : parentBranch.maxWidth; + } + /** * Calculate the forces based on the distance. * @param {number} distance @@ -38,12 +89,15 @@ class ForceAtlas2BasedRepulsionSolver extends BarnesHutSolver { ); } - // the dividing by the distance cubed instead of squared allows us to get the fx and fy components without sines and cosines - // it is shorthand for gravityforce with distance squared and fx = dx/distance * gravityForce + // dividing by the distance squared (rather than the cubed used by + // BarnesHutSolver) yields ForceAtlas2's linear 1/distance repulsion once the + // dx/dy components are factored in. Both parentBranch.mass and the node's + // mass are degree-weighted (see _getMass), so the repulsion is proportional + // to (degree+1) on both sides and stays reciprocal. const gravityForce = (this.options.gravitationalConstant * parentBranch.mass * - node.options.mass) / + this._getMass(node)) / Math.pow(distance, 2); const fx = dx * gravityForce; const fy = dy * gravityForce; diff --git a/test/physics/forceatlas2-rotation.test.ts b/test/physics/forceatlas2-rotation.test.ts index e88870a615..01e4c664ae 100644 --- a/test/physics/forceatlas2-rotation.test.ts +++ b/test/physics/forceatlas2-rotation.test.ts @@ -4,29 +4,36 @@ import { canvasMockify } from "../canvas-mock.js"; /** * Regression test for issue #2193: with the `forceAtlas2Based` solver an - * asymmetric graph never settles. After the layout has visually converged the - * whole graph keeps rotating about its center of mass forever. + * asymmetric graph converges to a stable *shape* but the whole layout then + * keeps rotating about its center of mass forever. * - * Root cause: `FA2BasedRepulsionSolver._calculateForces` scaled the pairwise - * repulsion by the *receiving* node's degree only, so the repulsion between two - * nodes of different degree was non-reciprocal (it violated Newton's third law). - * Each tick that injected a small net torque into the system; velocity damping - * caps but never removes it, so the converged layout spins at a terminal - * angular velocity. The bug only manifests on a degree-heterogeneous, - * asymmetric graph -- on a symmetric one the per-pair torques cancel, which is - * why it stayed undiagnosed. + * Root cause: `FA2BasedRepulsionSolver` scaled the pairwise repulsion by the + * *receiving* node's degree only, so the repulsion between two nodes of + * different degree was non-reciprocal -- it violated Newton's third law. + * ForceAtlas2's repulsion is defined as `k * (deg(a)+1) * (deg(b)+1) / dist` + * (Jacomy et al. 2014): symmetric in the two degrees, hence reciprocal. The + * one-sided factor injects a small net torque every tick; velocity damping + * removes relative motion but not a rigid rotation, so the converged layout + * spins at a terminal rate. The fix folds `(degree+1)` into the Barnes-Hut + * effective mass (both when building the tree and when evaluating the force), + * which keeps the degree-weighted "hubs to the periphery" behaviour while + * making the force reciprocal again. The bug only shows on a + * degree-heterogeneous, asymmetric graph -- on a symmetric one the per-pair + * torques cancel, which is why it stayed undiagnosed. * - * This test builds an asymmetric, strongly degree-heterogeneous graph (a few - * hubs in a cycle, each carrying very different numbers of leaves), converges - * it with a fixed timestep, and then measures the net angular velocity of the - * whole configuration about its center of mass over a window of steps, where - * omega = L / I, L = sum(rx*vy - ry*vx), I = sum(rx^2 + ry^2), and - * r = pos - centerOfMass. + * What this test measures: the *actual rotation of the node positions*, not a + * velocity proxy. It drives the engine with `minVelocity: 0` so the solver + * never reports "stabilized" (the regime in which the live network keeps + * ticking and the rotation is visible), converges the shape, then measures the + * best-fit rigid rotation angle of the whole configuration over a long window + * via the orthogonal Procrustes formula + * angle = atan2( Σ(x0*y1 - y0*x1), Σ(x0*x1 + y0*y1) ) + * where (x0,y0) and (x1,y1) are node positions relative to the centroid at the + * start and end of the window. * - * On the buggy code this is a *sustained* spin (|omega| ~= 1.4e-4 -- the window - * mean and window max coincide, i.e. it does not decay). With the fix the - * configuration is genuinely settled and |omega| drops by ~4 orders of - * magnitude (~1.1e-8). The threshold below sits well clear of both. + * On the buggy code the layout rotates ~122 deg over the 20000-step window and + * never settles; with the fix it rotates < 0.5 deg and settles. The threshold + * below (5 deg) sits more than an order of magnitude clear of both. */ describe("ForceAtlas2 rotation (issue #2193)", function (): void { beforeEach(function () { @@ -79,11 +86,15 @@ describe("ForceAtlas2 rotation (issue #2193)", function (): void { { physics: { solver: "forceAtlas2Based", - // Constant timestep -> deterministic, comparable angular-velocity - // measurements. Drive the engine ourselves below instead of relying - // on the asynchronous, render-loop-based stabilization. + // Constant timestep -> deterministic, comparable measurements. Drive + // the engine ourselves below instead of relying on the asynchronous, + // render-loop-based stabilization. adaptiveTimestep: false, stabilization: { enabled: false }, + // Never report "stabilized" so physicsTick keeps integrating; this is + // the regime in which the perpetual rotation is observable (a settled + // engine would simply freeze and hide it). + minVelocity: 0, }, }, ); @@ -91,43 +102,42 @@ describe("ForceAtlas2 rotation (issue #2193)", function (): void { const engine = network.physics; engine.updatePhysicsData(); - // Converge the layout. + // Converge the shape. for (let i = 0; i < 6000; i++) { engine.physicsTick(); } - // Measure the net rotation over a window of further steps. Averaging cancels - // the small transient settling motion and isolates *sustained* rotation; - // the per-step maximum is also tracked to confirm the value is steady (a - // terminal spin) rather than a decaying transient. - let omegaSum = 0; - let omegaAbsMax = 0; - const window = 200; - for (let i = 0; i < window; i++) { + // Snapshot positions relative to the centroid, advance a long window, then + // measure the net rigid rotation between the two snapshots. + const before = relativePositions(network); + const windowSteps = 20000; + for (let i = 0; i < windowSteps; i++) { engine.physicsTick(); - const omega = angularVelocityAboutCenterOfMass(network); - omegaSum += omega; - omegaAbsMax = Math.max(omegaAbsMax, Math.abs(omega)); } - const meanOmega = Math.abs(omegaSum / window); + const after = relativePositions(network); - // Buggy code: meanOmega ~= 1.39e-4 and omegaAbsMax ~= 1.39e-4 (sustained). - // Fixed code: meanOmega ~= 1.14e-8 and omegaAbsMax ~= 1.14e-8 (settled). - // The threshold sits cleanly between the two by orders of magnitude. - expect(meanOmega).to.be.below(1e-5); - expect(omegaAbsMax).to.be.below(1e-5); + const rotationDeg = Math.abs( + (bestFitRotation(before, after) * 180) / Math.PI, + ); + + // Stop the render loop before the JSDOM is torn down (the canvas mock + // shims requestAnimationFrame onto setTimeout, which would otherwise fire + // after teardown and throw "window is not defined"). + network.destroy(); + + // Buggy code: ~122 deg (and still climbing -- it never settles). + // Fixed code: < 0.5 deg. The threshold sits well clear of both. + expect(rotationDeg).to.be.below(5); }); }); /** - * Net angular velocity of the whole configuration about its center of mass, - * omega = L / I, using live node positions and physics velocities. + * Node positions relative to the configuration's centroid. * @param network - the vis Network to inspect - * @returns the (signed) angular velocity about the center of mass + * @returns parallel arrays of x and y offsets from the centroid */ -function angularVelocityAboutCenterOfMass(network: any): number { +function relativePositions(network: any): { x: number[]; y: number[] } { const nodes = network.body.nodes; - const velocities = network.physics.physicsBody.velocities; const ids = network.physics.physicsBody.physicsNodeIndices; let cx = 0; @@ -139,16 +149,34 @@ function angularVelocityAboutCenterOfMass(network: any): number { cx /= ids.length; cy /= ids.length; - let angularMomentum = 0; // L = sum( rx*vy - ry*vx ) - let momentOfInertia = 0; // I = sum( rx^2 + ry^2 ) + const x: number[] = []; + const y: number[] = []; for (const id of ids) { - const rx = nodes[id].x - cx; - const ry = nodes[id].y - cy; - const v = velocities[id]; - angularMomentum += rx * v.y - ry * v.x; - momentOfInertia += rx * rx + ry * ry; + x.push(nodes[id].x - cx); + y.push(nodes[id].y - cy); + } + return { x, y }; +} + +/** + * Best-fit rigid rotation angle taking configuration `a` onto configuration `b` + * (orthogonal Procrustes, 2D): angle = atan2(Σ a×b, Σ a·b). Both configurations + * must be centered (see relativePositions) and index-aligned. + * @param a - start positions relative to centroid + * @param b - end positions relative to centroid + * @returns the signed best-fit rotation angle in radians + */ +function bestFitRotation( + a: { x: number[]; y: number[] }, + b: { x: number[]; y: number[] }, +): number { + let cross = 0; // Σ (ax*by - ay*bx) + let dot = 0; // Σ (ax*bx + ay*by) + for (let i = 0; i < a.x.length; i++) { + cross += a.x[i] * b.y[i] - a.y[i] * b.x[i]; + dot += a.x[i] * b.x[i] + a.y[i] * b.y[i]; } - return angularMomentum / momentOfInertia; + return Math.atan2(cross, dot); } /**