nx_f32_cvt.nx source
↩ module page · 304 lines · 11950 B
1// nx_f32_cvt.nx -- IEEE 754 binary32 conversions (bits-up).
2//
3// L4.5 brick that bridges the L1 fixed-point substrate (Q10/Q14) and
4// the L2 f16 layer into L4 f32. All conversions are clean-room from
5// the IEEE 754-2019 spec; no libm, no compiler-builtin lowerings.
6//
7// Operations in v1:
8// nx_f16_to_f32(raw_f16) -> raw_f32 widening with subnormal-
9// to-normal renormalization
10// nx_i32_to_f32(value_i32) -> raw_f32 small-int (round-to-
11// nearest-even for >24 bits)
12// nx_q10_to_f32(value_q10) -> raw_f32 Q10 (* 1024) into f32
13// nx_q14_to_f32(value_q14) -> raw_f32 Q14 (* 16384) into f32
14//
15// Composes with nx_f32_mul / nx_f32_add / nx_f32_sub / nx_f32_div /
16// nx_f32_sqrt (all shipped previously in this L4 tower) to produce a
17// real-model dequant pipeline emitting IEEE 754 outputs that can be
18// bit-exact compared against ggml's reference dequant.
19//
20// Reference: IEEE 754-2019 binary16 / binary32 layouts; standard
21// widening + signed-int-to-float algorithms. No code borrowed.
22//
23// genealogy_id: ieee754_2019_widening + signed_int_to_float_canonical
24// lineage_id: substrate_f32_cvt_v1
25
26import "nx_syscalls.nx"
27import "nx_tier.nx"
28import "nx_f16.nx"
29import "nx_f32.nx"
30const NX_MAGIC_8388608: i64 = 8388608
31const NX_MAGIC_8388607: i64 = 8388607
32
33// ===== f16 -> f32 widening =========================================
34//
35// f16: 1 sign + 5 exp (bias 15) + 10 mantissa
36// f32: 1 sign + 8 exp (bias 127) + 23 mantissa
37//
38// Bias adjustment: f32_exp = f16_exp - 15 + 127 = f16_exp + 112
39// Mantissa: shift left 13 (10 -> 23 bits)
40//
41// Special cases:
42// zero -> signed zero
43// subnormal -> normal f32 (since 2^-14 well within f32 normal
44// range -126..127)
45// inf -> f32 inf
46// NaN -> f32 canonical qNaN
47
48func nx_f16_to_f32(raw: i64) -> i64 {
49 let sign: i64 = (raw >> 15) & 1
50 let exp_f16: i64 = (raw >> 10) & 31
51 let mant_f16: i64 = raw & 0x3FF
52
53 if exp_f16 == 0 {
54 if mant_f16 == 0 {
55 return sign << 31 // signed zero
56 }
57 // Subnormal f16: renormalize. Find leading 1 bit in mantissa,
58 // shift up to position 10 (implicit-1 in normal form), then
59 // drop leading 1 + widen.
60 var m: i64 = mant_f16
61 var exp_adj: i64 = 0
62 var iter: nx_int = 0
63 while iter < 11 {
64 if m < 0x400 {
65 m = m << 1
66 exp_adj = exp_adj + 1
67 }
68 iter = iter + 1
69 }
70 // Now m has implicit-1 at bit 10.
71 // Subnormal f16 value = mant * 2^-24; after exp_adj left-shifts to seat the
72 // implicit 1 at bit 10, the real exponent is -14 - exp_adj (NOT -13; the old
73 // "+1" here doubled every subnormal f16 -> the Q6_K ffn_down 2x forward bug).
74 // f32 biased = (-14 - exp_adj) + 127 = 113 - exp_adj
75 let f32_exp: i64 = 113 - exp_adj
76 let f32_mant: i64 = (m & 0x3FF) << 13
77 return (sign << 31) | (f32_exp << 23) | f32_mant
78 }
79 if exp_f16 == 31 {
80 if mant_f16 == 0 {
81 return (sign << 31) | 0x7F800000 // ±inf
82 }
83 return 0x7FC00000 // canonical qNaN
84 }
85
86 // Normal f16 -> normal f32.
87 let f32_exp: i64 = exp_f16 + 112
88 let f32_mant: i64 = mant_f16 << 13
89 return (sign << 31) | (f32_exp << 23) | f32_mant
90}
91
92// ===== i32 -> f32 ==================================================
93//
94// Convert signed i32 to f32. For |x| <= 2^24 the result is exact;
95// for larger magnitudes round-to-nearest-even on the dropped bits.
96//
97// Algorithm:
98// 1. Extract sign, take absolute value.
99// 2. Find position of leading 1 bit ("top").
100// 3. Biased exp = top + 127.
101// 4. Mantissa = bits below leading 1, shifted into 23-bit field.
102// 5. For top > 23: dropped bits become guard+sticky for rounding.
103//
104// For Q4_K dequant usage (sc, m in [0..63], q4 in [0..15]),
105// top <= 5 -- always exact.
106
107func nx_i32_to_f32(value: i64) -> i64 {
108 if value == 0 { return 0 }
109 var sign: i64 = 0
110 var v: i64 = value
111 if v < 0 {
112 sign = 1
113 v = 0 - v
114 }
115
116 // Find leading 1 bit position (0..31) via top-down scan.
117 var top: i64 = 31
118 var found: nx_int = 0
119 var ti: nx_int = 0
120 while ti < 32 {
121 if found == 0 {
122 if ((v >> top) & 1) == 1 { found = 1 }
123 else { top = top - 1 }
124 }
125 ti = ti + 1
126 }
127 // top is now the position of the leading 1 (0..31).
128
129 let f32_exp: i64 = top + 127
130 var f32_mant: i64 = 0
131
132 if top <= 23 {
133 // All bits fit below position 23; shift left to fill.
134 let shift_up: i64 = 23 - top
135 f32_mant = (v << shift_up) & 0x7FFFFF
136 return (sign << 31) | (f32_exp << 23) | f32_mant
137 }
138
139 // top > 23: need to round.
140 let shift_down: i64 = top - 23
141 let guard_mask: i64 = (1 << shift_down) - 1
142 let lost: i64 = v & guard_mask
143 f32_mant = (v >> shift_down) & 0x7FFFFF
144 let halfway: i64 = 1 << (shift_down - 1)
145 var round_up: i64 = 0
146 if lost > halfway { round_up = 1 }
147 if lost == halfway {
148 if (f32_mant & 1) == 1 { round_up = 1 }
149 }
150 if round_up == 1 {
151 f32_mant = f32_mant + 1
152 // Handle mantissa overflow from rounding.
153 if f32_mant >= 0x800000 {
154 f32_mant = 0
155 // exp_out + 1 -- but we're returning packed; rebuild:
156 let f32_exp_bumped: i64 = f32_exp + 1
157 return (sign << 31) | (f32_exp_bumped << 23) | f32_mant
158 }
159 }
160 return (sign << 31) | (f32_exp << 23) | f32_mant
161}
162
163// ===== Q10 -> f32 ==================================================
164//
165// Q10 = value_real * 1024. Convert by dividing by 1024 via the
166// i32->f32 path then multiplying by 2^-10 (= adjusting exponent
167// down by 10). Special-case zero.
168
169func nx_q10_to_f32(q10: i64) -> i64 {
170 if q10 == 0 { return 0 }
171 let as_f32: i64 = nx_i32_to_f32(q10)
172 // as_f32 represents q10 as a real number. We want q10 / 1024.
173 // Adjust the biased exponent by -10 (= divide by 2^10).
174 let exp_field: i64 = (as_f32 >> 23) & 0xFF
175 let new_exp: i64 = exp_field - 10
176 if new_exp <= 0 { return 0 } // underflow -> zero (v1)
177 let cleared: i64 = as_f32 & 0x807FFFFF // clear exp field
178 return cleared | (new_exp << 23)
179}
180
181// ===== f32 -> Q10 ==================================================
182//
183// The inverse of nx_q10_to_f32, and the primitive that was MISSING: the fused
184// integer Q4_K dot (nx_q4k_dot_row_col) consumes a dense Q10 column, but the
185// forward path carries activations as f32 bits, so nothing could feed it.
186// Without this the ~10x fused kernel is unreachable from the forward -- the
187// conversion, not the kernel, was the adoption blocker.
188//
189// value = (-1)^s * 2^(e-127) * (1 + m/2^23), and Q10 = round(value * 1024).
190// Let full_mant = 2^23 | m (a 24-bit integer). Then
191// value * 1024 = (-1)^s * full_mant * 2^(e - 127 + 10 - 23) = full_mant * 2^(e-140).
192// So a single shift by (e-140) does the whole conversion; no division, no
193// emulated-f32 arithmetic. Right shifts round HALF AWAY FROM ZERO on the
194// magnitude, matching nx_q4km_q20_to_q10's convention so a round-trip through
195// the fused dot does not drift.
196//
197// DECLARED BOUNDS (this is fixed point -- the limits are real, not theoretical):
198// |value| < 2^-11 -> returns 0 (below Q10 resolution; genuine underflow)
199// |value| > 2^48 -> SATURATES to +/-NX_Q10_SAT rather than wrapping i64.
200// Wrapping would silently flip the sign of a huge
201// activation; saturation keeps the sign and magnitude
202// ordering intact. Real activations are nowhere near
203// this, so hitting it means something upstream is wrong.
204// exp==0 (zero/denormal) and exp==255 (inf/nan) -> 0. This substrate has no
205// traps; a NaN reaching a dot product would poison a whole row silently.
206
207const NX_Q10_SAT: i64 = 4611686018427387904 // 2^62, headroom for the Q34 accumulate
208
209func nx_f32_to_q10(f32bits: i64) -> i64 {
210 let exp_field: i64 = (f32bits >> 23) & 0xFF
211 if exp_field == 0 { return 0 }
212 if exp_field == 255 { return 0 }
213 let sign: i64 = (f32bits >> 31) & 1
214 let full_mant: i64 = NX_MAGIC_8388608 | (f32bits & NX_MAGIC_8388607) // 2^23 | mantissa
215 let shift: i64 = exp_field - 140
216 var mag: i64 = 0
217 if shift >= 0 {
218 if shift > 38 {
219 if sign == 1 { return 0 - NX_Q10_SAT }
220 return NX_Q10_SAT
221 }
222 mag = full_mant << shift
223 } else {
224 let rs: i64 = 0 - shift
225 if rs >= 32 { return 0 }
226 mag = (full_mant + (1 << (rs - 1))) >> rs
227 }
228 if sign == 1 { return 0 - mag }
229 return mag
230}
231
232// ===== Q20 <-> f32 (the fused GEMM's activation format) ============
233//
234// WHY Q20 AND NOT Q10. The fused Q4_K GEMM's ONLY source of error versus the f32 path is the
235// resolution of the activation fixed-point format. Measured on real blk.0.attn_q weights, Q10
236// activations (1/1024) cost ~1.37% max relative deviation -- inside the 2% KAT band, but far too much
237// to compound across ~36 layers. Widening the activation to Q20 costs NOTHING: it is the same i64
238// multiply and the same instruction count. Only the accumulator range changes, and it has room --
239// Q24 weight x Q20 activation = Q44, and over k=4096 terms a realistic |value|<=8 stays near 2^52,
240// comfortably inside i64. u2605 PRECISION HERE WAS FREE AND WE WERE PAYING FOR IT ANYWAY.
241//
242// Recovery reuses nx_q4km_q20_to_q10, which is a rounded >>24 -- applied to a Q44 accumulator it
243// yields Q20 (that function's own header notes the name is historical and it is a >>24 normalisation).
244
245func nx_q20_to_f32(q20: i64) -> i64 {
246 if q20 == 0 { return 0 }
247 let as_f32: i64 = nx_i32_to_f32(q20)
248 let exp_field: i64 = (as_f32 >> 23) & 0xFF
249 let new_exp: i64 = exp_field - 20
250 if new_exp <= 0 { return 0 }
251 let cleared: i64 = as_f32 & 0x807FFFFF
252 return cleared | (new_exp << 23)
253}
254
255// value * 2^20 = full_mant * 2^(e - 127 + 20 - 23) = full_mant * 2^(e - 130).
256// Same construction and the same declared bounds as nx_f32_to_q10, shifted 10 bits finer.
257func nx_f32_to_q20(f32bits: i64) -> i64 {
258 let exp_field: i64 = (f32bits >> 23) & 0xFF
259 if exp_field == 0 { return 0 }
260 if exp_field == 255 { return 0 }
261 let sign: i64 = (f32bits >> 31) & 1
262 let full_mant: i64 = NX_MAGIC_8388608 | (f32bits & NX_MAGIC_8388607)
263 let shift: i64 = exp_field - 130
264 var mag: i64 = 0
265 if shift >= 0 {
266 if shift > 38 {
267 if sign == 1 { return 0 - NX_Q10_SAT }
268 return NX_Q10_SAT
269 }
270 mag = full_mant << shift
271 } else {
272 let rs: i64 = 0 - shift
273 if rs >= 32 { return 0 }
274 mag = (full_mant + (1 << (rs - 1))) >> rs
275 }
276 if sign == 1 { return 0 - mag }
277 return mag
278}
279
280// Q24 -> f32. The Q4_K iterator holds d/dmin as EXACT Q24 (its fields keep historical *_q10 names).
281// The i8-SIMD kernel needs them as f32 scalars to scale __f32_i8dot32's f32 result.
282func nx_q24_to_f32(q24: i64) -> i64 {
283 if q24 == 0 { return 0 }
284 let as_f32: i64 = nx_i32_to_f32(q24)
285 let exp_field: i64 = (as_f32 >> 23) & 0xFF
286 let new_exp: i64 = exp_field - 24
287 if new_exp <= 0 { return 0 }
288 let cleared: i64 = as_f32 & 0x807FFFFF
289 return cleared | (new_exp << 23)
290}
291
292// ===== Q14 -> f32 ==================================================
293//
294// Q14 = value_real * 16384. Adjust biased exponent down by 14.
295
296func nx_q14_to_f32(q14: i64) -> i64 {
297 if q14 == 0 { return 0 }
298 let as_f32: i64 = nx_i32_to_f32(q14)
299 let exp_field: i64 = (as_f32 >> 23) & 0xFF
300 let new_exp: i64 = exp_field - 14
301 if new_exp <= 0 { return 0 }
302 let cleared: i64 = as_f32 & 0x807FFFFF
303 return cleared | (new_exp << 23)
304}