code wiki / _hdl_build / nx_mathspec_tan.nx
nx_mathspec_tan.nx source
↩ module page · 171 lines · 6005 B
1// nx_mathspec_tan.nx -- SOVEREIGN SPEC EMITTER for the TAN kernel: every
2// constant from the 120-bit bigfloat (trig + tan oracle organs). No python,
3// no minimax: the small-residual polynomial uses the EXACT tangent-number
4// Taylor family tan d = d + sum_{n>=2} T_n d^(2n-1)/(2n-1)! with the integer
5// tangent numbers T_2..T_8 = 2 16 272 7936 353792 22368256 1903757312
6// (OEIS A000182, mathematical integers like k!); each coefficient divides by
7// 2..(2n-1) one factor at a time (bf_div_small bound). |d| <= 1/16 -> 7
8// correction terms truncate below 2^-57 relative.
9//
10// Kernel scheme (accurate-table + dd): pi/2 Cody-Waite 3-part reduce (same
11// 21-zeroed HI/MID construction as the sincos spec) -> r dd; anchors a_j = j/8
12// (j = 0..6 covers [0, pi/4]); tan(a+d) = T + t*sec2(a)/(1 - T*t) with T and
13// sec2 carried hi/lo from this spec and the quotient refined in dd; odd
14// quadrants reconstruct -1/tan via a dd reciprocal.
15//
16// Layout: c[0]=one c[1]=half c[2]=2/pi c[3]=HI c[4]=MID c[5]=LO c[6]=eight
17// c[7..13] anchors j/8 c[14..20] T_hi[j] = fl(tan(j/8))
18// c[21..27] T_lo[j] (signed) c[28..34] S_hi[j] = fl(1+tan^2(j/8))
19// c[35..41] S_lo[j] (signed) c[42..48] Q coeffs hi-first (T_8.. down to T_2)
20// license_tier: ORIGINAL
21import "nx_syscalls.nx"
22import "nx_pattern_emit.nx"
23import "nx_bigfloat120.nx"
24import "nx_bigfloat120_div.nx"
25import "nx_bigfloat120_trig.nx"
26const K_MAGIC_7936: i64 = 7936
27const K_MAGIC_353792: i64 = 353792
28const K_MAGIC_22368256: i64 = 22368256
29const K_MAGIC_1903757312: i64 = 1903757312
30
31func mtn_wlit(fd: i64, v: i64) -> i64 {
32 if v < 0 {
33 pe_w(fd, "0 - " as *u8)
34 pe_wn(fd, 0 - v)
35 return 0
36 }
37 pe_wn(fd, v)
38 return 0
39}
40
41func mtn_emit_row(fd: i64, idx: i64, v: i64) -> i64 {
42 pe_w(fd, " c[" as *u8); pe_wn(fd, idx); pe_w(fd, "] = " as *u8)
43 mtn_wlit(fd, v)
44 pe_w(fd, "\n" as *u8)
45 return 0
46}
47
48// split a positive bigfloat into f64 hi + signed f64 lo
49func mtn_emit_hilo(fd: i64, idxhi: i64, idxlo: i64, a: *i64) -> i64 {
50 if bf_is_zero(a) == 1 {
51 mtn_emit_row(fd, idxhi, 0)
52 mtn_emit_row(fd, idxlo, 0)
53 return 0
54 }
55 let hi: i64 = bf_to_f64(a, 0, 0)
56 let hib: *i64 = bf_new()
57 bf_set_f64(hib, hi)
58 let rem: *i64 = bf_new()
59 var losig: i64 = 0
60 if bf_cmp(a, hib) >= 0 {
61 bf_sub(rem, a, hib)
62 } else {
63 bf_sub(rem, hib, a)
64 losig = 1
65 }
66 mtn_emit_row(fd, idxhi, hi)
67 mtn_emit_row(fd, idxlo, bf_to_f64(rem, losig, 0))
68 return 0
69}
70
71func main() -> i64 {
72 let fd: i64 = sys_openat_wr("runtime/_hdl_build/_pm_tan_spec.nx" as *u8, 0x1a4)
73 if fd < 0 { sys_exit(2) }
74 pe_w(fd, "// _pm_tan_spec.nx -- GENERATED by nx_mathspec_tan (SOVEREIGN bigfloat,\n" as *u8)
75 pe_w(fd, "// no python). DO NOT HAND-EDIT. Layout: see nx_mathspec_tan.nx header.\n" as *u8)
76 pe_w(fd, "import \"nx_syscalls.nx\"\n\nfunc pm_tan_spec_fill(c: *i64) -> i64 {\n" as *u8)
77
78 let one: *i64 = bf_new()
79 bf_set_int(one, 1)
80 mtn_emit_row(fd, 0, bf_to_f64(one, 0, 0))
81 let half: *i64 = bf_new()
82 bf_copy(half, one)
83 half[0] = half[0] - 1
84 mtn_emit_row(fd, 1, bf_to_f64(half, 0, 0))
85
86 // 2/pi + 3-part pi/2 split (same construction as nx_mathspec_sincos)
87 let pi: *i64 = bf_new()
88 bf_pi(pi)
89 let pio2: *i64 = bf_new()
90 bf_copy(pio2, pi)
91 pio2[0] = pio2[0] - 1
92 let two: *i64 = bf_new()
93 bf_set_int(two, 2)
94 let topi: *i64 = bf_new()
95 bf_div(topi, two, pi)
96 mtn_emit_row(fd, 2, bf_to_f64(topi, 0, 0))
97 let m21: i64 = (1 << 21) - 1
98 let hi_bits: i64 = bf_to_f64(pio2, 0, 0) & (0 - 1 - m21)
99 let hib: *i64 = bf_new()
100 bf_set_f64(hib, hi_bits)
101 let rem1: *i64 = bf_new()
102 bf_sub(rem1, pio2, hib)
103 let mid_bits: i64 = bf_to_f64(rem1, 0, 0) & (0 - 1 - m21)
104 let midb: *i64 = bf_new()
105 bf_set_f64(midb, mid_bits)
106 let rem2: *i64 = bf_new()
107 bf_sub(rem2, rem1, midb)
108 mtn_emit_row(fd, 3, hi_bits)
109 mtn_emit_row(fd, 4, mid_bits)
110 mtn_emit_row(fd, 5, bf_to_f64(rem2, 0, 0))
111 let eightb: *i64 = bf_new()
112 bf_set_int(eightb, 8)
113 mtn_emit_row(fd, 6, bf_to_f64(eightb, 0, 0))
114
115 // anchors + tan/sec2 hi-lo tables: a_j = j/8 <= 0.75 < pi/4, so tan(a_j)
116 // comes straight from the blessed first-octant series organs.
117 var j: i64 = 0
118 while j <= 6 {
119 if j == 0 {
120 mtn_emit_row(fd, 7, 0)
121 mtn_emit_row(fd, 14, 0)
122 mtn_emit_row(fd, 21, 0)
123 mtn_emit_hilo(fd, 28, 35, one) // sec^2(0) = 1
124 } else {
125 let aj: *i64 = bf_new()
126 bf_set_int(aj, j)
127 aj[0] = aj[0] - 3
128 mtn_emit_row(fd, 7 + j, bf_to_f64(aj, 0, 0))
129 let sv: *i64 = bf_new()
130 let cv: *i64 = bf_new()
131 bf_sin_r(sv, aj)
132 bf_cos_r(cv, aj)
133 let tj: *i64 = bf_new()
134 bf_div(tj, sv, cv)
135 mtn_emit_hilo(fd, 14 + j, 21 + j, tj)
136 let t2: *i64 = bf_new()
137 bf_mul(t2, tj, tj)
138 let sj: *i64 = bf_new()
139 bf_add(sj, one, t2)
140 mtn_emit_hilo(fd, 28 + j, 35 + j, sj)
141 }
142 j = j + 1
143 }
144
145 // Q coeffs hi-first: coeff_n = T_n/(2n-1)!, n = 8 down to 2
146 let tn: *i64 = sys_mmap(8 * 10) as *i64
147 tn[2] = 2; tn[3] = 16; tn[4] = 272; tn[5] = K_MAGIC_7936
148 tn[6] = K_MAGIC_353792; tn[7] = K_MAGIC_22368256; tn[8] = K_MAGIC_1903757312
149 var n: i64 = 8
150 var idx: i64 = 42
151 while n >= 2 {
152 let cf: *i64 = bf_new()
153 bf_set_int(cf, tn[n])
154 let tmp: *i64 = bf_new()
155 var dvr: i64 = 2
156 while dvr <= 2 * n - 1 {
157 bf_div_small(tmp, cf, dvr)
158 bf_copy(cf, tmp)
159 dvr = dvr + 1
160 }
161 mtn_emit_row(fd, idx, bf_to_f64(cf, 0, 0))
162 idx = idx + 1
163 n = n - 1
164 }
165
166 pe_w(fd, " return 49\n}\n" as *u8)
167 sys_close(fd)
168 pe_w(1, "SPEC EMITTED: _pm_tan_spec.nx (49 consts, all from sovereign bigfloat)\n" as *u8)
169 sys_exit(0)
170 return 0
171}