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}