Finds where a 3D point of an edge lies on a free-form curve or surface, so the face it bounds can be laid out flat and triangulated.
How the precision is kept
These are searches and iterations in 64-bit whose stopping rule sets the accuracy. Surface derivatives are numerical central differences. The parameters found only place points in the (u, v) layout: the tessellator keeps each boundary point's own 3D position, so an inexact parameter changes the shape of triangles, not where the edges of the solid are.
Calculations under this policy
5 calculations in 1 source file. The code is read from the released app when this page loads. It is the calculation people run, not a copy kept for this page.
The ellipse parameter of the point on the ellipse nearest a given point.
f(t) = (E(t) − q)·E′(t) = 0, solved by t ← t − f ÷ f′
Starts from atan2(l_y ÷ b, l_x ÷ a) and takes exactly eight Newton steps; there is no test that the steps have become small. The loop leaves early only when the second derivative is 1e-30 or less in size.
step/geom.js
lines 59 to 67
BimCraft 0.8.1
let t = Math.atan2(l[1] / b, l[0] / a);
for (let i = 0; i < 8; i++) { // Newton on the distance (the eccentric angle is not the closest point when a ≠ b)
const x = a * Math.cos(t), y = b * Math.sin(t), dx = -a * Math.sin(t), dy = b * Math.cos(t);
const f = (x - l[0]) * dx + (y - l[1]) * dy, df = dx * dx + dy * dy + (x - l[0]) * -x + (y - l[1]) * -y;
if (!(Math.abs(df) > 1e-30)) break;
t -= f / df;
}
if (t < 0) t += TAU;
return t % TAU;
Parameter of a point on a B-spline curve
64-bit
The curve parameter nearest a given point.
A scan of max(32, 8 × control points) + 1 even samples picks a start, then the best of three points is kept and the step is halved when neither neighbor is closer. It stops when the step falls below 1e-13 of the parameter range or after 200 rounds. Distances are compared, not derivatives, so near the minimum two parameters can tie within rounding.
step/geom.js
lines 123 to 136
BimCraft 0.8.1
const m = Math.max(32, 8 * (n + 1));
const span0 = range[1] - range[0] || 1;
let best = range[0], bd = Infinity;
for (let k = 0; k <= m; k++) { const t = range[0] + span0 * k / m; const d = len(sub(at(t), q)); if (d < bd) { bd = d; best = t; } }
const lim = (t) => (per ? t : Math.max(range[0], Math.min(range[1], t)));
let step = span0 / m;
for (let it = 0; it < 200 && step > 1e-13 * span0; it++) {
const a = lim(best - step), b = lim(best + step);
const da = len(sub(at(a), q)), db = len(sub(at(b), q));
if (da < bd && da <= db) { bd = da; best = a; }
else if (db < bd) { bd = db; best = b; }
else step /= 2;
}
return per ? range[0] + ((((best - range[0]) % per) + per) % per) : best;
Parameters of a point on a B-spline surface
64-bit
The (u, v) of the surface point nearest a given point.
[S_u·S_u S_u·S_v; S_u·S_v S_v·S_v] d = −[r·S_u; r·S_v], r = S(u, v) − q
Gauss–Newton with derivatives from central differences of step 1e-6 × the parameter range, so each derivative is a difference of two nearby surface points divided by a small number. At most 30 rounds; it stops when both steps are below 1e-13 of the range. Where one derivative vanishes (a pole) it steps along the other direction alone.
step/geom.js
lines 183 to 207
BimCraft 0.8.1
const project = (q, hint) => {
let u, v;
if (hint) [u, v] = hint;
else { let bd = Infinity; for (const [a, b, p] of samples()) { const d = len(sub(p, q)); if (d < bd) { bd = d; u = a; v = b; } } }
const hu = 1e-6 * ((ru[1] - ru[0]) || 1), hv = 1e-6 * ((rv[1] - rv[0]) || 1);
for (let it = 0; it < 30; it++) {
const S = at(u, v);
const Su = mul(sub(at(u + hu, v), at(u - hu, v)), 1 / (2 * hu)), Sv = mul(sub(at(u, v + hv), at(u, v - hv)), 1 / (2 * hv));
const r = sub(S, q);
// Gauss–Newton on |S − q|²: [Su·Su Su·Sv; Su·Sv Sv·Sv] d = −[r·Su; r·Sv]
const a11 = dot(Su, Su), a12 = dot(Su, Sv), a22 = dot(Sv, Sv), b1 = -dot(r, Su), b2 = -dot(r, Sv);
const det = a11 * a22 - a12 * a12;
let du, dv;
if (Math.abs(det) > 1e-12 * Math.max(1e-30, a11 * a22)) { du = (b1 * a22 - b2 * a12) / det; dv = (a11 * b2 - a12 * b1) / det; }
// at a pole one derivative vanishes (the row collapses to a point): step along the other direction alone
else if (a22 > 1e-30 && a22 >= a11) { du = 0; dv = b2 / a22; }
else if (a11 > 1e-30) { du = b1 / a11; dv = 0; }
else break;
const nu2 = periodU ? u + du : Math.max(ru[0], Math.min(ru[1], u + du)), nv2 = periodV ? v + dv : Math.max(rv[0], Math.min(rv[1], v + dv));
du = nu2 - u; dv = nv2 - v; u = nu2; v = nv2;
if (Math.abs(du) < 1e-13 * ((ru[1] - ru[0]) || 1) && Math.abs(dv) < 1e-13 * ((rv[1] - rv[0]) || 1)) break;
}
if (periodU) u = ru[0] + ((((u - ru[0]) % periodU) + periodU) % periodU);
if (periodV) v = rv[0] + ((((v - rv[0]) % periodV) + periodV) % periodV);
return [u, v];
Surface of linear extrusion
64-bit
The point of a curve swept along a vector, and the (u, v) of a point.
S(u, v) = C(u) + v · vec; v = (q − C(u))·d ÷ |vec|, u = project(q − v · vec)
The point itself is one 64-bit sum. The inverse alternates between v and the curve's own projection for at most 10 rounds and stops when u changes by less than 1e-12, an absolute figure in the curve's parameter units. Its mesh scale comes from 16 samples of the curve.
step/geom.js
lines 294 to 306
BimCraft 0.8.1
export function extrusion(curve, vec) {
const d = norm(vec), m = len(vec);
const span = curve.range ? curve.range[1] - curve.range[0] : 1;
let L = 0; { const k = 16; let prev = curve.at(curve.range ? curve.range[0] : 0); for (let i = 1; i <= k; i++) { const q = curve.at((curve.range ? curve.range[0] : 0) + span * i / k); L += len(sub(q, prev)); prev = q; } }
return { kind: 'extrusion', at: (u, v) => add(curve.at(u), mul(vec, v)), periodU: curve.period || 0, periodV: 0, scale: [L / (span || 1), m || 1], flat: curve.kind === 'line' || curve.kind === 'polyline',
project: (q) => {
// the point moved back along the extrusion onto the curve's plane is not known in general: iterate
let v = 0, u = curve.project(q);
for (let it = 0; it < 10; it++) { const c = curve.at(u); v = dot(sub(q, c), d) / (m || 1); const nu = curve.project(sub(q, mul(vec, v))); if (Math.abs(nu - u) < 1e-12) break; u = nu; }
return [u, v];
},
singular: () => null };
}
Offset surface
64-bit
The point of a surface moved a distance along its normal.
S_off(u, v) = S(u, v) + d · n, n = (ΔS_u × ΔS_v) ÷ |…|
The normal is the cross product of two central differences taken with a fixed step of 1e-6 in parameter units, whatever the surface's parameter range is. Where that cross product is exactly zero the normal falls back to (0, 0, 1). The inverse is the basis surface's, so it ignores the offset.
step/geom.js
lines 307 to 316
BimCraft 0.8.1
/** An offset surface, approximated: the basis moved along its numerical normal by d (exact for planes, cylinders, spheres). */
export function offset(basis, d) {
const nAt = (u, v) => {
const hu = 1e-6, hv = 1e-6;
const Su = sub(basis.at(u + hu, v), basis.at(u - hu, v)), Sv = sub(basis.at(u, v + hv), basis.at(u, v - hv));
return norm(cross(Su, Sv));
};
return { ...basis, kind: 'offset', flat: basis.flat, at: (u, v) => add(basis.at(u, v), mul(nAt(u, v), d)),
project: (q, hint) => basis.project(q, hint) };
}
In detail
Method
Ellipse: Newton on the squared distance from the eccentric angle. B-spline curve: a coarse scan then a shrinking three-point search. B-spline surface: a sample grid for the first guess, then Gauss–Newton with central-difference derivatives. Extrusion: alternating projection. Offset surface: the basis surface moved along a numerical normal.
Number types
64-bit numbers throughout.
Tolerances
Ellipse: 8 Newton steps, no convergence test; stops early only when the second derivative is under 1e-30. Curve search: starts at range ÷ max(32, 8 × control points), halves until the step is under 1e-13 × range, at most 200 rounds. Surface: difference step 1e-6 × parameter range; at most 30 rounds; stops when both steps are under 1e-13 × range; the 2 × 2 system counts as singular when its determinant is under 1e-12 × the product of its diagonal. Extrusion: at most 10 rounds, stops when u changes by less than 1e-12 (absolute, in parameter units). Offset normal: difference step 1e-6 in parameter units (absolute).
Limits
The surface search can stop at a local closest point when the first guess from the sample grid is on another part of the surface. An offset surface is an approximation: its points use a numerical normal, and its inverse is the basis surface's.
Approach: 64-bit, stated limit
The arithmetic is 64-bit, but the method or the data sets the accuracy: a solver tolerance, a mesh, a step size, a datum shift, a table. The limit is stated where the result is shown.