function solveEquation(terms) {
const parsed = terms
.slice()
.sort()
.map(t => parseTerm(t));
const [x1, y1, c1] = parsed[0];
const [x2, y2, c2] = parsed[1];
const [x3, y3, c3] = parsed[2];
return gauss([[x2 - x1, y2 - y1], [x3 - x2, y3 - y2]], [c1 - c2, c2 - c3]);
}
function gauss(A, x) {
// thanks to https://github.com/itsravenous/gaussian-elimination/blob/master/gauss.js
var i, k, j;
// Just make a single matrix
for (i = 0; i < A.length; i++) {
A[i].push(x[i]);
}
var n = A.length;
for (i = 0; i < n; i++) {
// Search for maximum in this column
var maxEl = abs(A[i][i]),
maxRow = i;
for (k = i + 1; k < n; k++) {
if (abs(A[k][i]) > maxEl) {
maxEl = abs(A[k][i]);
maxRow = k;
}
}
// Swap maximum row with current row (column by column)
for (k = i; k < n + 1; k++) {
var tmp = A[maxRow][k];
A[maxRow][k] = A[i][k];
A[i][k] = tmp;
}
// Make all rows below this one 0 in current column
for (k = i + 1; k < n; k++) {
var c = -A[k][i] / A[i][i];
for (j = i; j < n + 1; j++) {
if (i === j) {
A[k][j] = 0;
} else {
A[k][j] += c * A[i][j];
}
}
}
}
// Solve equation Ax=b for an upper triangular matrix A
x = array_fill(0, n, 0);
for (i = n - 1; i > -1; i--) {
x[i] = A[i][n] / A[i][i];
for (k = i - 1; k > -1; k--) {
A[k][n] -= A[k][i] * x[i];
}
}
return x;
}
function array_fill(i, n, v) {
var a = [];
for (; i < n; i++) {
a.push(v);
}
return a;
}
function drawNewtonTriangulation3D(triangulation, scene) {
triangulation.forEach(face => {
const projectedUp = face.map(v =>
v.map((k, i) => (i === 2 ? maxCoeff + 0.5 : k))
);
addFace(projectedUp, scene, { wireframe: true });
});
}
function drawNewtonTriangulation2D(triangulation, context, options) {
const opts = Object.assign({ gridStep: gridStep }, options);
for (let i = 0; i < triangulation.length; ++i) {
const t = triangulation[i];
triangle(t, context, { color: "gray", lineWidth: 3 });
}
points.forEach(p =>
circle([p[0] * opts.gridStep, p[1] * opts.gridStep, 4], context)
);
}
function drawTropicalCurve(triangulation, polynom, context, options) {
const opts = Object.assign({ gridStep: gridStep, r: 4 }, options);
const tropicalPoints = findTropicalPoints(
triangulation,
polynom,
opts.gridStep
);
const angles = findAngles(triangulation);
const neighbors = findNeighbors(triangulation);
tropicalPoints.forEach((p, i) => {
circle([...p, opts.r], context);
neighbors[i].forEach((n, j) => {
if (n !== 0 && !n) {
const inf = coordsFromDeg(angles[i][j] - 90, 1000, p);
line(p, inf, context);
} else {
const neigbor = tropicalPoints[n];
line(p, neigbor, context);
}
});
});
}
function draw3DHull(hull, scene) {
for (let [i1, i2, i3] of hull) {
const face = [points[i1], points[i2], points[i3]];
const color = randomColor();
addFace(face, scene, { color });
}
}
function drawGrid2D(context, gs = gridStep) {
for (let i = -gridSize; i < gridSize; ++i) {
for (let j = -gridSize; j < gridSize; ++j) {
const color = "rgba(200, 200, 200, 1)";
line([i * gs, 0], [i * gs, j * gs], context, { color });
line([0, i * gs], [j * gs, i * gs], context, { color });
}
}
}
function drawExamplePolynomialCurve() {
const context = DOM.context2d(300, 300);
context.translate(150, 150);
circle([0, 0, 3], context);
textAt([5, 15], "(0,0)", 11, context);
line([0, 0], [-300, 0], context);
line([0, 0], [0, 300], context);
line([0, 0], [300, -300], context);
// arrows
context.moveTo(-149, 0);
context.lineTo(-139, 5);
context.lineTo(-139, -5);
context.fill();
context.moveTo(0, 149);
context.lineTo(5, 139);
context.lineTo(-5, 139);
context.fill();
context.moveTo(149, -149);
context.lineTo(135, -141);
context.lineTo(142, -133);
context.fill();
return context.canvas;
}
function findNeighbor([l1, l2], faceIndex, triangulation) {
let neighbor = false;
for (let i = 0; i < triangulation.length; ++i) {
if (i === faceIndex) continue;
const [p1, p2, p3] = triangulation[i];
const lines = [[p1, p2], [p2, p3], [p3, p1]];
lines.forEach(([p1, p2]) => {
const p1Match = dist2D(p1, l1) === 0 || dist2D(p1, l2) === 0;
const p2Match = dist2D(p2, l1) === 0 || dist2D(p2, l2) === 0;
if (p1Match && p2Match) {
neighbor = i;
}
});
}
return neighbor;
}
function findAngles(triangulation) {
const res = [];
triangulation.forEach(([p1, p2, p3]) => {
res.push([
getAlphaBetweenPoints(p1, p2),
getAlphaBetweenPoints(p2, p3),
getAlphaBetweenPoints(p3, p1)
]);
});
return res;
}
function findNeighbors(triangulation) {
const res = [];
for (let i = 0; i < triangulation.length; ++i) {
const [p1, p2, p3] = triangulation[i];
const lines = [[p1, p2], [p2, p3], [p3, p1]];
const neighbors = [];
for (let l of lines) {
neighbors.push(findNeighbor(l, i, triangulation));
}
res.push(neighbors);
}
return res;
}
function findTropicalPoints(triangulation, poly, gs = gridStep) {
const res = [];
triangulation.forEach((t, i) => {
const term = findTermsFromPoints(t, poly);
const p = solveEquation(term).map(v => v * gs);
res.push(p);
});
return res;
}
function findTermsFromPoints(points, poly) {
const res = [];
for (let [x, y, z] of points) {
poly.forEach(({ point, term }, i) => {
const [px, py, pz] = point;
if (px === x && py === y && pz === z) res.push(term);
});
}
return res;
}
function findNewtonTriangulation(p) {
const hull = qh3D(p);
const hullUpperFaces = [];
for (let [i1, i2, i3] of hull) {
const face = [p[i1], p[i2], p[i3]];
const [f, _] = makeFace(face);
const normals = f.normal;
const isUpper = normals.angleTo(new THREE.Vector3(0, 1, 0)) >= Math.PI / 2;
const areaNotNull = areaTriangle(...face) > 0;
if (isUpper && areaNotNull) {
hullUpperFaces.push(face);
}
}
return hullUpperFaces;
}
function newtonTriangulation2D(points, options) {
const opts = Object.assign({ gridStep: gridStep }, options);
const triangulation = findNewtonTriangulation(points);
const projected = projectTo2D(triangulation, opts);
// filter triangles that are colinear
const res = [];
for (let i = 0; i < projected.length; ++i) {
if (areaTriangle(...projected[i]) > 0) {
res.push(projected[i]);
}
}
return res;
}
function projectTo2D(faces, options) {
const opts = Object.assign({ gridStep: gridStep }, options);
return faces.map(f =>
f.map(v => [v[0] * opts.gridStep, v[1] * opts.gridStep])
);
}
function addCube([x, z, y], scene, options) {
const opts = Object.assign({ color: 0x000000, r: 2 }, options);
const geometry = new THREE.BoxBufferGeometry(opts.r, opts.r, opts.r);
const material = new THREE.MeshBasicMaterial({ color: opts.color });
const cube = new THREE.Mesh(geometry, material);
cube.position.set(x * gridStep, y * gridStep, z * gridStep);
scene.add(cube);
}
function addLine([[x1, z1, y1], [x2, z2, y2]], scene, options) {
const opts = Object.assign(
{ color: "rgba(100,100,100,1)", linewidth: 2 },
options
);
const geometry = new THREE.Geometry();
geometry.vertices.push(
new THREE.Vector3(x1 * gridStep, y1 * gridStep, z1 * gridStep),
new THREE.Vector3(x2 * gridStep, y2 * gridStep, z2 * gridStep)
);
const material = new THREE.LineBasicMaterial({
color: opts.color,
linewidth: opts.linewidth
});
const line = new THREE.Line(geometry, material);
scene.add(line);
}
function addPolygon(p, scene) {
const points = p.slice().sort(([x1, y1], [x2, y2]) => x1 - x2);
const geometry = new THREE.Geometry();
// const pl = points.length;
for (let i = 0; i < points.length; ++i) {
const [x, z, y] = points[i];
geometry.vertices.push(
new THREE.Vector3(x * gridStep, y * gridStep, z * gridStep)
);
if (i >= 2) {
geometry.faces.push(new THREE.Face3(i - 2, i - 1, i));
}
}
const material = new THREE.MeshBasicMaterial({
color: "rgba(200, 200, 200, 1)",
side: THREE.DoubleSide
});
const mesh = new THREE.Mesh(geometry, material);
scene.add(mesh);
}
function addFace(points, scene, options) {
const opts = Object.assign({ color: 0x000000, wireframe: false }, options);
const [face, geometry] = makeFace(points);
const material = new THREE.MeshBasicMaterial({
color: opts.color,
wireframe: opts.wireframe,
side: THREE.DoubleSide
});
const mesh = new THREE.Mesh(geometry, material);
scene.add(mesh);
}
function addGrid([x, y], scene, options) {
const opts = Object.assign({ color: 0x000000 }, options);
const points = [];
let counter = 0;
for (let i = -x; i < x; ++i) {
for (let j = -y; j < y; ++j) {
const c = [i, j, 0];
points.push(c);
if (points.length > 1 && j !== -x) {
addLine([points[counter - 1], c], scene);
}
counter += 1;
if (i > -x) {
addLine([points[counter - 1], points[counter - y - x - 1]], scene);
}
}
}
}
function makeFace(points) {
const geometry = new THREE.Geometry();
for (let i = 0; i < 3; ++i) {
geometry.vertices.push(
new THREE.Vector3(
points[i][0] * gridStep,
points[i][2] * gridStep,
points[i][1] * gridStep
)
);
}
const face = new THREE.Face3(0, 1, 2);
geometry.faces.push(face);
geometry.computeBoundingSphere();
geometry.computeFaceNormals();
return [face, geometry];
}
scene = {
const scene = new THREE.Scene();
scene.background = new THREE.Color(0xffffff);
addGrid([gridSize, gridSize], scene);
const hull2D = convexHull2D(poly.map(({ point }) => [point[0], point[1]]));
addPolygon(hull2D, scene);
poly.forEach(({ point }) => {
addCube(point, scene, { color: "red", r: 5 });
addLine([point, [point[0], point[1], 0]], scene);
addCube([point[0], point[1], 0], scene, { r: 4 });
});
const hull = qh3D(points);
draw3DHull(hull, scene);
const triangulation = findNewtonTriangulation(points);
drawNewtonTriangulation3D(triangulation, scene);
yield scene;
}
camera = {
const fov = 50;
const near = 0.1;
const far = 10000;
const camera = new THREE.PerspectiveCamera(fov, 1, near, far);
camera.lookAt(new THREE.Vector3(0, 0, 0));
camera.position.set(500, 1800, 1500);
yield camera;
}
function randomColor() {
return `rgba(${Math.floor(Math.random() * 255)},
${Math.floor(Math.random() * 255)},
${Math.floor(Math.random() * 255)},1)`;
}
function buildPoly(polynom) {
const terms = toConventionalNotation(polynom);
const points = terms.map(t => parseTerm(t));
const res = [];
for (let i = 0; i < points.length; ++i) {
res.push({ point: points[i], term: terms[i] });
}
return res;
}
function toConventionalNotation(polynom) {
const splited = polynom.split("+");
const prodAsSum = splited.map(term => {
const asSum = term.replace(/(\(|\))/g, "").replace(/\*/g, "+");
const sumTerms = asSum.split("+");
const withExponent = sumTerms.map(t => {
const [base, exp] = t.split("^");
return exp ? `${exp}${base}` : t;
});
return withExponent.join("+").replace(/\s/g, "");
});
return prodAsSum;
}
function parseTerm(t) {
const term = t.replace(/\s/g, "");
let xTerm = term.match(/\(?-?\d*\/?\d*\)?\*?x/g) || ["0"];
let yTerm = term.match(/\(?-?\d*\/?\d*\)?\*?y/g) || ["0"];
let constantTerm = term.split("+").filter(v => !v.match(/(x|y)/g));
// defaults
if (constantTerm.length === 0) constantTerm = ["0"];
if (xTerm[0] === "x") xTerm = ["1"];
if (xTerm[0] === "-x") xTerm = ["-1"];
if (yTerm[0] === "y") yTerm = ["1"];
if (yTerm[0] === "-y") yTerm = ["-1"];
// cleanup
xTerm = xTerm.map(v => cleanUpTerm(v))[0];
yTerm = yTerm.map(v => cleanUpTerm(v))[0];
constantTerm = constantTerm.map(v => cleanUpTerm(v))[0];
return [xTerm, yTerm, constantTerm];
}
function cleanUpTerm(term) {
// first time ever i use eval ... not sure if evil or right...
return eval(term.replace(/(\(|\)|\*|x|y)/g, ""));
}
maxCoeff = points.reduce(
(acc, curr) => (curr[2] >= acc ? curr[2] : acc),
points[0][2]
)
points = poly.map(({ point }) => point)
function header() {
const context = DOM.context2d(width, width / 2);
context.translate(width / 2 - 5, width / 3 - 30);
context.scale(1, -1);
const headerPoly = buildPoly(headerPolynom);
const headerPoints = headerPoly.map(({ point }) => point);
const triangulation = findNewtonTriangulation(headerPoints);
drawTropicalCurve(triangulation, headerPoly, context, {
gridStep: 30,
r: 2.3
});
return context.canvas;
}
import {
line,
circle,
triangle,
areaTriangle,
textAt,
dist2D,
getAlphaBetweenPoints,
coordsFromDeg,
} from '@timhau/geometry'
One platform to build and deploy the best data apps
Experiment and prototype by building visualizations in live JavaScript notebooks. Collaborate with your team and decide which concepts to build out.
Use Observable Framework to build data apps locally. Use data loaders to build in any language or library, including Python, SQL, and R.
Seamlessly deploy to Observable. Test before you ship, use automatic deploy-on-commit, and ensure your projects are always up-to-date.