Flow study: fluid cells, openings and wall distances
Finds which grid cells are fluid, closes the openings, measures each opening's area, and cuts every wall link against the part for its wall distance.
How the precision is kept
Cell classification is whole-number work on the voxel grid. Opening areas, wall distances and the frontal area are 64-bit calculations on the part's triangles. Some of them are rasterized: an opening's holes are found on pixels of half a cell and the frontal area on a 1024-pixel grid, so those values are accurate to the pixel, not to the arithmetic.
Calculations under this policy
8 calculations in 2 source files. 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 area-weighted outward normal of a face, its area and how flat it is.
area = Σ |nₜ| ÷ 2; flatness = |Σ nₜ| ÷ (2 × area)
Cross products of triangle edge vectors, summed in triangle order in 64-bit. A face with no triangles returns the normal (0, 0, 1) and area 0 instead of dividing by zero.
flow-domain.js
lines 30 to 42
PartCraft 0.8.1
/** Area-weighted outward normal and area of a face from the mesh's triangles. */
export function faceNormal(mesh, fi) {
const P = mesh.positions, I = mesh.indices, tf = mesh.triFace, s = [0, 0, 0];
let area = 0;
for (let t = 0; t < I.length / 3; t++) {
if ((tf ? tf[t] : 0) !== fi) continue;
const a = 3 * I[3 * t], b = 3 * I[3 * t + 1], c = 3 * I[3 * t + 2];
const n = cross([P[b] - P[a], P[b + 1] - P[a + 1], P[b + 2] - P[a + 2]], [P[c] - P[a], P[c + 1] - P[a + 1], P[c + 2] - P[a + 2]]);
s[0] += n[0]; s[1] += n[1]; s[2] += n[2]; area += Math.hypot(n[0], n[1], n[2]) / 2;
}
const l = Math.hypot(s[0], s[1], s[2]);
return { n: l > 0 ? [s[0] / l, s[1] / l, s[2] / l] : [0, 0, 1], flat: area > 0 ? l / 2 / area : 0, area };
}
Area of an opening's holes
64-bit
The area of the holes in a flat face, which is the opening's true area.
loop area = ½ Σ (pₐ × p_b) · n; holes = Σ −area over clockwise loops
Each loop's area is summed from cross products of position vectors measured from the model origin, not from a point on the face. For an opening far from the origin the terms are large and nearly cancel, which costs digits in proportion to (distance ÷ opening size)².
flow-domain.js
lines 44 to 70
PartCraft 0.8.1
/**
* Area of the holes of a flat face (mm²), exact for its triangles: the face's boundary edges (used by one of its
* triangles) chained into loops; a loop turning clockwise about the outward normal bounds a hole. → { holes, loops }
*/
export function faceHoleArea(mesh, fi, n) {
const P = mesh.positions, I = mesh.indices, tf = mesh.triFace, cnt = new Map(), next = new Map();
const key = (a, b) => a * 4294967296 + b;
for (let t = 0; t < I.length / 3; t++) {
if ((tf ? tf[t] : 0) !== fi) continue;
for (let k = 0; k < 3; k++) { const a = I[3 * t + k], b = I[3 * t + (k + 1) % 3]; cnt.set(key(a, b), (cnt.get(key(a, b)) || 0) + 1); }
}
for (const [kk] of cnt) { const a = Math.floor(kk / 4294967296), b = kk % 4294967296; if (!cnt.has(key(b, a))) next.set(a, b); }
let holes = 0, loops = 0;
const seen = new Set();
for (const a0 of next.keys()) {
if (seen.has(a0)) continue;
let a = a0, area = 0, guard = 0;
do {
seen.add(a); const b = next.get(a); if (b === undefined) break;
const pa = [P[3 * a], P[3 * a + 1], P[3 * a + 2]], pb = [P[3 * b], P[3 * b + 1], P[3 * b + 2]];
area += dot(cross(pa, pb), n) / 2; a = b;
} while (a !== a0 && ++guard < 1e6);
loops++;
if (area < 0) holes -= area;
}
return { holes, loops };
}
Lid cells over an opening
64-bit
Finds the cells just behind a flat face that cover its holes.
raster area = hole pixels × (h ÷ 2)²
The face is rasterized on pixels of half a cell; a pixel is covered when its center is inside a triangle within 1e-9 in barycentric coordinates, and edges are sampled every sixth of a cell. The returned area is a pixel count, used only as a fallback for the loop area.
flow-domain.js
lines 72 to 128
PartCraft 0.8.1
/**
* The lid cells of a flat face: cells (not solid) just behind its plane whose centres project onto a hole of the face.
* → { cells: [grid index], holes: hole pixel count } or null when the face has no hole.
*/
export function lidCells(grid, solid, mesh, fi, n, o) {
const { h, nx, ny, nz, origin } = grid;
const P = mesh.positions, I = mesh.indices, tf = mesh.triFace;
const ref = Math.abs(n[0]) < 0.9 ? [1, 0, 0] : [0, 1, 0], u = unit(cross(n, ref)), v = cross(n, u);
const tris = [];
let umin = Infinity, umax = -Infinity, vmin = Infinity, vmax = -Infinity;
const lo = [Infinity, Infinity, Infinity], hi = [-Infinity, -Infinity, -Infinity];
for (let t = 0; t < I.length / 3; t++) {
if ((tf ? tf[t] : 0) !== fi) continue;
const pts = [0, 1, 2].map((k) => { const a = 3 * I[3 * t + k], p = [P[a], P[a + 1], P[a + 2]]; for (let c = 0; c < 3; c++) { lo[c] = Math.min(lo[c], p[c]); hi[c] = Math.max(hi[c], p[c]); } const d = sub(p, o); return [dot(d, u), dot(d, v)]; });
for (const [a, b] of pts) { umin = Math.min(umin, a); umax = Math.max(umax, a); vmin = Math.min(vmin, b); vmax = Math.max(vmax, b); }
tris.push(pts);
}
if (!tris.length) return null;
const ps = h / 2, W = Math.ceil((umax - umin) / ps) + 5, H = Math.ceil((vmax - vmin) / ps) + 5;
if (W * H > 16e6) return null;
const u0 = umin - 2 * ps, v0 = vmin - 2 * ps, cov = new Uint8Array(W * H);
const mark = (a, b) => { const i = Math.floor((a - u0) / ps), j = Math.floor((b - v0) / ps); if (i >= 0 && j >= 0 && i < W && j < H) cov[i + W * j] = 1; };
for (const [A, B, C] of tris) {
// pixel centres inside the triangle, and its edges sampled finely (conservative: thin rims stay closed)
const i0 = Math.max(0, Math.floor((Math.min(A[0], B[0], C[0]) - u0) / ps)), i1 = Math.min(W - 1, Math.floor((Math.max(A[0], B[0], C[0]) - u0) / ps));
const j0 = Math.max(0, Math.floor((Math.min(A[1], B[1], C[1]) - v0) / ps)), j1 = Math.min(H - 1, Math.floor((Math.max(A[1], B[1], C[1]) - v0) / ps));
const d = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]);
if (Math.abs(d) > 1e-18) for (let j = j0; j <= j1; j++) for (let i = i0; i <= i1; i++) {
const x = u0 + (i + 0.5) * ps, y = v0 + (j + 0.5) * ps;
const w1 = ((B[0] - x) * (C[1] - y) - (B[1] - y) * (C[0] - x)) / d, w2 = ((C[0] - x) * (A[1] - y) - (C[1] - y) * (A[0] - x)) / d, w3 = 1 - w1 - w2;
if (w1 >= -1e-9 && w2 >= -1e-9 && w3 >= -1e-9) cov[i + W * j] = 1;
}
for (const [p, q] of [[A, B], [B, C], [C, A]]) { const L = Math.hypot(q[0] - p[0], q[1] - p[1]), k = Math.max(1, Math.ceil(L / (ps / 3))); for (let s = 0; s <= k; s++) mark(p[0] + (q[0] - p[0]) * s / k, p[1] + (q[1] - p[1]) * s / k); }
}
// flood the uncovered pixels from the border: what it cannot reach is a hole
const out = new Uint8Array(W * H), stack = [];
const push = (i) => { if (!cov[i] && !out[i]) { out[i] = 1; stack.push(i); } };
for (let i = 0; i < W; i++) { push(i); push(i + W * (H - 1)); }
for (let j = 0; j < H; j++) { push(W * j); push(W - 1 + W * j); }
while (stack.length) { const p = stack.pop(), i = p % W, j = (p - i) / W; if (i > 0) push(p - 1); if (i < W - 1) push(p + 1); if (j > 0) push(p - W); if (j < H - 1) push(p + W); }
let holes = 0;
for (let p = 0; p < W * H; p++) if (!cov[p] && !out[p]) holes++;
if (!holes) return null;
// cells behind the plane (signed distance in (−T, 0], T the L1 thickness) over the face or its holes
const T = h * (Math.abs(n[0]) + Math.abs(n[1]) + Math.abs(n[2])), cells = [];
const box = [0, 1, 2].map((c) => [Math.max(0, Math.floor((lo[c] - T - h - origin[c]) / h)), Math.min([nx, ny, nz][c] - 1, Math.floor((hi[c] + T + h - origin[c]) / h))]);
for (let k = box[2][0]; k <= box[2][1]; k++) for (let j = box[1][0]; j <= box[1][1]; j++) for (let i = box[0][0]; i <= box[0][1]; i++) {
const e = i + nx * (j + ny * k);
if (solid[e]) continue;
const c = [origin[0] + (i + 0.5) * h, origin[1] + (j + 0.5) * h, origin[2] + (k + 0.5) * h], d = sub(c, o), s = dot(d, n);
if (!(s > -T - 1e-9 * h && s <= 1e-9 * h)) continue;
const pi = Math.floor((dot(d, u) - u0) / ps), pj = Math.floor((dot(d, v) - v0) / ps);
if (pi < 0 || pj < 0 || pi >= W || pj >= H || out[pi + W * pj]) continue;
cells.push(e);
}
return { cells, holes, area: holes * ps * ps };
}
Wall distance along a link
64-bit
Where the link from a fluid cell to a solid cell crosses the part's surface.
t = (e₂ · (tv × e₁)) ÷ (e₁ · (d × e₂))
The Möller-Trumbore crossing in 64-bit. A triangle is skipped when the determinant is under 1e-20 in absolute value, which is not scaled to the triangle's size; the barycentric and distance tests allow 1e-9. The fraction is clamped to [0, 1].
flow-domain.js
lines 130 to 149
PartCraft 0.8.1
/** Segment p → p + d against the triangles near it: the first crossing as { t (0…1), face, n (unit triangle normal along d) } or null. */
function firstHit(mesh, hash, p, d) {
const P = mesh.positions, I = mesh.indices, tf = mesh.triFace;
const mid = [p[0] + d[0] / 2, p[1] + d[1] / 2, p[2] + d[2] / 2];
let best = null;
for (const t of hash.near(mid)) {
const a = 3 * I[3 * t], b = 3 * I[3 * t + 1], c = 3 * I[3 * t + 2];
const e1 = [P[b] - P[a], P[b + 1] - P[a + 1], P[b + 2] - P[a + 2]], e2 = [P[c] - P[a], P[c + 1] - P[a + 1], P[c + 2] - P[a + 2]];
const pv = cross(d, e2), det = dot(e1, pv);
if (Math.abs(det) < 1e-20) continue;
const inv = 1 / det, tv = [p[0] - P[a], p[1] - P[a + 1], p[2] - P[a + 2]], uu = dot(tv, pv) * inv;
if (uu < -1e-9 || uu > 1 + 1e-9) continue;
const qv = cross(tv, e1), vv = dot(d, qv) * inv;
if (vv < -1e-9 || uu + vv > 1 + 1e-9) continue;
const s = dot(e2, qv) * inv;
if (s < -1e-9 || s > 1 + 1e-9) continue;
if (!best || s < best.t) { const nn = cross(e1, e2), l = Math.hypot(nn[0], nn[1], nn[2]) || 1, sg = dot(nn, d) >= 0 ? 1 / l : -1 / l; best = { t: Math.max(0, Math.min(1, s)), face: tf ? tf[t] : 0, n: [nn[0] * sg, nn[1] * sg, nn[2] * sg] }; }
}
return best;
}
Frontal area
64-bit
The area of the part seen along the flow, used as the reference area for drag and lift coefficients.
A pixel count: a pixel is covered when its center lies inside a projected triangle. The area is therefore accurate to about the outline length times one pixel, not to 15 digits, and the drag and lift coefficients divide by it.
flow-domain.js
lines 308 to 329
PartCraft 0.8.1
/** Frontal area of a mesh seen along `dir` (projected triangles rasterized on a res × res grid over their extent), mm². */
export function frontalArea(mesh, dir, res = 1024) {
const d = unit(dir), ref = Math.abs(d[0]) < 0.9 ? [1, 0, 0] : [0, 1, 0], u = unit(cross(d, ref)), v = cross(d, u);
const P = mesh.positions, I = mesh.indices, nv = P.length / 3, U = new Float64Array(nv), V = new Float64Array(nv);
let u0 = Infinity, u1 = -Infinity, v0 = Infinity, v1 = -Infinity;
for (let i = 0; i < nv; i++) { const p = [P[3 * i], P[3 * i + 1], P[3 * i + 2]]; U[i] = dot(p, u); V[i] = dot(p, v); u0 = Math.min(u0, U[i]); u1 = Math.max(u1, U[i]); v0 = Math.min(v0, V[i]); v1 = Math.max(v1, V[i]); }
const ps = Math.max(u1 - u0, v1 - v0) / res || 1, W = Math.ceil((u1 - u0) / ps) + 1, H = Math.ceil((v1 - v0) / ps) + 1, cov = new Uint8Array(W * H);
for (let t = 0; t < I.length / 3; t++) {
const a = I[3 * t], b = I[3 * t + 1], c = I[3 * t + 2];
const A = [U[a], V[a]], B = [U[b], V[b]], C = [U[c], V[c]], det = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]);
if (Math.abs(det) < 1e-18) continue;
const i0 = Math.max(0, Math.floor((Math.min(A[0], B[0], C[0]) - u0) / ps)), i1 = Math.min(W - 1, Math.floor((Math.max(A[0], B[0], C[0]) - u0) / ps));
const j0 = Math.max(0, Math.floor((Math.min(A[1], B[1], C[1]) - v0) / ps)), j1 = Math.min(H - 1, Math.floor((Math.max(A[1], B[1], C[1]) - v0) / ps));
for (let j = j0; j <= j1; j++) for (let i = i0; i <= i1; i++) {
const x = u0 + (i + 0.5) * ps, y = v0 + (j + 0.5) * ps;
const w1 = ((B[0] - x) * (C[1] - y) - (B[1] - y) * (C[0] - x)) / det, w2 = ((C[0] - x) * (A[1] - y) - (C[1] - y) * (A[0] - x)) / det;
if (w1 >= 0 && w2 >= 0 && w1 + w2 <= 1) cov[i + W * j] = 1;
}
}
let n = 0; for (const c of cov) n += c;
return n * ps * ps;
}
Picks the cell size that gives about a target number of cells in the grid box.
h = ∛(volume ÷ target), rounded up to 1, 1.25, 1.5, 2, 2.5, 3, 4, 5, 6, 8 or 10 times a power of ten
The size is rounded up onto a preferred series with 1e-9 allowance and then rounded to 3 significant digits with toPrecision(3). This is a mesh choice; it sets the resolution of the result, not its digits.
flow-study.js
lines 23 to 27
PartCraft 0.8.1
function autoSize(volume, target) {
const raw = Math.cbrt(volume / target), e = Math.floor(Math.log10(raw)), b = raw / 10 ** e;
const s = [1, 1.25, 1.5, 2, 2.5, 3, 4, 5, 6, 8, 10].find((x) => x >= b - 1e-9) || 10;
return +(s * 10 ** e).toPrecision(3);
}
Cell count estimate
64-bit
The number of grid cells a study will use.
n = ⌈extent ÷ h − 1e-6⌉ + margin cells, per axis
The 1e-6 allowance stops an extent that is a whole number of cells within rounding from gaining one more cell.
The direction lift is measured along, made square to the flow.
l = l₀ − (l₀ · d) d, then divided by its length
One Gram-Schmidt step in 64-bit. If less than 1e-6 is left after removing the flow direction, the default +Y (or +Z) is returned instead of dividing by a near-zero length.
flow-study.js
lines 111 to 118
PartCraft 0.8.1
/** The lift direction: the user's (made square to the flow) or +Y (+Z when the flow is mostly along Y). */
export function liftDir(ext) {
const d = unitv(ext.dir);
let l = ext.lift ? unitv(ext.lift) : Math.abs(d[1]) > 0.9 ? [0, 0, 1] : [0, 1, 0];
const k = l[0] * d[0] + l[1] * d[1] + l[2] * d[2];
l = l.map((c, i) => c - k * d[i]);
return Math.hypot(...l) > 1e-6 ? unitv(l) : (Math.abs(d[1]) > 0.9 ? [0, 0, 1] : [0, 1, 0]);
}
In detail
Method
Voxelization and flood fills; opening hole area from the face's boundary loops; ray and triangle crossing (the Möller-Trumbore form) for wall distances; pixel counting for the frontal area.
Number types
64-bit numbers; Int32Array and Uint8Array for cells and pixel masks. The mesh handed to these functions has 32-bit positions (see the 32-bit policy for the flow study).
Tolerances
Ray crossing: determinant under 1e-20 skipped, barycentric and distance margins 1e-9. A face counts as flat when its summed normal is over 0.999 of its area. Cell size: the automatic size is rounded to 3 significant digits on a 1, 1.25, 1.5, 2 … series.
Limits
The frontal area is a pixel count at 1024 pixels across the larger extent. An opening's raster area is used only when its boundary loops give no hole area.
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.