Resources Numerical accuracy Calculation policy
Places sprinkler heads on a grid in a room, builds the branch and main network, sizes the pipes and finds the pressure and flow the design area needs against the water supply.
The calculation runs in SI units in 64-bit. Its accuracy is set by the method: pipe friction is the empirical Hazen-Williams formula with a coefficient and exponents of 4 digits, the network pressures are found by repeated substitution to within 1 Pa, and the source pressure by 40 bisection steps. Water density is the constant 998.2 kg/m³ and gravity the standard 9.80665 m/s².
Open item. Report when the network iteration reaches 300 passes without meeting the 1 Pa tolerance.
16 calculations in 3 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.
Tests tools/bimcraft-test/eg-mep-test.mjs
Converts a head's K-factor from L/min per √bar to m³/s per √Pa.
k = K ÷ 1000 ÷ 60 ÷ √100000
The factors 1000, 60 and 100000 are exact; the square root of 100000 is the one rounded step.
/** K in L/min per √bar → m³/s per √Pa. */
export const kSi = (K) => (K / 1000 / 60) / Math.sqrt(1e5);The friction loss in pascals per meter of pipe for a flow, an inside diameter and a roughness coefficient C.
Δp per m = 10.67 × Q^1.852 ÷ (C^1.852 × D^4.871) × 998.2 × 9.80665
The coefficient 10.67 and the exponents 1.852 and 4.871 are the rounded constants of the empirical formula, and 998.2 is a rounded water density; they, not the 64-bit arithmetic, limit the result. A flow or diameter of zero or less returns 0.
export const hw = (Q, D, C) => (Q > 0 && D > 0 ? (10.67 * Q ** 1.852 / (C ** 1.852 * D ** 4.871)) * RHO * G : 0); // Pa per metreChooses for each pipe the smallest steel schedule 40 size that keeps the velocity within the limit.
Q = heads below × flow per head ÷ 1000 × 1.15; v = Q ÷ (π d² ÷ 4) ≤ vmax
The inside diameter is the table's inches times the exact 0.0254. The flow carries a fixed factor of 1.15 and the default velocity limit is 6 m/s; if no size meets it the largest allowed size is used.
export function sizeNetwork(net, perHeadLps, opts = {}) {
const sizes = sizesList(opts.minSize, opts.maxSize), vmax = opts.vmax ?? 6, B = below(net);
for (const e of net.edges) {
const Q = B.count(e.b) * perHeadLps / 1000 * 1.15;
e.size = sizes.find((s) => Q / (Math.PI * idOf(s) ** 2 / 4) <= vmax) || sizes[sizes.length - 1];
}
return net;
}The groups of heads tried as the design area: blocks of rows at the far and near ends of the branch lines.
heads = ceil(area ÷ coverage); per row m = ceil(1.2 × √area ÷ spacing); rows = ceil(heads ÷ m)
Counts are ceilings of 64-bit quotients, held within the grid size. More than 60 candidates are thinned to every n-th one, so the governing area is the worst of those tested.
/** Design areas to try: k rows × m heads, at the ends of the branch lines farthest from and nearest to the main, on several rows. */
export function candidates(layout, net, hazard) {
const per = layout.coveragePerHead, nHeads = Math.max(1, Math.ceil(hazard.area / per));
const m = Math.max(1, Math.min(Math.ceil(1.2 * Math.sqrt(hazard.area) / (layout.sx / 1000)), layout.nx));
const k = Math.max(1, Math.min(Math.ceil(nHeads / m), layout.ny));
const out = [];
const byRow = new Map();
for (const h of layout.heads) { if (!byRow.has(h.row)) byRow.set(h.row, []); byRow.get(h.row).push(h); }
const rows = [...byRow.keys()].sort((a, b) => a - b);
const c0 = net.mainCol;
for (let r0 = 0; r0 + k <= rows.length; r0++) for (const mode of ['far', 'near']) for (const side of mode === 'far' ? [-1, 1] : [0]) {
const set = [];
for (let r = r0; r < r0 + k; r++) {
const hs = (byRow.get(rows[r]) || []).filter((h) => (side < 0 ? h.col < c0 : side > 0 ? h.col >= c0 : true));
const dist = (h) => Math.abs(h.col - (c0 - 0.5));
hs.sort((a, b) => (mode === 'far' ? dist(b) - dist(a) : dist(a) - dist(b)));
set.push(...hs.slice(0, m));
}
if (set.length) out.push({ id: `${mode}${side}@${r0}`, heads: set.map((h) => h.id) });
}
return out.length > 60 ? out.filter((_, i) => i % Math.ceil(out.length / 60) === 0) : out;
}The least flow and pressure every operating head must have.
q = density ÷ 1000 ÷ 60 × coverage; P = max((q ÷ k)², minimum head pressure)
The design density in mm/min times the coverage in m² gives m³/s with exact factors; the pressure is the square of one quotient.
const qMin = hz.density / 1000 / 60 * (opts.coverageM2 ?? 1); // m³/s per head
const pMinHead = Math.max((qMin / k) ** 2, (o.minHeadBar ?? 0) * 1e5);The pressure at every node and the flow in every pipe when the source is at a given pressure and each operating head discharges K √P.
P(b) = P(a) − friction × (length + equivalent length) − ρ g Δz
Head flows are computed from the current pressures, summed up the tree, and the pressures recomputed from the source down; this repeats until the largest change is under 1 Pa or 300 passes are done. A negative pressure is treated as 0 under the square root.
const pressures = (Ps) => {
const P = new Map([[net.source, Ps]]);
let q = new Map();
for (let it = 0; it < 300; it++) {
// flows from the heads' present pressures
const flow = new Map(edges.map((e) => [e.id, 0]));
const headQ = new Map();
for (const id of set) { const p = P.get(id) ?? Ps; const qh = k * sqrtP(p); headQ.set(id, qh); }
for (let i = order.length - 1; i >= 0; i--) {
const e = order[i];
let Q = headQ.get(e.b) || 0;
for (const c of kidsOf.get(e.b) || []) Q += flow.get(c.id);
flow.set(e.id, Q);
}
let delta = 0;
for (const e of order) {
const D = idOf(e.size), Q = flow.get(e.id);
const dp = hw(Q, D, o.C) * (e.length + eqLen(e, D, o.LD)) + RHO * G * e.dz;
const next = P.get(e.a) - dp;
delta = Math.max(delta, Math.abs(next - (P.get(e.b) ?? Ps)));
P.set(e.b, next);
}
q = flow;
if (delta < 1) break;
}
return { P, flow: q };
};The lowest source pressure at which every operating head has its minimum.
upper bound doubled from 200000 Pa until it works, then halved 40 times
Forty halvings narrow the bracket by a factor of 2⁴⁰, far below the 1 Pa tolerance of the inner iteration, so the result is good to about 1 Pa. If no pressure up to 50 MPa works the area is reported as failed.
const ok = (Ps) => { const { P } = pressures(Ps); for (const id of set) if ((P.get(id) ?? 0) < pMinHead - 1) return false; return true; };
let lo = 0, hi = 2e5;
while (!ok(hi) && hi < 5e7) hi *= 2;
if (hi >= 5e7) return { failed: true };
for (let i = 0; i < 40; i++) { const mid = (lo + hi) / 2; if (ok(mid)) hi = mid; else lo = mid; }
const { P, flow } = pressures(hi);Each head's flow, pressure and density, each pipe's flow, velocity and pressure drop, and the total flow.
density = q × 1000 × 60 ÷ coverage; v = Q ÷ (π D² ÷ 4); total = Σ head flows
The values are recomputed in 64-bit from the converged pressures with exact unit factors; nothing is rounded here. The total is the sum of the head flows in order.
const heads = [...set].map((id) => { const p = P.get(id), q = k * sqrtP(p); return { id, qLps: q * 1000, pBar: p / 1e5, density: q * 1000 * 60 / (opts.coverageM2 || 1) }; });
const eo = edges.map((e) => { const D = idOf(e.size), Q = flow.get(e.id); return { id: e.id, qLps: Q * 1000, v: D ? Q / (Math.PI * D * D / 4) : 0, dpPa: hw(Q, D, o.C) * (e.length + eqLen(e, D, o.LD)) + RHO * G * e.dz, size: e.size, kind: e.kind }; });
const qTot = heads.reduce((s, h) => s + h.qLps, 0);
return { psrcPa: hi, qLps: qTot, heads, edges: eo, maxV: Math.max(0, ...eo.filter((x) => x.qLps > 0).map((x) => x.v)), pMinHeadBar: pMinHead / 1e5 };The pressure the water supply gives at a flow, from the static and residual pressures of a flow test.
P = Ps − (Ps − Pr) × (Q ÷ Qr)^1.85
One 64-bit expression with the exponent 1.85. Without a test flow, or when the static pressure is not above the residual, the static pressure is returned unchanged.
/** Supply pressure (bar) at a flow (L/min) from the test: static Ps, residual Pr (bar) at Qr (L/min). */
export function supplyAt(sup, qLmin) {
if (!sup || !(sup.qr > 0) || !(sup.ps > sup.pr)) return sup?.ps ?? 0;
return sup.ps - (sup.ps - sup.pr) * (qLmin / sup.qr) ** 1.85;
}Takes the design area needing the most pressure, adds the hose allowance and margin, and compares with the supply.
Q total = Q heads × 60 + hose; required = P source ÷ 100000 + margin; ok if available ≥ required − 10⁻⁹
Flows and pressures are converted with exact factors in 64-bit. The comparison with the supply allows 1e-9 bar, and the least density among the heads is reported unrounded.
/** The whole calculation: size, try the candidate areas, take the one needing most pressure, add the hose allowance, compare with the supply. */
export function design(layout, net, hazard, supply, opts = {}) {
const o = { ...opts };
const cover = layout.coveragePerHead;
const qHead = hazard.density * cover / 60; // L/s at the design density
sizeNetwork(net, qHead, o);
let gov = null, tried = 0;
for (const c of candidates(layout, net, hazard)) {
const r = solveArea(net, c.heads, hazard, { ...o, coverageM2: cover });
tried++;
if (r.failed) continue;
if (!gov || r.psrcPa > gov.psrcPa) gov = { ...r, candidate: c };
}
if (!gov) return { failed: true, tried };
const hose = hazard.hose || 0;
const qLmin = gov.qLps * 60 + hose;
const need = gov.psrcPa / 1e5 + (o.margin || 0);
const have = supply ? supplyAt(supply, qLmin) : null;
return { gov, tried, qHeadsLmin: gov.qLps * 60, hoseLmin: hose, qTotalLmin: qLmin, requiredBar: gov.psrcPa / 1e5, requiredWithMarginBar: need, availableBar: have, ok: have == null ? null : have >= need - 1e-9,
minDensity: Math.min(...gov.heads.map((h) => h.density)), designDensity: hazard.density, maxV: gov.maxV, sizes: Object.fromEntries(net.edges.map((e) => [e.id, e.size])) };
}Tests tools/bimcraft-test/eg-mep-test.mjs
The number of heads each way and their spacing, the largest the hazard's maximum spacing and coverage allow.
nx = ceil(W ÷ maxS − 10⁻⁹); ny likewise; raised until (W ÷ nx)(D ÷ ny) ≤ max coverage; sx = W ÷ nx
The room is first turned so its longest wall lies along X, using the wall's atan2 angle. The slack of 1e-9 keeps a width that is an exact multiple of the spacing from gaining a row through rounding; the coverage test allows 1e-6 mm².
let nx = Math.max(1, Math.ceil(W / maxS - 1e-9)), ny = Math.max(1, Math.ceil(D / maxS - 1e-9));
for (let guard = 0; guard < 200 && (W / nx) * (D / ny) > maxA + 1e-6; guard++) { if (W / nx >= D / ny) nx++; else ny++; }
const sx = W / nx, sy = D / ny;Whether a point lies inside the room outline.
The even-odd ray test with one 64-bit division per crossing edge. A point exactly on an edge may fall on either side; no tolerance is applied.
function inside(poly, x, y) {
let c = false;
for (let i = 0, j = poly.length - 1; i < poly.length; j = i++) {
const a = poly[i], b = poly[j];
if ((a[1] > y) !== (b[1] > y) && x < ((b[0] - a[0]) * (y - a[1])) / (b[1] - a[1]) + a[0]) c = !c;
}
return c;
}Moves the head of a grid cell the room only partly fills to the middle of the part inside.
position = mean of the cell's 6 × 6 sample points that lie inside
The position is the average of up to 36 sample points, so it is an estimate of the centroid of the part inside. A cell with fewer than 3 samples inside gets no head.
// a cell the room only partly fills: the head goes to the middle of the part inside
let sxs = 0, sys = 0, n = 0;
for (let a = 0; a < 6; a++) for (let b = 0; b < 6; b++) { const px = x0 + sx * (i + (a + 0.5) / 6), py = y0 + sy * (j + (b + 0.5) / 6); if (inside(P, px, py)) { sxs += px; sys += py; n++; } }
if (n >= 3) { heads.push({ row: j, col: i, q: [sxs / n, sys / n], offGrid: true }); offGrid++; }The share of the room farther than 0.75 × the allowed spacing from every head.
uncovered = samples with no head within 0.75 × maxS ÷ samples inside the room
The room is sampled on a grid whose step is the larger of 150 mm and a quarter of the smaller spacing, so the share is a sampled estimate. It is reported as a percentage rounded to 1 decimal in the issue text and kept unrounded in the result.
let samples = 0, miss = 0;
const step = Math.max(150, Math.min(sx, sy) / 4);
for (let x = x0 + step / 2; x < x1; x += step) for (let y = y0 + step / 2; y < y1; y += step) {
if (!inside(P, x, y)) continue;
samples++;
if (!heads.some((h) => Math.hypot(h.q[0] - x, h.q[1] - y) <= 0.75 * maxS)) miss++; // no point farther than 0.75 × the allowed spacing from a head
}
const uncovered = samples ? miss / samples : 0;The length in meters and the rise of each pipe of the network from its two end points.
length = √(Δx² + Δy² + Δz²) ÷ 1000; dz = Δz ÷ 1000
One Math.hypot of the three differences in millimeters, divided by the exact 1000. Head positions pass through a rotation to the grid frame and back, each with its own rounding.
const edge = (a, b, kind, fit = {}) => {
const A = nodes.get(a).p, B = nodes.get(b).p;
edges.push({ id: a + '>' + b, a, b, length: Math.hypot(A[0] - B[0], A[1] - B[1], A[2] - B[2]) / 1000, dz: (B[2] - A[2]) / 1000, kind, fit });
};Tests None found by the review
The K-factor, coverage, design density and hazard class written on each placed head.
CoverageM2 = round(coverage × 100) ÷ 100
The coverage per head is rounded to 0.01 m² before it is stored; the K-factor and the design density are written as given.
w.upsertPset(id, 'CadCraft_Sprinkler', [['Head', 'IFCLABEL', h.id], ['KFactor', 'IFCREAL', k], ['CoverageM2', 'IFCREAL', Math.round(coverageM2 * 100) / 100], ['DesignDensity', 'IFCREAL', hazard.density], ['Hazard', 'IFCLABEL', hazard.id]]);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.