Force directed graph with layered gravity and physically modeled collisions

Force Directed Graph with Custom Gravity

Features

Metrics display and inputs
The metrics panel across the top of the svg element gives a live display of the layout state. The inputs on the left allow for the number of nodes, the maximum recovery velocity for escaped nodes and the windup factor for clearing overlaps to be adjusted live. The current alpha value for the layout is displayed along with instantaneous and averaged tick time and the average frame rate of the layout. Changing the number of nodes re-starts the layout.
Custom gravity (layers)
The nodes have a cy member which is its target y value, the gravity function, which is called from the force.on("tick"...) callback has behaviour to move the nodes toward their cy position on each tick.
Collisions between nodes Based on this example but enhanced to accurately model relative "mass" the rebound velocities are solved for using the equations for conservation of momentum and conservation of energy. The "efficiency" of the collisions - as defined by the ratio of the actual final velocities divided by the ideal, perfectly elastic system - is set internally.
Boundary constraints and collisions
This behaviour is also included in the gravity function and uses basic geometry to reflect the incident velocity in the plane of the boundary. It does this by using the Node velocity API to manipulate the d.px and d.py values of the nodes. It is possible, however, for nodes to penetrate the boundaries due to limitations in the temporal resolution of the layout.
Recovering escaped nodes
If a node escapes the boundaries of the bounding box, the velocity is still reflected in the plane of the penetrated boundary. If the escaped node has no velocity (p.q.x/y == p.q.px/py) then the node is steered back toward the boundary, with constant speed, or a random speed (between 0 and recover speed) if it has come to rest. After it is fully recovered, the velocity of the node will be naturally determined by the distance it was last moved by the recovery behaviour.
Node velocity API
Behaviour is added to the nodes to allow easy setting and getting of node velocity vector. This is accomplished by manipulating the previous value of the node position vector.
Node velocity indicators
Velocity vector tracers are toggled by a button at the bottom of the graph. A line is simply drawn from the present position to the present position minus the velocity, times the magnification factor, to indicate the node velocity.
Quantisation of screen position
It is possible for nodes with different values to be stationary due to quantisation of d.x and d.px by the rendering process. This means that decisions in the code based on position may not reflect the rendered state properly. This is fixed by adding a getter to each data point that returns a rounded version of d.x/y and d.px/y: d.q.x/y and d.q.px/y. The quantised versions are used for the decision making in the boundary collisions and recovery behaviour.
Guaranteed un-mixing of nodes
It's possible for nodes to be blocked from reaching they're cy positions by other, more massive nodes in adjacent layers, in order to overcome this, a frustration factor is applied to the isolated nodes which has the effect of increasing they're effective mass with each tick that they are outside they're designated band. The increase in mass factor per tick is the windup which is internally set to 0.01. As soon as the nodes make it to they're layer they're anxiety is switched back to 1.
Clearing node overlaps
Smaller nodes can be squashed in the matrix of bigger nodes, because they have mass is insignificant, the matrix of bigger nodes is not moved and the smallest nodes remain overlapped. In order to clear this issue, the overlaps are monitored and the effective mass of the smaller node in the overlapping pair is increased by fear factor, every collision where the overlap increases. The increase per collision is set by the windUp control in the metrics panel. Nodes with a fear greater than one are highlighted with a lighter stroke color. If the overlap decreases, the fear factor is also decreased by the wind-up value.
Extended cooling time
The cooling time for the force layout is determined by the evolution of its alpha value. This is set to 0.1 when the layout is started and is normally reduced by 1% each tick until it reaches 0.005, at which point the layout will stop updating. It is possible to manipulate the cooling rate so that it cools at x% by dividing the current alpha value by the standard factor of 0.99 and multiplying it by (1 - x%) at the end of the tick callback:

  force.alpha(a/0.99*(1 - x))

d3 features used

<!DOCTYPE html>
<meta charset="utf-8">
<link rel="stylesheet" type="text/css"
      href="https://gitcdn.xyz/repo/cool-Blue/d3-lib/1dcef5868b3e774e8303d0c6eeaaef2528867893/inputs/button/style.css">
<style>
    body {
        /*margin: 200px 500px 100px 500px;*/
    }

#inputs {
        display: inline-block;
        margin: 0 0 0 100px;
        border: none;
        padding: 0 0 0 1em;
        box-sizing: border-box;
        background-color: black;
    }

#metrics {
        display: inline-block;
    }

label, input {
        text-align: left;
        width: 3.5em;
        color: orange;
        /*padding-left: 1em;*/
        background-color: black;
        outline: none;
        border: none;
    }

circle {
        stroke: black;
    }

svg {
        display: block;
        overflow: visible;
        border: none;
        background: black;
        margin: 0 100px 0 100px;
    }

text {
        text-anchor: middle;
    }

rect {
        stroke: #ccc;
    }
    .g-button {
        color: #804700;
        background: black;
        border-color: orange;
    }
    .g-button.g-active {
        color: orange;
        background: #333333;
        border-color: orange;
    }

</style>

width = 960 - 200,
        height = 500 - elapsedTime.selection.node().clientHeight - buttons.height,
        r0 = 5.5,

n = 2000, // total number of nodes
        m = 10; // number of distinct layers

nodeCount
        .property("value", n)
        .on("change", function() {
            viz = update(force, this.value);
            this.blur();
        });
    reEntrySpeed
        .property("value", 5)
        .value = function() { return this.property("value")};
    windUp
        .property("value", 1)
        .value = function() { return +this.property("value")};

elapsedTime.selection.style({
        width: (width - parseFloat(window.getComputedStyle(d3.select("#inputs").node()).getPropertyValue("width"))) + "px"
    });

var color = d3.scale.category10()
        .domain(d3.range(m)),

y = d3.scale.ordinal()
            .domain(d3.range(m))
            .rangePoints([height, 0], 1),
        w = d3.scale.ordinal()
            .domain(d3.range(m))
            .rangeBands([height, 0]),
        wRange = w.range(),

bubble = Bubble(svg);

svg.attr("width", width)
        .attr("height", height)
        .append("g");
    svg.selectAll("wells").data(wRange).enter().append("rect")
        .attr({width: width, height: w.rangeBand(), y: function(d) {return d}});

var force = d3.layout.force()
        .size([width, height])
        .gravity(0)
        .charge(0)
        .friction(0.5)
        .on("tick", tick)
        .on("start", function() {
            elapsedTime.start(1000);
        })
        .on("end", function() {
            window.requestAnimationFrame(tick)
        }),

viz = update(force, n);

// Move nodes toward cluster focus.
    function gravity(alpha) {
        var moreThan;
        return function(d) {
            //reflect off the edges of the container
            // check for boundary collisions and reverse velocity if necessary
            if((moreThan = d.x >= (width - d.radius)) || d.x <= d.radius) {
                // if the object is outside the boundaries
                //   manage the sign of its x velocity component to ensure it is moving back into the bounds
                if(~~d.v.x) d.v.x *= moreThan && d.v.x > 0 || !moreThan && d.v.x < 0 ? -1 : 1;
                //   if vx is too small, then steer it back in
                else d.sx = (~~Math.abs(d.v.y) || Math.random() * reEntrySpeed.value()) * (moreThan ? -1 : 1)

}
            if((moreThan = d.y >= (height - d.radius)) || d.y <= d.radius) {
                if(~~d.v.y) d.v.y *= moreThan && d.v.y > 0 || !moreThan && d.v.y < 0 ? -1 : 1;
                else d.sy = (~~Math.abs(d.v.x) || Math.random() * reEntrySpeed.value()) * (moreThan ? -1 : 1)
            }

//find the layers

d.y += d.fixed ? 0 : (d.cy - d.y) /** d.frustration()*/ * alpha;

};
    }

// collision detection
    // physically accurate: conservation of energy and momentum
    function Collide(data, eff) {
        var maxRadius = d3.max(data, function(d) {
            return d.radius
        });
        return function collide(alpha) {
            var quadtree = Math.random() > 0 && d3.geom.quadtree(data) || null;
            return quadtree ? function(d) {
                var d_r  = d.radius,
                    d_x  = d.x,
                    d_y  = d.y,
                    cell = d_r + maxRadius,
                    nx1  = d_x - cell,
                    nx2  = d_x + cell,
                    ny1  = d_y - cell,
                    ny2  = d_y + cell;
                quadtree.visit(function(quad, x1, y1, x2, y2) {
                    var possible = !(x1 > nx2 || x2 < nx1 || y1 > ny2 || y2 < ny1);
                    var q = quad.point;
                    if(q && (q !== d) && possible) {
                        var x = d_x - q.x,
                            y = d_y - q.y,
                            l = Math.sqrt(x * x + y * y),
                            r = d_r + q.radius,
                            margin = l - r;
                        if(margin < 0) {
                            var m    = d.m,
                                mq   = q.m,
                                quadIsBigger = mq > m,
                                mT   = m + mq,
                                mr   = m / mT * d.frustration() * d.anxiety(margin, quadIsBigger),
                                mqr  = mq / mT * q.frustration() * q.anxiety(margin, !quadIsBigger),

vels = bounce(d, q, x, y, m, mq),
                                _l   = (margin) / l * (1 + alpha);

d.x -= (x *= _l) * mqr;
                            d.y -= (y *= _l) * mqr;
                            q.x += x * mr;
                            q.y += y * mr;

if(false && (d.s > 5 || q.s > 5)) {
                                var fmt = " >8,.3f";
                                console.log("hit: " + [f(fmt, (l - r)), f(fmt, _l)]
                                    + "\td: " + f(fmt, [d.s, vels.d.x, vels.d.y])
                                    + "\tq: " + f(fmt, [q.s, vels.q.x, vels.q.y]));
                            }

d.v = vels.d;
                            q.v = vels.q;

}
                    }
                    return !possible;
                });
                function bounce(d, q, x, y, m, mq) {
                    // Note: confirmed that lookup tables for sin and cos are an order of magnitude slower
                    var collision_angle             = Math.atan2(y, x),
                        magnitude_d                 = Math.sqrt(d.v.x * d.v.x + d.v.y * d.v.y),
                        magnitude_q                 = Math.sqrt(q.v.x * q.v.x + q.v.y * q.v.y),
                        direction_d                 = Math.atan2(d.v.y, d.v.x),
                        direction_q                 = Math.atan2(q.v.y, q.v.x),
                        new_vx_d                    = magnitude_d * Math.cos(direction_d - collision_angle) * eff,
                        new_vy_d                    = magnitude_d * Math.sin(direction_d - collision_angle) * eff,
                        new_vx_q                    = magnitude_q * Math.cos(direction_q - collision_angle) * eff,
                        new_vy_q                    = magnitude_q * Math.sin(direction_q - collision_angle) * eff,
                        final_vx_d                  = ((m - mq) * new_vx_d + (mq + mq) * new_vx_q) / (m + mq),
                        final_vx_q                  = ((m + m) * new_vx_d + (mq - m) * new_vx_q) / (m + mq),
                        final_vy_d                  = new_vy_d,
                        final_vy_q                  = new_vy_q,
                        cos_collision_angle         = Math.cos(collision_angle),
                        sin_collision_angle         = Math.sin(collision_angle),
                        cos_collision_angle_plus_90 = Math.cos(collision_angle + Math.PI / 2),
                        sin_collision_angle_plus_90 = Math.sin(collision_angle + Math.PI / 2);

d.v = {
                        x: cos_collision_angle * final_vx_d + cos_collision_angle_plus_90 * final_vy_d,
                        y: sin_collision_angle * final_vx_d + sin_collision_angle_plus_90 * final_vy_d
                    };
                    q.v = {
                        x: cos_collision_angle * final_vx_q + cos_collision_angle_plus_90 * final_vy_q,
                        y: sin_collision_angle * final_vx_q + sin_collision_angle_plus_90 * final_vy_q
                    };
                    return {d: d.v, q: q.v}
                }
            } : function() {
            };
        }
    }

/*    function Collide(nodes, padding) {
     // Resolve collisions between nodes.
     var r0 = d3.max(nodes, function(d) {return d.radius});
     return function collide(alpha) {
     var quadtree = d3.geom.quadtree(nodes);
     return function(d) {
     var r = d.radius + r0 + padding,
     nx1 = d.x - r,
     nx2 = d.x + r,
     ny1 = d.y - r,
     ny2 = d.y + r;
     quadtree.visit(function(quad, x1, y1, x2, y2) {
     var possible = !(x1 > nx2 || x2 < nx1 || y1 > ny2 || y2 < ny1);
     if (quad.point && (quad.point !== d) && possible) {
     var x = d.x - quad.point.x,
     y = d.y - quad.point.y,
     l = Math.sqrt(x * x + y * y),
     r = d.radius + quad.point.radius + padding,
     m = Math.pow(quad.point.radius, 3),
     mq = Math.pow(d.radius, 3),
     mT = m + mq;
     if (l < r) {
     //move the nodes away from each other along the radial (normal) vector
     //taking relative mass into consideration, the sign is already established
     //in calculating x and y and the nodes are modelled as spheres for calculating mass
     l = (r - l) / l * (1 + alpha);
     d.x += (x *= l) * m/mT;
     d.y += (y *= l) * m/mT;
     quad.point.x -= x * mq/mT;
     quad.point.y -= y * mq/mT;
     }
     }
     return !possible;
     });
     };
     }
     }*/
    function initNodes(force, n) {
        force.nodes(d3.range(n).map(function() {
            var layer  = Math.floor(Math.random() * m),
                //      v = (layer + 1) / m * -Math.log(Math.random());
                v      = -Math.log(Math.random()),
                radius = Math.sqrt(v) * r0;
            return {
                radius: radius,
                m: Math.pow(radius, 3),
                color: layer,
                cy: y(layer),
                get v() {
                    var d = this;
                    return {x: d.x - d.px || 0, y: d.y - d.py || 0}
                },
                set v(v) {
                    var d = this;
                    d.px = d.x - v.x;
                    d.py = d.y - v.y;
                },
                set sx(s) {
                    this.v = {x: s, y: this.v.y}
                },
                set sy(s) {
                    this.v = {y: s, x: this.v.x}
                },
                get s() {
                    var v = this.v;
                    return Math.sqrt(v.x * v.x + v.y * v.y)
                },
                frustration: (function() {
                    //if they can't get home, they get angry, but, as soon as they're home, they're fine
                    var anger = 1, windUp = 0.01;
                    return function() {
                        // adjust frustration level based on context and windup rate
                        var d = this, anxious = (Math.abs(d.cy - d.y) > w.rangeBand()
                        / 2);
                        return anger = anxious ? anger + windUp : 1;
                    }
                })(),
                anxiety: (function() {
                    // get agitated if overlaps keep increasing
                    var fear = 1, pervOverlap;
                    return function(overlap, runt) {
                        if(typeof overlap == "undefined") return fear;
                        // adjust anxiety level based on context and windup rate
                        var afraid = ((overlap = (-overlap)) > pervOverlap);
                        pervOverlap = overlap;
                        return fear += afraid && runt ? windUp.value() : fear - windUp.value() < 1 ? fear -1 : -windUp.value();
                    }
                })()
            };
        }));
        //      var collide = Collide(force.nodes(), padding);
        force.start();
        // add a quantiser object that returns a quantised version of all numerical properties
        force.nodes().forEach(function(d) {
            d.q = {};
            Object.keys(d).forEach(function(p) {
                if(!isNaN(d[p])) Object.defineProperty(d.q, p, {
                    get: function() {return Math.round(d[p])}
                });
            })
        });
        return Collide(force.nodes(), .8);
    }
    function update(force, n) {
        return {
            Collide: initNodes(force, n),
            circle: (function() {
                var update = svg.selectAll("circle")
                    .data(force.nodes());
                update.enter().append("circle");
                update.exit().remove();
                update.attr("r", function(d) {
                    return d.radius;
                })
                    //            .style("fill", function (d) {
                    //              return d.color;
                    //            })
                    .call(bubble.call)
                    .call(force.drag)
                return update;
            })(),
            velocities: function() {
                var update = svg.selectAll("line")
                    .data(buttons.showVectors ? force.nodes() : []);
                update.enter().append("line");
                update.exit().remove();
                update
                    .attr("x1", function(d) {
                        return d.px;
                    })
                    .attr("y1", function(d) {
                        return d.py;
                    })
                    .attr("x2", function(d) {
                        return d.x;
                    })
                    .attr("y2", function(d) {
                        return d.y;
                    })
                    .style({"stroke": "red",
                        "stroke-width": "2px",
                        "stroke-linecap": "round",
                        opacity: 0.5
                    });
                return update;
            }
        };
    }
    function Bubble(svg) {
        var colors = d3.range(20).map(d3.scale.category10()).map(function(d) {
            return filters.sphere(svg, d, 1)
        });
        return {
            call: function(selection) {
                selection.style("fill", function(d) {
                    return colors[d.color]
                })
            },
            map: function(d, i, data) {
                d.fill = colors[~~(Math.random() * 20)];
            },
            fill: function(d) {return d.fill}
        }
    };
    function f(_fmt, x) {
        return Array.isArray(x) ? x.map(f.bind(null, _fmt)) : d3.format(_fmt)(x);
    }

</script>
</body>

Updated missing url https://gitcdn.xyz/repo/cool-Blue/d3-lib/master/filters/shadow.js to https://cdn.jsdelivr.net/gh/repo/cool-blue/d3-lib/filters/shadow.js

Updated missing url https://gitcdn.xyz/repo/cool-Blue/d3-lib/master/elapsedTime/elapsed-time-1.0.js to https://cdn.jsdelivr.net/gh/repo/cool-blue/d3-lib/elapsedtime/elapsed-time-1.0.js

Updated missing url https://gitcdn.xyz/repo/cool-Blue/d3-lib/master/inputs/button/2.0.0/button.js to https://cdn.jsdelivr.net/gh/repo/cool-blue/d3-lib/inputs/button/2.0.0/button.js

https://cdnjs.cloudflare.com/ajax/libs/d3/3.5.5/d3.min.js

https://cdnjs.cloudflare.com/ajax/libs/pixi.js/3.0.7/pixi.js

https://cdnjs.cloudflare.com/ajax/libs/tinycolor/1.1.2/tinycolor.min.js

https://gitcdn.xyz/repo/cool-Blue/d3-lib/master/filters/shadow.js

https://gitcdn.xyz/repo/cool-Blue/d3-lib/master/elapsedTime/elapsed-time-1.0.js

https://gitcdn.xyz/repo/cool-Blue/d3-lib/master/inputs/button/2.0.0/button.js