code wiki / _hdl_build / _f64_sincos_gate_authored.nx
_f64_sincos_gate_authored.nx source
↩ module page · 147 lines · 5876 B
1// _f64_sincos_gate_authored.nx -- ME1 sin/cos gate: TEAM-AUTHORED kernel vs the
2// LIVE sovereign bigfloat oracle (no stored vector files; expected values are
3// computed in-process by nx_bigfloat120_trig). Pure Nishi end to end.
4// 300 xorshift points across the v1 domain |x| < 2^20 + edge raws; ULP measured
5// via the signed ordered map. GREEN iff max ULP <= 1 on both functions.
6// Markers: SCBAD fn= x= exp= got= ulp= / SIN-ULP .../ COS-ULP ...
7// SINCOS-GATE verdict=GREEN|RED
8
9import "nx_syscalls.nx"
10import "nx_f64.nx"
11import "nx_f64_div.nx"
12import "nx_f64_sqrt.nx"
13import "nx_f64_cvt.nx"
14import "nx_bigfloat120.nx"
15import "nx_bigfloat120_div.nx"
16import "nx_bigfloat120_trig.nx"
17import "_pe_f64sincos.nx"
18
19func scg_puts(s: *u8) -> i64 { var n: i64 = 0; while s[n] != (0 as u8) { n = n + 1 } sys_write(1, s, n); return 0 }
20func scg_putn(v: i64) -> i64 { let bb: *u8 = sys_mmap(28); var m: i64 = v; if m < 0 { m = 0 - m; sys_write(1, "-" as *u8, 1) }; let t: *u8 = sys_mmap(28); var k: i64 = 0; if m == 0 { t[0] = 48; k = 1 }; while m > 0 { t[k] = 48 + (m % 10); m = m / 10; k = k + 1 }; var i: i64 = 0; while i < k { bb[i] = t[k-1-i]; i = i + 1 }; sys_write(1, bb, k); return 0 }
21func scg_puthex(v: i64) -> i64 { let bb: *u8 = sys_mmap(20); var i: i64 = 0; while i < 16 { let nib: i64 = (v >> ((15 - i) * 4)) & 15; if nib < 10 { bb[i] = 48 + nib } else { bb[i] = 55 + nib } i = i + 1 } sys_write(1, bb, 16); return 0 }
22
23func scg_rng(state: *i64) -> i64 {
24 var s: i64 = state[0]
25 s = s ^ ((s >> 12) & 0x000FFFFFFFFFFFFF)
26 s = s ^ (s << 25)
27 s = s ^ ((s >> 27) & 0x0000001FFFFFFFFF)
28 state[0] = s
29 return s * 2685821657736338717
30}
31
32func scg_ord(v: i64) -> i64 {
33 if v >= 0 { return v }
34 return (1 << 63) - v
35}
36
37// returns ulp distance; NaN-expected handled by caller
38func scg_ulp(got: i64, want: i64) -> i64 {
39 if got == want { return 0 }
40 var u: i64 = scg_ord(got) - scg_ord(want)
41 if u < 0 { u = 0 - u }
42 return u
43}
44
45func main() -> i64 {
46 let xs: *i64 = sys_mmap(8 * 512) as *i64
47 var n: i64 = 0
48 // edges: tiny, near pi/2 multiples, near domain edge, subnormals
49 xs[n] = 0x3FF0000000000000; n = n + 1 // 1.0
50 xs[n] = 0xBFF0000000000000; n = n + 1 // -1.0
51 xs[n] = 0x3FF921FB54442D18; n = n + 1 // ~pi/2
52 xs[n] = 0x400921FB54442D18; n = n + 1 // ~pi
53 xs[n] = 0x401921FB54442D18; n = n + 1 // ~2pi
54 xs[n] = 0xC012D97C7F3321D2; n = n + 1 // ~-3pi/2
55 xs[n] = 0x3E50000000000000; n = n + 1 // 2^-26
56 xs[n] = 0x0000000000000001; n = n + 1 // min subnormal
57 xs[n] = 0x800FFFFFFFFFFFFF; n = n + 1 // -max subnormal
58 xs[n] = 0x412FFFFFFFFFFFFF; n = n + 1 // just under 2^20... (2^18.99)
59 xs[n] = 0x4130000000000000; n = n + 1 // 2^18
60 xs[n] = 0x40F869F000000000; n = n + 1 // ~99999
61 let st: *i64 = sys_mmap(16) as *i64
62 st[0] = 31415926535897932
63 var j: i64 = 0
64 while j < 300 {
65 var raw: i64 = scg_rng(st)
66 var efr: i64 = 0
67 let band: i64 = raw & 3
68 if band == 0 { efr = 1015 + (raw & 7) } // near 1
69 if band == 1 { efr = 1023 + (raw & 15) } // 1 .. 2^16
70 if band == 2 { efr = 950 + (raw & 63) } // small
71 if band == 3 { efr = 1033 + (raw & 7) } // 2^10 .. 2^17
72 raw = (efr << 52) | (scg_rng(st) & 0x000FFFFFFFFFFFFF)
73 if (scg_rng(st) & 1) == 1 { raw = raw | (1 << 63) }
74 xs[n] = raw
75 n = n + 1
76 j = j + 1
77 }
78
79 var s_exact: i64 = 0
80 var s_u1: i64 = 0
81 var s_worse: i64 = 0
82 var s_max: i64 = 0
83 var c_exact: i64 = 0
84 var c_u1: i64 = 0
85 var c_worse: i64 = 0
86 var c_max: i64 = 0
87 var i: i64 = 0
88 while i < n {
89 let x: i64 = xs[i]
90 let sw: i64 = bf_sin_f64(x)
91 let sg: i64 = nx_f64_sin(x)
92 let su: i64 = scg_ulp(sg, sw)
93 if su == 0 { s_exact = s_exact + 1 }
94 if su == 1 { s_u1 = s_u1 + 1 }
95 if su > 1 {
96 s_worse = s_worse + 1
97 if s_worse <= 20 {
98 scg_puts("SCBAD fn=sin x=" as *u8); scg_puthex(x)
99 scg_puts(" exp=" as *u8); scg_puthex(sw)
100 scg_puts(" got=" as *u8); scg_puthex(sg)
101 scg_puts(" ulp=" as *u8); scg_putn(su)
102 scg_puts("\n" as *u8)
103 }
104 }
105 if su > s_max { s_max = su }
106 let cw: i64 = bf_cos_f64(x)
107 let cg: i64 = nx_f64_cos(x)
108 let cu: i64 = scg_ulp(cg, cw)
109 if cu == 0 { c_exact = c_exact + 1 }
110 if cu == 1 { c_u1 = c_u1 + 1 }
111 if cu > 1 {
112 c_worse = c_worse + 1
113 if c_worse <= 20 {
114 scg_puts("SCBAD fn=cos x=" as *u8); scg_puthex(x)
115 scg_puts(" exp=" as *u8); scg_puthex(cw)
116 scg_puts(" got=" as *u8); scg_puthex(cg)
117 scg_puts(" ulp=" as *u8); scg_putn(cu)
118 scg_puts("\n" as *u8)
119 }
120 }
121 if cu > c_max { c_max = cu }
122 i = i + 1
123 }
124
125 scg_puts("SIN-ULP total=" as *u8); scg_putn(n)
126 scg_puts(" exact=" as *u8); scg_putn(s_exact)
127 scg_puts(" ulp1=" as *u8); scg_putn(s_u1)
128 scg_puts(" worse=" as *u8); scg_putn(s_worse)
129 scg_puts(" max=" as *u8); scg_putn(s_max)
130 scg_puts("\n" as *u8)
131 scg_puts("COS-ULP total=" as *u8); scg_putn(n)
132 scg_puts(" exact=" as *u8); scg_putn(c_exact)
133 scg_puts(" ulp1=" as *u8); scg_putn(c_u1)
134 scg_puts(" worse=" as *u8); scg_putn(c_worse)
135 scg_puts(" max=" as *u8); scg_putn(c_max)
136 scg_puts("\n" as *u8)
137 if s_worse == 0 {
138 if c_worse == 0 {
139 scg_puts("SINCOS-GATE verdict=GREEN\n" as *u8)
140 return 0
141 }
142 }
143 scg_puts("SINCOS-GATE verdict=RED\n" as *u8)
144 var rc: i64 = s_worse + c_worse
145 if rc > 100 { rc = 100 }
146 return rc
147}