thresh2 = Array.from({ length: 4 }, (_, i) => Math.pow(2.5, i-1))
Getting actual contours
https://en.wikipedia.org/wiki/Marching_squares
{
const img_width = 600;
const height = 300;
let ranges = {x: [-100, 100], y: [-50, 50]};
let graph = createStandardGrid(ranges.x, ranges.y);
// Calculate field
let nodes = [
// earth
{
x: -5, y: -5,
//m: 7.342e22,
m: 1000,
vx: 0, vy: 0
}
];
graph.svg.append('circle')
.attr('r', 4)
.attr('cx', graph.x(nodes[0].x))
.attr('cy', graph.y(nodes[0].y))
.attr('fill', 'gray');
// const G = 6.674e-11 / (unitDist ** 3);
// Draw field
const F = gravityField(nodes, 1);
let s0 = makeStateVec(nodes);
const makeGrid = () => {
const q = pixelsPerSample; // The level of detail, e.g., sample every 4 pixels in x and y.
const x0 = -q / 2, x1 = img_width + q; // not sure where 28 comes from here...
const y0 = -q / 2, y1 = height + q;
const n = Math.ceil((x1 - x0) / q);
const m = Math.ceil((y1 - y0) / q);
const grid = new Array(n * m);
let pts = [];
for (let j = 0; j < m; ++j) {
for (let i = 0; i < n; ++i) {
let pt = [graph.x.invert(i*q + x0), graph.y.invert(j*q + y0)];
pts.push(pt);
let gvr = F(graph.x.invert(i*q + x0), graph.y.invert(j*q + y0), s0);
grid[j * n + i] = Math.sqrt(gvr.xpp**2 + gvr.ypp**2);
}
}
grid.x = -q;
grid.y = -q;
grid.k = q;
grid.n = n;
grid.m = m;
return [grid, pts];
}
const [grid, pts] = makeGrid();
// https://github.com/d3/d3-geo/blob/main/README.md#path_projection
// https://github.com/d3/d3-geo/blob/main/README.md#transforms
let geoTrans = d3.geoTransform({
point: function(x, y) {
// translate from "sample domain" into "pixel domain"
// note: we don't return to the original domain
this.stream.point(grid.x + grid.k*x, grid.y + grid.k*y);
}
});
const contours = d3.contours()
.size([grid.n, grid.m])
.thresholds(thresh2)
(grid)
//.map(transform);
//color = d3.scaleSequentialLog(d3.extent(thresholds), d3.interpolateMagma)
pts.map(([x, y]) => {
graph.svg.append('circle')
.attr('r', 1.5)
.attr('cx', graph.x(x))
.attr('cy', graph.y(y))
.attr('fill', 'lightgray');
});
graph.svg.append("g")
.attr("fill", "none")
.attr("stroke", "#000")
.attr("stroke-opacity", 0.5)
.selectAll("path")
.data(contours)
.join("path")
//.attr("fill", d => color(d.value))
.attr("d", d3.geoPath(geoTrans));
return graph.svg.node();
}
thresh1 = Array.from({ length: 30 }, (_, i) => Math.pow(2, i-20))
Calculating gravity
makeStateVec = (nodes) => {
return nodes.map(n => n.x)
.concat(nodes.map(n => n.y))
.concat(nodes.map(n => n.vx))
.concat(nodes.map(n => n.vy));
}
gravityField = (nodes, G) => {
// nodes are the nodes exhibiting gravitational forces, generally celestial bodies
// the s input below is supposed to be the values affecting each other, but since we're just doing fields,
// s will always be the initial state.
G = G || 1;
const n = nodes.length;
// Borrowed from https://observablehq.com/@mcmcclur/gravitation-and-the-n-body-problem
let x_accel = (x1, y1, s, j) => {
let m2 = nodes[j].m;
// In the n-body simulation, all the nodes are moving, but in this modified version
// nothing is moving, we're just grabbing the point in time acceleration.
// let x1 = s[i];
// let y1 = s[i + n];
let x2 = s[j];
let y2 = s[j + n];
return (m2 * (x2 - x1)) / ((x2 - x1) ** 2 + (y2 - y1) ** 2) ** (3 / 2);
};
let y_accel = (x1, y1, s, j) => {
let m2 = nodes[j].m;
// let x1 = s[i];
// let y1 = s[i + n];
let x2 = s[j];
let y2 = s[j + n];
return (m2 * (y2 - y1)) / ((x2 - x1) ** 2 + (y2 - y1) ** 2) ** (3 / 2);
};
let gr = (x, y, s) => {
return {
xpp: G*d3.sum(d3.range(nodes.length).map((i) => x_accel(x, y, s, i))),
ypp: G*d3.sum(d3.range(nodes.length).map((i) => y_accel(x, y, s, i))),
};
};
return gr
}
import {createStandardGrid, drawVector, addArrowHead, vectorAdd, vectorScale} from "@pcarleton/range-kutta-supporting-functions"
d3 = require('d3@7.3.0')
Purpose-built for displays of data
Observable is your go-to platform for exploring data and creating expressive data visualizations. Use reactive JavaScript notebooks for prototyping and a collaborative canvas for visual data exploration and dashboard creation.
