The Barnes-Hut Approximation The Barnes-Hut approximation is a method for calculating the forces exerted by a collection of interacting particles (also known as the "n-body problem"), for example due to gravity or electromagnetism. The method is useful for physical simulation, including force-directed graph layout. Instead of a costly comparison of all individual particles, the approximation uses a quadtree for efficient comparison with centers of charge. Click or drag on the plot to probe the force exerted by a set of "charged" points and see which centers of charge are used for the force calculation. Adjust the slider to change the level of approximation (the Barnes-Hut "theta" parameter); a value of zero corresponds to naïve comparison with all points. For more details, see the interactive explanation of the Barnes-Hut approximation.
chart = { const svg = d3.select(DOM.svg(width, height)); // visualize quadtree boxes svg.append('g') .selectAll('rect').data(quadBoxes(quad)) .enter().append('rect') .attr('x', q => q[0][0] + 0.5) .attr('y', q => q[0][1] + 0.5) .attr('width', q => q[1][0] - q[0][0]) .attr('height', q => q[1][1] - q[0][1]) .style('fill', 'none') .style('stroke', '#ddd') .style('line-width', 0.5); // add Barnes-Hut annotation layer svg.append('g') .attr('class', 'barneshut'); // visualize input points svg.append('g') .selectAll('circle') .data(points) .enter().append('circle') .attr('cx', d => d[0]) .attr('cy', d => d[1]) .attr('r', pointRadius) .style('fill', '#666') .style('fill-opacity', 0.5); // drag interaction indicator svg.append("path") .attr("class", "path-move") .attr("transform", `translate(25,${height - 22})scale(0.05)`) .attr("d", dragPath) .append("title") .text("Click or drag within the plot!"); return svg.node(); }
viewof theta = DOM.range(0, 2, 0.1)
`theta = ${theta}`
probe = Generators.observe(notify => { // set initial probe coordinates notify({x: width / 2, y: height / 2}); // update probe coordinate on mousedown or drag const onDown = function(event) { notify(domPoint(chart, event)); }; const onMove = function(event) { if (event.buttons) notify(domPoint(chart, event)); }; chart.addEventListener("mousedown", onDown); chart.addEventListener("mousemove", onMove); return function() { chart.removeEventListener("mousedown", onDown); chart.removeEventListener("mousemove", onMove); }; })
estimate = { let g = d3.select(chart).select('.barneshut').html(''), fg = g.append('g'), lg = g.append('g'), cg = g.append('g'); const p = probe; if (!p) return; const theta2 = theta * theta, charges = [], boxes = []; let fx = 0, fy = 0; quad.visit(function(quad, x1, y1, x2, y2) { if (!quad.value) return true; let x = quad.x - p.x, y = quad.y - p.y, w = x2 - x1, l = x * x + y * y; // Apply the Barnes-Hut approximation. if (quad.length && w * w / theta2 < l) { let c = { x: quad.x, y: quad.y, v: quad.value, s: 5e3 * quad.value / l }; charges.push(c); boxes.push({x: x1, y: y1, w: w, h: y2 - y1}); fx += x * quad.value / l; fy += y * quad.value / l; return true; } // Otherwise, process points directly. else if (quad.length || !l) return; do if (quad.data !== p) { charges.push({ x: quad.data[0], y: quad.data[1], v: 1, s: 5e3 / l }); fx += x / l; fy += y / l; } while (quad = quad.next); }); fg.selectAll('rect').data(boxes) .enter().append('rect') .attr('x', d => d.x) .attr('y', d => d.y) .attr('width', d => d.w) .attr('height', d => d.h) .style('pointer-events', 'none') .style('fill', 'none') .style('stroke', 'd5aaaa') .style('stroke-width', 2); lg.selectAll('path').data(charges) .enter().append('path') .style('pointer-events', 'none') .style('stroke', '#991151') .style('stroke-opacity', 0.3) .style('stroke-dasharray', [5, 5]) .style('stroke-linecap', 'round') .style('stroke-width', d => Math.max(1, Math.min(5, d.s || 1))) .attr('d', d => 'M' + d.x + ',' + d.y + 'L' + p.x + ',' + p.y); cg.selectAll('circle').data(charges) .enter().append('circle') .attr('cx', d => d.x) .attr('cy', d => d.y) .attr('r', d => pointRadius * (Math.sqrt(d.v) || 1)) .style('pointer-events', 'none') .style('fill', 'firebrick'); cg.append('circle') .attr('cx', p.x) .attr('cy', p.y) .attr('r', pointRadius) .style('pointer-events', 'none') .style('fill', 'white') .style('stroke-width', 1.5) .style('stroke-linecap', 'round') .style('stroke', '#800080'); fx *= forceScalar; fy *= forceScalar; cg.append('path') .attr('d', 'M' + p.x + ',' + p.y + 'L' + (p.x - fx) + ',' + (p.y - fy)) .style('pointer-events', 'none') .style('fill', 'none') .style('stroke-width', 1.5) .style('stroke-linecap', 'round') .style('stroke', 'purple'); return {fx: fx, fy: fy}; }
height = 514
pointRadius = 4
forceScalar = 90
quad = d3.quadtree() .extent([[1, 1], [width - 1, height - 1]]) .addAll(points) .visitAfter(accumulate)
function accumulate(quad, x1, y1, x2, y2) { let strength = 0, q, c, x, y, i; // For internal nodes, accumulate forces from child quadrants. if (quad.length) { let pid = [x1,y1,x2,y2].join(','); for (x = y = i = 0; i < 4; ++i) { if ((q = quad[i]) && (c = q.value)) { strength += c, x += c * q.x, y += c * q.y; q.parent = quad; q.pid = pid; } } quad.x = x / strength; quad.y = y / strength; quad.x1 = x1; quad.y1 = y1; quad.w = x2 - x1; quad.h = y2 - y1; } // For leaf nodes, accumulate forces from coincident quadrants. else { q = quad; q.x = q.data[0]; q.y = q.data[1]; do strength += 1; while (q = q.next); } quad.value = strength; }
function quadBoxes(quad) { let boxes = []; function processNode(node, extent, depth) { if (Array.isArray(node)) { processSplit(node, extent, depth); } } function processSplit(node, extent, depth) { let lo = extent[0], hi = extent[1], mp = [(lo[0] + hi[0]) >> 1, (lo[1] + hi[1]) >> 1]; let e = [ [lo, mp], [[mp[0], lo[1]], [hi[0], mp[1]]], [[lo[0], mp[1]], [mp[0], hi[1]]], [mp, hi] ]; for (let i=0; i<4; ++i) { boxes.push(e[i]); e[i].depth = depth; if (node[i]) processNode(node[i], e[i], depth + 1); } } if (quad.root()) { // add quadtree root extents boxes.push(quad.extent().slice()); // recurse to process quadtree content processNode(quad.root(), quad.extent(), 1); } return boxes; }
function toPoint(xy) { return {x: xy[0], y: xy[1]}; }
points = [[292,318], [248,330], [348,290], [354,311], [251,310], [270,353], [256,337], [255,299], [259,348], [247,320], [443,314], [440,276], [412,196], [428,312], [404,287], [403,272], [451,80], [398,50], [420,15], [419,43], [400,18], [437,28], [385,34], [417,123], [462,206], [500,237], [478,188], [458,248], [389,243], [399,299], [385,139], [418,188], [415,306], [417,249], [421,362], [382,340], [406,347], [399,370], [378,358], [551,207], [476,211], [544,257], [508,208], [430,223], [385,274], [345,232], [564,425], [566,381], [567,320], [520,218], [492,174], [519,191], [554,160], [526,150], [548,187], [574,245], [571,208], [625,274], [584,300], [633,311], [647,345], [623,325], [626,296], [653,322], [595,323], [654,303], [635,361], [660,257], [505,279], [487,284], [502,262], [499,304], [441,222], [583,365], [553,363], [530,295], [622,344]]
dragPath = "M0,-300l100,100h-50v150h150v-50L300,0l-100,100v-50h-150v150h50L0,300l-100,-100h50v-150h-150v50L-300,0l100,-100v50h150v-150h-50z"
function domPoint(node, event) { if (event.changedTouches) event = event.changedTouches[0]; var svg = node.ownerSVGElement || node; if (svg.createSVGPoint) { var point = svg.createSVGPoint(); point.x = event.clientX, point.y = event.clientY; point = point.matrixTransform(node.getScreenCTM().inverse()); return {x: point.x, y: point.y}; } var rect = node.getBoundingClientRect(); return { x: event.clientX - rect.left - node.clientLeft, y: event.clientY - rect.top - node.clientTop }; }
d3 = require("https://d3js.org/d3.v5.min.js")