code wiki / (root) / nx_icbrt.nx

nx_icbrt.nx source

↩ module page · 98 lines · 3794 B

1// nx_icbrt.nx -- bits-up integer cube root + Q14 fixed-point wrapper. 2// 3// Returns floor(cube_root(n)) for n >= 0; 0 for n < 0. 4// 5// Algorithm: bit-by-bit, analogous to nx_isqrt. Iterate from the 6// highest possible bit of the cube root down; at each bit, try setting 7// it and check whether candidate^3 still fits within n. 8// 9// Why this exists (substrate-level motivation, 2026-05-19): 10// Q14 fixed-point cube root is needed for unit-exact bending mechanics 11// in nx_pillar_physics v2.1 (the EXCEED axis vs Orca/Bambu/Cura/Prusa 12// which use constant pillar diameter regardless of load). Square 13// pillar bending stress: σ = M × (s/2) / I where I = s⁴/12 yields 14// the textbook formula s = ∛(6M / σ_yield). No production slicer 15// implements this; no academic system implements it at slice-time. 16// 17// Cube root has uses beyond pillar sizing -- volumetric flow at the 18// extruder (radius scales as cbrt of flow per area), surface-area / 19// volume ratios for cooling models, etc. Ships as a general bits-up 20// math primitive that any future arc can compose. 21// 22// Overflow safety: for i64 inputs the cube root fits in 22 bits (since 23// cube_root(2^63) ≈ 2.1e21 < 2^22 = 4.2e21). At candidate ≈ 2^21, 24// candidate^2 ≈ 2^42 and the test candidate^3 ≤ n is rewritten as 25// candidate^2 ≤ n / candidate to keep all intermediates in i64. 26// 27// license_tier: ORIGINAL 28 29import "nx_syscalls.nx" 30import "nx_tier.nx" 31 32const NX_ICBRT_TOP_BIT: i64 = 21 33 34// Q14 scaling constant: 16384² = 268435456 = 2^28. 35// Used to lift a Q14 value into a regime where icbrt produces a Q14 36// result. See nx_icbrt_q14 derivation below. 37const NX_ICBRT_Q14_SHIFT: i64 = 268435456 38 39// Safe input cap for nx_icbrt_q14 before the internal × 2^28 overflows. 40// Signed i64 max = 2^63 - 1, so the largest x_q14 with x_q14 × 2^28 41// still positive is floor((2^63 - 1) / 2^28) = 2^35 - 1. 42// Using 2^35 exactly would produce 2^63 which overflows to negative. 43const NX_ICBRT_Q14_LIMIT: i64 = 34359738367 44 45// ===== integer cube root ========================================== 46// 47// Returns floor(cube_root(n)) for n >= 0. Negative input -> 0. 48// 49// Closed-form invariants: 50// icbrt(0) = 0 51// icbrt(1) = 1 52// icbrt(k³) = k for any nonneg k whose cube fits in i64 53// icbrt(k³ - 1) = k - 1 54// icbrt(k³ + 1) = k 55 56func nx_icbrt(n: nx_int) -> nx_int { 57 if n <= 0 { return 0 } 58 if n == 1 { return 1 } 59 var x: i64 = 0 60 var b: i64 = NX_ICBRT_TOP_BIT 61 while b >= 0 { 62 let one: i64 = 1 63 let bit: i64 = one << b 64 let candidate: i64 = x | bit 65 // candidate^3 ≤ n iff candidate^2 ≤ n / candidate (candidate > 0). 66 // Rewritten form keeps all intermediates within i64. 67 let csq: i64 = candidate * candidate 68 let lim: i64 = n / candidate 69 if csq <= lim { x = candidate } 70 b = b - 1 71 } 72 return x 73} 74 75// ===== Q14 cube root =============================================== 76// 77// Input: x_q14 = x_real × 16384 (Q14 representation of x_real) 78// Output: r_q14 = cube_root(x_real) × 16384 (Q14 of cube root) 79// 80// Derivation: 81// r_real = cube_root(x_real) 82// r_q14 = r_real × 16384 83// r_q14³ = r_real³ × 16384³ 84// = x_real × 16384³ 85// = (x_q14 / 16384) × 16384³ 86// = x_q14 × 16384² 87// So: r_q14 = icbrt(x_q14 × 16384²) = icbrt(x_q14 × NX_ICBRT_Q14_SHIFT). 88// 89// Overflow: x_q14 × 2^28 must fit in i64 -> x_q14 ≤ 2^35 = NX_ICBRT_Q14_LIMIT. 90// Inputs beyond this saturate at the limit (cube root of saturated value). 91 92func nx_icbrt_q14(x_q14: nx_int) -> nx_int { 93 if x_q14 <= 0 { return 0 } 94 var capped: i64 = x_q14 95 if capped > NX_ICBRT_Q14_LIMIT { capped = NX_ICBRT_Q14_LIMIT } 96 let scaled: i64 = capped * NX_ICBRT_Q14_SHIFT 97 return nx_icbrt(scaled) 98}