nx_pow10.nx source
↩ module page · 142 lines · 5927 B
1// nx_pow10.nx -- SHARED NUMERIC rung: base-10 exponential and logarithm in
2// pure integer arithmetic. Extracted because two unrelated lanes needed the
3// same primitive -- peptide ionisation fractions (Henderson-Hasselbalch) and
4// thermal death kinetics (D/z values) are both decade-log physics -- and a
5// physical primitive with two copies is a drift waiting to happen.
6//
7// HOW IT WORKS. A 10-entry mantissa table holds 10^0.0 .. 10^0.9 scaled by
8// 1000; any exponent splits into whole decades (repeated x10) plus a
9// fractional part read from the table. ilog10_milli is the same table run
10// backwards: normalise into [1, 10), find the bracketing pair, interpolate.
11//
12// INTERPOLATED, AND THE ERROR IS BOUNDED AND SIGNED. Nearest-entry lookup
13// alone is 0.1 log units granular = up to 12% error, which is far too coarse
14// for a food-safety lethality calculation. Linear interpolation between
15// table entries drops that to <=0.7%. Because 10^x is CONVEX, a chord lies
16// above the curve, so interpolation always returns a value >= the true one.
17// That sign is not incidental -- callers computing a lethality with a
18// NEGATIVE exponent (process cooler than reference) get a slightly
19// UNDER-stated result, which is the safe direction. A caller relying on a
20// positive exponent should carry margin; see nx_thermal_process, which does.
21//
22// ⚠THAT BOUND WAS WRONG UNTIL 2026-07-25 AND SAID <=0.3%. This file had no
23// gate of its own -- it was exercised only through its callers, so nothing
24// ever checked its self-description. nx_pow10_test now MEASURES the bound
25// via the functional equation 10^a*10^b = 10^(a+b), scanning the decade for
26// the worst case: 13 permil on the doubled identity, hence ~6.5 permil on a
27// single evaluation. That matches the closed form for a chord over a convex
28// function, (1+r)/2/sqrt(r) - 1 with r = 10^0.1, which is 0.664%. The old
29// figure understated the error by a factor of two.
30//
31// The one caller that reasons about this numerically (nx_thermal_process,
32// botulinum cook adequacy) was re-checked against the corrected figure and
33// its conclusion is UNCHANGED: 0.7% is still an order of magnitude inside
34// the ~19% margin the 3.0 min industry F0 target carries over the 2.52 min
35// 12D floor. The claim was wrong; the safety verdict it supported was not.
36//
37// Units: _q3 = value x 1000, _milli = value x 1000.
38//
39// Grounding: standard_base10_logarithm_identities (no citation needed beyond
40// the definition; the gate checks against exact known values).
41//
42// genealogy_id: integer_numerics + nishi_shared
43
44import "nx_syscalls.nx"
45const IP10_MAGIC_1259: i64 = 1259
46const IP10_MAGIC_1585: i64 = 1585
47const IP10_MAGIC_1995: i64 = 1995
48const IP10_MAGIC_2512: i64 = 2512
49const IP10_MAGIC_3162: i64 = 3162
50const IP10_MAGIC_3981: i64 = 3981
51const IP10_MAGIC_5012: i64 = 5012
52const IP10_MAGIC_6310: i64 = 6310
53const IP10_MAGIC_7943: i64 = 7943
54const IP10_MAGIC_10000: i64 = 10000
55const IP10_MAGIC_1000000: i64 = 1000000
56
57const IP10_INVALID: i64 = 0 - 1
58
59// Beyond +/-6 decades the caller is out of any physical range this serves,
60// and 10^6 x 1000 already sits comfortably inside i64. Saturate.
61const IP10_CLAMP_MILLI: i64 = 6000
62
63// 10^(t/10) x 1000 for t = 0..10. Entry 10 is the next decade, present so
64// interpolation has a right-hand bracket at the top of the range.
65func ipow10_mantissa_q3(tenth: i64) -> i64 {
66 if tenth == 0 { return 1000 }
67 if tenth == 1 { return IP10_MAGIC_1259 }
68 if tenth == 2 { return IP10_MAGIC_1585 }
69 if tenth == 3 { return IP10_MAGIC_1995 }
70 if tenth == 4 { return IP10_MAGIC_2512 }
71 if tenth == 5 { return IP10_MAGIC_3162 }
72 if tenth == 6 { return IP10_MAGIC_3981 }
73 if tenth == 7 { return IP10_MAGIC_5012 }
74 if tenth == 8 { return IP10_MAGIC_6310 }
75 if tenth == 9 { return IP10_MAGIC_7943 }
76 if tenth == 10 { return IP10_MAGIC_10000 }
77 return 1000
78}
79
80// 10^(d_milli/1000) x 1000, interpolated.
81func ipow10_q3(d_milli: i64) -> i64 {
82 var neg: i64 = 0
83 var k: i64 = 0
84 var d: i64 = d_milli
85 if d < 0 { neg = 1 }
86 if neg == 1 { d = 0 - d }
87 if d > IP10_CLAMP_MILLI { d = IP10_CLAMP_MILLI }
88 let decades: i64 = d / 1000
89 let rem: i64 = d - decades * 1000
90 let tenth: i64 = rem / 100
91 let frac: i64 = rem - tenth * 100
92 let lo: i64 = ipow10_mantissa_q3(tenth)
93 let hi: i64 = ipow10_mantissa_q3(tenth + 1)
94 let span: i64 = hi - lo
95 var v: i64 = lo + span * frac / 100
96 while k < decades {
97 v = v * 10
98 k = k + 1
99 }
100 if neg == 0 { return v }
101 if v <= 0 { return 0 }
102 return IP10_MAGIC_1000000 / v
103}
104
105// log10(v_q3/1000) x 1000. Input is the value scaled by 1000, so v_q3=1000
106// means 1.0 and returns 0. Non-positive input is REFUSED, not clamped: the
107// logarithm of zero or a negative is undefined and returning 0 would be a
108// lie that propagates silently into a divisor.
109func ilog10_milli(v_q3: i64) -> i64 {
110 if v_q3 <= 0 { return IP10_INVALID }
111 var decade: i64 = 0
112 var m: i64 = v_q3
113 while m >= IP10_MAGIC_10000 {
114 m = m / 10
115 decade = decade + 1
116 }
117 while m < 1000 {
118 m = m * 10
119 decade = decade - 1
120 }
121 // Find the bracketing pair. NishiLang has no break, and assigning the
122 // loop variable to escape would destroy the index the next lines need,
123 // so this runs all 10 iterations and latches the FIRST match.
124 var t: i64 = 0
125 var found: i64 = 0
126 var tenth: i64 = 9
127 while t < 10 {
128 if found == 0 {
129 let hi: i64 = ipow10_mantissa_q3(t + 1)
130 if m < hi { found = 1; tenth = t }
131 }
132 t = t + 1
133 }
134 let lo: i64 = ipow10_mantissa_q3(tenth)
135 let hi2: i64 = ipow10_mantissa_q3(tenth + 1)
136 let span: i64 = hi2 - lo
137 var frac: i64 = 0
138 if span > 0 { frac = (m - lo) * 100 / span }
139 let base: i64 = decade * 1000
140 let tpart: i64 = tenth * 100
141 return base + tpart + frac
142}