c/
math_log.rs

1// Translated from musl: log.c, logf.c, log10.c, log10f.c, log2.c, log2f.c
2// and their associated *_data.c files.
3
4// ============================================================
5// log (double precision) data
6// ============================================================
7
8const LOG_TABLE_BITS: u32 = 7;
9const LOG_N: usize = 1 << LOG_TABLE_BITS;
10const LN2HI: f64 = asdouble(0x3fe62e42fefa3800);
11const LN2LO: f64 = asdouble(0x3d2ef35793c76730);
12
13const LOG_POLY1: [f64; 11] = [
14    asdouble(0xbfe0000000000000),
15    asdouble(0x3fd5555555555577),
16    asdouble(0xbfcffffffffffdcb),
17    asdouble(0x3fc999999995dd0c),
18    asdouble(0xbfc55555556745a7),
19    asdouble(0x3fc24924a344de30),
20    asdouble(0xbfbfffffa4423d65),
21    asdouble(0x3fbc7184282ad6ca),
22    asdouble(0xbfb999eb43b068ff),
23    asdouble(0x3fb78182f7afd085),
24    asdouble(0xbfb5521375d145cd),
25];
26
27const LOG_POLY: [f64; 5] = [
28    asdouble(0xbfe0000000000001),
29    asdouble(0x3fd555555551305b),
30    asdouble(0xbfcfffffffeb4590),
31    asdouble(0x3fc999b324f10111),
32    asdouble(0xbfc55575e506c89f),
33];
34
35const LOG_TAB: [(u64, u64); 128] = [
36    (0x3ff734f0c3e0de9f, 0xbfd7cc7f79e69000),
37    (0x3ff713786a2ce91f, 0xbfd76feec20d0000),
38    (0x3ff6f26008fab5a0, 0xbfd713e31351e000),
39    (0x3ff6d1a61f138c7d, 0xbfd6b85b38287800),
40    (0x3ff6b1490bc5b4d1, 0xbfd65d5590807800),
41    (0x3ff69147332f0cba, 0xbfd602d076180000),
42    (0x3ff6719f18224223, 0xbfd5a8ca86909000),
43    (0x3ff6524f99a51ed9, 0xbfd54f4356035000),
44    (0x3ff63356aa8f24c4, 0xbfd4f637c36b4000),
45    (0x3ff614b36b9ddc14, 0xbfd49da7fda85000),
46    (0x3ff5f66452c65c4c, 0xbfd445923989a800),
47    (0x3ff5d867b5912c4f, 0xbfd3edf439b0b800),
48    (0x3ff5babccb5b90de, 0xbfd396ce448f7000),
49    (0x3ff59d61f2d91a78, 0xbfd3401e17bda000),
50    (0x3ff5805612465687, 0xbfd2e9e2ef468000),
51    (0x3ff56397cee76bd3, 0xbfd2941b3830e000),
52    (0x3ff54725e2a77f93, 0xbfd23ec58cda8800),
53    (0x3ff52aff42064583, 0xbfd1e9e129279000),
54    (0x3ff50f22dbb2bddf, 0xbfd1956d2b48f800),
55    (0x3ff4f38f4734ded7, 0xbfd141679ab9f800),
56    (0x3ff4d843cfde2840, 0xbfd0edd094ef9800),
57    (0x3ff4bd3ec078a3c8, 0xbfd09aa518db1000),
58    (0x3ff4a27fc3e0258a, 0xbfd047e65263b800),
59    (0x3ff4880524d48434, 0xbfcfeb224586f000),
60    (0x3ff46dce1b192d0b, 0xbfcf474a7517b000),
61    (0x3ff453d9d3391854, 0xbfcea4443d103000),
62    (0x3ff43a2744b4845a, 0xbfce020d44e9b000),
63    (0x3ff420b54115f8fb, 0xbfcd60a22977f000),
64    (0x3ff40782da3ef4b1, 0xbfccc00104959000),
65    (0x3ff3ee8f5d57fe8f, 0xbfcc202956891000),
66    (0x3ff3d5d9a00b4ce9, 0xbfcb81178d811000),
67    (0x3ff3bd60c010c12b, 0xbfcae2c9ccd3d000),
68    (0x3ff3a5242b75dab8, 0xbfca45402e129000),
69    (0x3ff38d22cd9fd002, 0xbfc9a877681df000),
70    (0x3ff3755bc5847a1c, 0xbfc90c6d69483000),
71    (0x3ff35dce49ad36e2, 0xbfc87120a645c000),
72    (0x3ff34679984dd440, 0xbfc7d68fb4143000),
73    (0x3ff32f5cceffcb24, 0xbfc73cb83c627000),
74    (0x3ff3187775a10d49, 0xbfc6a39a9b376000),
75    (0x3ff301c8373e3990, 0xbfc60b3154b7a000),
76    (0x3ff2eb4ebb95f841, 0xbfc5737d76243000),
77    (0x3ff2d50a0219a9d1, 0xbfc4dc7b8fc23000),
78    (0x3ff2bef9a8b7fd2a, 0xbfc4462c51d20000),
79    (0x3ff2a91c7a0c1bab, 0xbfc3b08abc830000),
80    (0x3ff293726014b530, 0xbfc31b996b490000),
81    (0x3ff27dfa5757a1f5, 0xbfc2875490a44000),
82    (0x3ff268b39b1d3bbf, 0xbfc1f3b9f879a000),
83    (0x3ff2539d838ff5bd, 0xbfc160c8252ca000),
84    (0x3ff23eb7aac9083b, 0xbfc0ce7f57f72000),
85    (0x3ff22a012ba940b6, 0xbfc03cdc49fea000),
86    (0x3ff2157996cc4132, 0xbfbf57bdbc4b8000),
87    (0x3ff201201dd2fc9b, 0xbfbe370896404000),
88    (0x3ff1ecf4494d480b, 0xbfbd17983ef94000),
89    (0x3ff1d8f5528f6569, 0xbfbbf9674ed8a000),
90    (0x3ff1c52311577e7c, 0xbfbadc79202f6000),
91    (0x3ff1b17c74cb26e9, 0xbfb9c0c3e7288000),
92    (0x3ff19e010c2c1ab6, 0xbfb8a646b372c000),
93    (0x3ff18ab07bb670bd, 0xbfb78d01b3ac0000),
94    (0x3ff1778a25efbcb6, 0xbfb674f145380000),
95    (0x3ff1648d354c31da, 0xbfb55e0e6d878000),
96    (0x3ff151b990275fdd, 0xbfb4485cdea1e000),
97    (0x3ff13f0ea432d24c, 0xbfb333d94d6aa000),
98    (0x3ff12c8b7210f9da, 0xbfb22079f8c56000),
99    (0x3ff11a3028ecb531, 0xbfb10e4698622000),
100    (0x3ff107fbda8434af, 0xbfaffa6c6ad20000),
101    (0x3ff0f5ee0f4e6bb3, 0xbfadda8d4a774000),
102    (0x3ff0e4065d2a9fce, 0xbfabbcece4850000),
103    (0x3ff0d244632ca521, 0xbfa9a1894012c000),
104    (0x3ff0c0a77ce2981a, 0xbfa788583302c000),
105    (0x3ff0af2f83c636d1, 0xbfa5715e67d68000),
106    (0x3ff09ddb98a01339, 0xbfa35c8a49658000),
107    (0x3ff08cabaf52e7df, 0xbfa149e364154000),
108    (0x3ff07b9f2f4e28fb, 0xbf9e72c082eb8000),
109    (0x3ff06ab58c358f19, 0xbf9a55f152528000),
110    (0x3ff059eea5ecf92c, 0xbf963d62cf818000),
111    (0x3ff04949cdd12c90, 0xbf9228fb8caa0000),
112    (0x3ff038c6c6f0ada9, 0xbf8c317b20f90000),
113    (0x3ff02865137932a9, 0xbf8419355daa0000),
114    (0x3ff0182427ea7348, 0xbf781203c2ec0000),
115    (0x3ff008040614b195, 0xbf60040979240000),
116    (0x3fefe01ff726fa1a, 0x3f6feff384900000),
117    (0x3fefa11cc261ea74, 0x3f87dc41353d0000),
118    (0x3fef6310b081992e, 0x3f93cea3c4c28000),
119    (0x3fef25f63ceeadcd, 0x3f9b9fc114890000),
120    (0x3feee9c8039113e7, 0x3fa1b0d8ce110000),
121    (0x3feeae8078cbb1ab, 0x3fa58a5bd001c000),
122    (0x3fee741aa29d0c9b, 0x3fa95c8340d88000),
123    (0x3fee3a91830a99b5, 0x3fad276aef578000),
124    (0x3fee01e009609a56, 0x3fb07598e598c000),
125    (0x3fedca01e577bb98, 0x3fb253f5e30d2000),
126    (0x3fed92f20b7c9103, 0x3fb42edd8b380000),
127    (0x3fed5cac66fb5cce, 0x3fb606598757c000),
128    (0x3fed272caa5ede9d, 0x3fb7da76356a0000),
129    (0x3fecf26e3e6b2ccd, 0x3fb9ab434e1c6000),
130    (0x3fecbe6da2a77902, 0x3fbb78c7bb0d6000),
131    (0x3fec8b266d37086d, 0x3fbd431332e72000),
132    (0x3fec5894bd5d5804, 0x3fbf0a3171de6000),
133    (0x3fec26b533bb9f8c, 0x3fc067152b914000),
134    (0x3febf583eeece73f, 0x3fc147858292b000),
135    (0x3febc4fd75db96c1, 0x3fc2266ecdca3000),
136    (0x3feb951e0c864a28, 0x3fc303d7a6c55000),
137    (0x3feb65e2c5ef3e2c, 0x3fc3dfc33c331000),
138    (0x3feb374867c9888b, 0x3fc4ba366b7a8000),
139    (0x3feb094b211d304a, 0x3fc5933928d1f000),
140    (0x3feadbe885f2ef7e, 0x3fc66acd2418f000),
141    (0x3feaaf1d31603da2, 0x3fc740f8ec669000),
142    (0x3fea82e63fd358a7, 0x3fc815c0f51af000),
143    (0x3fea5740ef09738b, 0x3fc8e92954f68000),
144    (0x3fea2c2a90ab4b27, 0x3fc9bb3602f84000),
145    (0x3fea01a01393f2d1, 0x3fca8bed1c2c0000),
146    (0x3fe9d79f24db3c1b, 0x3fcb5b515c01d000),
147    (0x3fe9ae2505c7b190, 0x3fcc2967ccbcc000),
148    (0x3fe9852ef297ce2f, 0x3fccf635d5486000),
149    (0x3fe95cbaeea44b75, 0x3fcdc1bd3446c000),
150    (0x3fe934c69de74838, 0x3fce8c01b8cfe000),
151    (0x3fe90d4f2f6752e6, 0x3fcf5509c0179000),
152    (0x3fe8e6528effd79d, 0x3fd00e6c121fb800),
153    (0x3fe8bfce9fcc007c, 0x3fd071b80e93d000),
154    (0x3fe899c0dabec30e, 0x3fd0d46b9e867000),
155    (0x3fe87427aa2317fb, 0x3fd13687334bd000),
156    (0x3fe84f00acb39a08, 0x3fd1980d67234800),
157    (0x3fe82a49e8653e55, 0x3fd1f8ffe0cc8000),
158    (0x3fe8060195f40260, 0x3fd2595fd7636800),
159    (0x3fe7e22563e0a329, 0x3fd2b9300914a800),
160    (0x3fe7beb377dcb5ad, 0x3fd3187210436000),
161    (0x3fe79baa679725c2, 0x3fd377266dec1800),
162    (0x3fe77907f2170657, 0x3fd3d54ffbaf3000),
163    (0x3fe756cadbd6130c, 0x3fd432eee32fe000),
164];
165
166const LOG_TAB2: [(u64, u64); 128] = [
167    (0x3fe61000014fb66b, 0x3c7e026c91425b3c),
168    (0x3fe63000034db495, 0x3c8dbfea48005d41),
169    (0x3fe650000d94d478, 0x3c8e7fa786d6a5b7),
170    (0x3fe67000074e6fad, 0x3c61fcea6b54254c),
171    (0x3fe68ffffedf0fae, 0xbc7c7e274c590efd),
172    (0x3fe6b0000763c5bc, 0xbc8ac16848dcda01),
173    (0x3fe6d0001e5cc1f6, 0x3c833f1c9d499311),
174    (0x3fe6efffeb05f63e, 0xbc7e80041ae22d53),
175    (0x3fe710000e869780, 0x3c7bff6671097952),
176    (0x3fe72ffffc67e912, 0x3c8c00e226bd8724),
177    (0x3fe74fffdf81116a, 0xbc6e02916ef101d2),
178    (0x3fe770000f679c90, 0xbc67fc71cd549c74),
179    (0x3fe78ffffa7ec835, 0x3c81bec19ef50483),
180    (0x3fe7affffe20c2e6, 0xbc707e1729cc6465),
181    (0x3fe7cfffed3fc900, 0xbc808072087b8b1c),
182    (0x3fe7efffe9261a76, 0x3c8dc0286d9df9ae),
183    (0x3fe81000049ca3e8, 0x3c897fd251e54c33),
184    (0x3fe8300017932c8f, 0xbc8afee9b630f381),
185    (0x3fe850000633739c, 0x3c89bfbf6b6535bc),
186    (0x3fe87000204289c6, 0xbc8bbf65f3117b75),
187    (0x3fe88fffebf57904, 0xbc89006ea23dcb57),
188    (0x3fe8b00022bc04df, 0xbc7d00df38e04b0a),
189    (0x3fe8cfffe50c1b8a, 0xbc88007146ff9f05),
190    (0x3fe8effffc918e43, 0x3c83817bd07a7038),
191    (0x3fe910001efa5fc7, 0x3c893e9176dfb403),
192    (0x3fe9300013467bb9, 0x3c7f804e4b980276),
193    (0x3fe94fffe6ee076f, 0xbc8f7ef0d9ff622e),
194    (0x3fe96fffde3c12d1, 0xbc7082aa962638ba),
195    (0x3fe98ffff4458a0d, 0xbc87801b9164a8ef),
196    (0x3fe9afffdd982e3e, 0xbc8740e08a5a9337),
197    (0x3fe9cfffed49fb66, 0x3c3fce08c19be000),
198    (0x3fe9f00020f19c51, 0xbc8a3faa27885b0a),
199    (0x3fea10001145b006, 0x3c74ff489958da56),
200    (0x3fea300007bbf6fa, 0x3c8cbeab8a2b6d18),
201    (0x3fea500010971d79, 0x3c88fecadd787930),
202    (0x3fea70001df52e48, 0xbc8f41763dd8abdb),
203    (0x3fea90001c593352, 0xbc8ebf0284c27612),
204    (0x3feab0002a4f3e4b, 0xbc69fd043cff3f5f),
205    (0x3feacfffd7ae1ed1, 0xbc823ee7129070b4),
206    (0x3feaefffee510478, 0x3c6a063ee00edea3),
207    (0x3feb0fffdb650d5b, 0x3c5a06c8381f0ab9),
208    (0x3feb2ffffeaaca57, 0xbc79011e74233c1d),
209    (0x3feb4fffd995badc, 0xbc79ff1068862a9f),
210    (0x3feb7000249e659c, 0x3c8aff45d0864f3e),
211    (0x3feb8ffff9871640, 0x3c7cfe7796c2c3f9),
212    (0x3febafffd204cb4f, 0xbc63ff27eef22bc4),
213    (0x3febcfffd2415c45, 0xbc6cffb7ee3bea21),
214    (0x3febeffff86309df, 0xbc814103972e0b5c),
215    (0x3fec0fffe1b57653, 0x3c8bc16494b76a19),
216    (0x3fec2ffff1fa57e3, 0xbc64feef8d30c6ed),
217    (0x3fec4fffdcbfe424, 0xbc843f68bcec4775),
218    (0x3fec6fffed54b9f7, 0x3c847ea3f053e0ec),
219    (0x3fec8fffeb998fd5, 0x3c7383068df992f1),
220    (0x3fecb0002125219a, 0xbc68fd8e64180e04),
221    (0x3feccfffdd94469c, 0x3c8e7ebe1cc7ea72),
222    (0x3fecefffeafdc476, 0x3c8ebe39ad9f88fe),
223    (0x3fed1000169af82b, 0x3c757d91a8b95a71),
224    (0x3fed30000d0ff71d, 0x3c89c1906970c7da),
225    (0x3fed4fffea790fc4, 0xbc580e37c558fe0c),
226    (0x3fed70002edc87e5, 0xbc7f80d64dc10f44),
227    (0x3fed900021dc82aa, 0xbc747c8f94fd5c5c),
228    (0x3fedafffd86b0283, 0x3c8c7f1dc521617e),
229    (0x3fedd000296c4739, 0x3c88019eb2ffb153),
230    (0x3fedefffe54490f5, 0x3c6e00d2c652cc89),
231    (0x3fee0fffcdabf694, 0xbc7f8340202d69d2),
232    (0x3fee2fffdb52c8dd, 0x3c7b00c1ca1b0864),
233    (0x3fee4ffff24216ef, 0x3c72ffa8b094ab51),
234    (0x3fee6fffe88a5e11, 0xbc57f673b1efbe59),
235    (0x3fee9000119eff0d, 0xbc84808d5e0bc801),
236    (0x3feeafffdfa51744, 0x3c780006d54320b5),
237    (0x3feed0001a127fa1, 0xbc5002f860565c92),
238    (0x3feef00007babcc4, 0xbc8540445d35e611),
239    (0x3fef0ffff57a8d02, 0xbc4ffb3139ef9105),
240    (0x3fef30001ee58ac7, 0x3c8a81acf2731155),
241    (0x3fef4ffff5823494, 0x3c8a3f41d4d7c743),
242    (0x3fef6ffffca94c6b, 0xbc6202f41c987875),
243    (0x3fef8fffe1f9c441, 0x3c777dd1f477e74b),
244    (0x3fefafffd2e0e37e, 0xbc6f01199a7ca331),
245    (0x3fefd0001c77e49e, 0x3c7181ee4bceacb1),
246    (0x3fefeffff7e0c331, 0xbc6e05370170875a),
247    (0x3ff00ffff465606e, 0xbc8a7ead491c0ada),
248    (0x3ff02ffff3867a58, 0xbc977f69c3fcb2e0),
249    (0x3ff04ffffdfc0d17, 0x3c97bffe34cb945b),
250    (0x3ff0700003cd4d82, 0x3c820083c0e456cb),
251    (0x3ff08ffff9f2cbe8, 0xbc6dffdfbe37751a),
252    (0x3ff0b000010cda65, 0xbc913f7faee626eb),
253    (0x3ff0d00001a4d338, 0x3c807dfa79489ff7),
254    (0x3ff0effffadafdfd, 0xbc77040570d66bc0),
255    (0x3ff110000bbafd96, 0x3c8e80d4846d0b62),
256    (0x3ff12ffffae5f45d, 0x3c9dbffa64fd36ef),
257    (0x3ff150000dd59ad9, 0x3c9a0077701250ae),
258    (0x3ff170000f21559a, 0x3c8dfdf9e2e3deee),
259    (0x3ff18ffffc275426, 0x3c910030dc3b7273),
260    (0x3ff1b000123d3c59, 0x3c997f7980030188),
261    (0x3ff1cffff8299eb7, 0xbc65f932ab9f8c67),
262    (0x3ff1effff48ad400, 0x3c937fbf9da75beb),
263    (0x3ff210000c8b86a4, 0x3c9f806b91fd5b22),
264    (0x3ff2300003854303, 0x3c93ffc2eb9fbf33),
265    (0x3ff24fffffbcf684, 0x3c7601e77e2e2e72),
266    (0x3ff26ffff52921d9, 0x3c7ffcbb767f0c61),
267    (0x3ff2900014933a3c, 0xbc7202ca3c02412b),
268    (0x3ff2b00014556313, 0xbc92808233f21f02),
269    (0x3ff2cfffebfe523b, 0xbc88ff7e384fdcf2),
270    (0x3ff2f0000bb8ad96, 0xbc85ff51503041c5),
271    (0x3ff30ffffb7ae2af, 0xbc810071885e289d),
272    (0x3ff32ffffeac5f7f, 0xbc91ff5d3fb7b715),
273    (0x3ff350000ca66756, 0x3c957f82228b82bd),
274    (0x3ff3700011fbf721, 0x3c8000bac40dd5cc),
275    (0x3ff38ffff9592fb9, 0xbc943f9d2db2a751),
276    (0x3ff3b00004ddd242, 0x3c857f6b707638e1),
277    (0x3ff3cffff5b2c957, 0x3c7a023a10bf1231),
278    (0x3ff3efffeab0b418, 0x3c987f6d66b152b0),
279    (0x3ff410001532aff4, 0x3c67f8375f198524),
280    (0x3ff4300017478b29, 0x3c8301e672dc5143),
281    (0x3ff44fffe795b463, 0x3c89ff69b8b2895a),
282    (0x3ff46fffe80475e0, 0xbc95c0b19bc2f254),
283    (0x3ff48fffef6fc1e7, 0x3c9b4009f23a2a72),
284    (0x3ff4afffe5bea704, 0xbc94ffb7bf0d7d45),
285    (0x3ff4d000171027de, 0xbc99c06471dc6a3d),
286    (0x3ff4f0000ff03ee2, 0x3c977f890b85531c),
287    (0x3ff5100012dc4bd1, 0x3c6004657166a436),
288    (0x3ff530001605277a, 0xbc96bfcece233209),
289    (0x3ff54fffecdb704c, 0xbc8902720505a1d7),
290    (0x3ff56fffef5f54a9, 0x3c9bbfe60ec96412),
291    (0x3ff5900017e61012, 0x3c887ec581afef90),
292    (0x3ff5b00003c93e92, 0xbc9f41080abf0cc0),
293    (0x3ff5d0001d4919bc, 0xbc98812afb254729),
294    (0x3ff5efffe7b87a89, 0xbc947eb780ed6904),
295];
296
297// ============================================================
298// logf (single precision) data
299// ============================================================
300
301const LOGF_TABLE_BITS: u32 = 4;
302const LOGF_N: usize = 1 << LOGF_TABLE_BITS;
303const LOGF_LN2: f64 = asdouble(0x3fe62e42fefa39ef);
304
305const LOGF_TAB: [(f64, f64); 16] = [
306    (asdouble(0x3ff661ec79f8f3be), asdouble(0xbfd57bf7808caade)),
307    (asdouble(0x3ff571ed4aaf883d), asdouble(0xbfd2bef0a7c06ddb)),
308    (asdouble(0x3ff49539f0f010b0), asdouble(0xbfd01eae7f513a67)),
309    (asdouble(0x3ff3c995b0b80385), asdouble(0xbfcb31d8a68224e9)),
310    (asdouble(0x3ff30d190c8864a5), asdouble(0xbfc6574f0ac07758)),
311    (asdouble(0x3ff25e227b0b8ea0), asdouble(0xbfc1aa2bc79c8100)),
312    (asdouble(0x3ff1bb4a4a1a343f), asdouble(0xbfba4e76ce8c0e5e)),
313    (asdouble(0x3ff12358f08ae5ba), asdouble(0xbfb1973c5a611ccc)),
314    (asdouble(0x3ff0953f419900a7), asdouble(0xbfa252f438e10c1e)),
315    (asdouble(0x3ff0000000000000), asdouble(0x0000000000000000)),
316    (asdouble(0x3fee608cfd9a47ac), asdouble(0x3faaa5aa5df25984)),
317    (asdouble(0x3feca4b31f026aa0), asdouble(0x3fbc5e53aa362eb4)),
318    (asdouble(0x3feb2036576afce6), asdouble(0x3fc526e57720db08)),
319    (asdouble(0x3fe9c2d163a1aa2d), asdouble(0x3fcbc2860d224770)),
320    (asdouble(0x3fe886e6037841ed), asdouble(0x3fd1058bc8a07ee1)),
321    (asdouble(0x3fe767dcf5534862), asdouble(0x3fd4043057b6ee09)),
322];
323
324const LOGF_POLY: [f64; 3] = [
325    asdouble(0xbfd00ea348b88334),
326    asdouble(0x3fd5575b0be00b6a),
327    asdouble(0xbfdffffef20a4123),
328];
329
330// ============================================================
331// log2 (double precision) data
332// ============================================================
333
334const LOG2_TABLE_BITS: u32 = 6;
335const LOG2_N: usize = 1 << LOG2_TABLE_BITS;
336const INVLN2HI: f64 = asdouble(0x3ff7154765200000);
337const INVLN2LO: f64 = asdouble(0x3de705fc2eefa200);
338
339const LOG2_POLY1: [f64; 10] = [
340    asdouble(0xbfe71547652b82fe),
341    asdouble(0x3fdec709dc3a03f7),
342    asdouble(0xbfd71547652b7c3f),
343    asdouble(0x3fd2776c50f05be4),
344    asdouble(0xbfcec709dd768fe5),
345    asdouble(0x3fca61761ec4e736),
346    asdouble(0xbfc7153fbc64a79b),
347    asdouble(0x3fc484d154f01b4a),
348    asdouble(0xbfc289e4a72c383c),
349    asdouble(0x3fc0b32f285aee66),
350];
351
352const LOG2_POLY: [f64; 6] = [
353    asdouble(0xbfe71547652b8339),
354    asdouble(0x3fdec709dc3a04be),
355    asdouble(0xbfd7154764702ffb),
356    asdouble(0x3fd2776c50034c48),
357    asdouble(0xbfcec7b328ea92bc),
358    asdouble(0x3fca6225e117f92e),
359];
360
361const LOG2_TAB: [(u64, u64); 64] = [
362    (0x3ff724286bb1acf8, 0xbfe1095feecdb000),
363    (0x3ff6e1f766d2cca1, 0xbfe08494bd76d000),
364    (0x3ff6a13d0e30d48a, 0xbfe00143aee8f800),
365    (0x3ff661ec32d06c85, 0xbfdefec5360b4000),
366    (0x3ff623fa951198f8, 0xbfddfdd91ab7e000),
367    (0x3ff5e75ba4cf026c, 0xbfdcffae0cc79000),
368    (0x3ff5ac055a214fb8, 0xbfdc043811fda000),
369    (0x3ff571ed0f166e1e, 0xbfdb0b67323ae000),
370    (0x3ff53909590bf835, 0xbfda152f5a2db000),
371    (0x3ff5014fed61addd, 0xbfd9217f5af86000),
372    (0x3ff4cab88e487bd0, 0xbfd8304db0719000),
373    (0x3ff49539b4334fee, 0xbfd74189f9a9e000),
374    (0x3ff460cbdfafd569, 0xbfd6552bb5199000),
375    (0x3ff42d664ee4b953, 0xbfd56b23a29b1000),
376    (0x3ff3fb01111dd8a6, 0xbfd483650f5fa000),
377    (0x3ff3c995b70c5836, 0xbfd39de937f6a000),
378    (0x3ff3991c4ab6fd4a, 0xbfd2baa1538d6000),
379    (0x3ff3698e0ce099b5, 0xbfd1d98340ca4000),
380    (0x3ff33ae48213e7b2, 0xbfd0fa853a40e000),
381    (0x3ff30d191985bdb1, 0xbfd01d9c32e73000),
382    (0x3ff2e025cab271d7, 0xbfce857da2fa6000),
383    (0x3ff2b404cf13cd82, 0xbfccd3c8633d8000),
384    (0x3ff288b02c7ccb50, 0xbfcb26034c14a000),
385    (0x3ff25e2263944de5, 0xbfc97c1c2f4fe000),
386    (0x3ff234563d8615b1, 0xbfc7d6023f800000),
387    (0x3ff20b46e33eaf38, 0xbfc633a71a05e000),
388    (0x3ff1e2eefdcda3dd, 0xbfc494f5e9570000),
389    (0x3ff1bb4a580b3930, 0xbfc2f9e424e0a000),
390    (0x3ff19453847f2200, 0xbfc162595afdc000),
391    (0x3ff16e06c0d5d73c, 0xbfbf9c9a75bd8000),
392    (0x3ff1485f47b7e4c2, 0xbfbc7b575bf9c000),
393    (0x3ff12358ad0085d1, 0xbfb960c60ff48000),
394    (0x3ff0fef00f532227, 0xbfb64ce247b60000),
395    (0x3ff0db2077d03a8f, 0xbfb33f78b2014000),
396    (0x3ff0b7e6d65980d9, 0xbfb0387d1a42c000),
397    (0x3ff0953efe7b408d, 0xbfaa6f9208b50000),
398    (0x3ff07325cac53b83, 0xbfa47a954f770000),
399    (0x3ff05197e40d1b5c, 0xbf9d23a8c50c0000),
400    (0x3ff03091c1208ea2, 0xbf916a2629780000),
401    (0x3ff0101025b37e21, 0xbf7720f8d8e80000),
402    (0x3fefc07ef9caa76b, 0x3f86fe53b1500000),
403    (0x3fef4465d3f6f184, 0x3fa11ccce10f8000),
404    (0x3feecc079f84107f, 0x3fac4dfc8c8b8000),
405    (0x3fee573a99975ae8, 0x3fb3aa321e574000),
406    (0x3fede5d6f0bd3de6, 0x3fb918a0d08b8000),
407    (0x3fed77b681ff38b3, 0x3fbe72e9da044000),
408    (0x3fed0cb5724de943, 0x3fc1dcd2507f6000),
409    (0x3feca4b2dc0e7563, 0x3fc476ab03dea000),
410    (0x3fec3f8ee8d6cb51, 0x3fc7074377e22000),
411    (0x3febdd2b4f020c4c, 0x3fc98ede8ba94000),
412    (0x3feb7d6c006015ca, 0x3fcc0db86ad2e000),
413    (0x3feb20366e2e338f, 0x3fce840aafcee000),
414    (0x3feac57026295039, 0x3fd0790ab4678000),
415    (0x3fea6d01bc2731dd, 0x3fd1ac056801c000),
416    (0x3fea16d3bc3ff18b, 0x3fd2db11d4fee000),
417    (0x3fe9c2d14967fead, 0x3fd406464ec58000),
418    (0x3fe970e4f47c9902, 0x3fd52dbe093af000),
419    (0x3fe920fb3982bcf2, 0x3fd651902050d000),
420    (0x3fe8d30187f759f1, 0x3fd771d2cdeaf000),
421    (0x3fe886e5ebb9f66d, 0x3fd88e9c857d9000),
422    (0x3fe83c97b658b994, 0x3fd9a80155e16000),
423    (0x3fe7f405ffc61022, 0x3fdabe186ed3d000),
424    (0x3fe7ad22181415ca, 0x3fdbd0f2aea0e000),
425    (0x3fe767dcf99eff8c, 0x3fdce0a43dbf4000),
426];
427
428const LOG2_TAB2: [(u64, u64); 64] = [
429    (0x3fe6200012b90a8e, 0x3c8904ab0644b605),
430    (0x3fe66000045734a6, 0x3c61ff9bea62f7a9),
431    (0x3fe69fffc325f2c5, 0x3c827ecfcb3c90ba),
432    (0x3fe6e00038b95a04, 0x3c88ff8856739326),
433    (0x3fe71fffe09994e3, 0x3c8afd40275f82b1),
434    (0x3fe7600015590e10, 0xbc72fd75b4238341),
435    (0x3fe7a00012655bd5, 0x3c7808e67c242b76),
436    (0x3fe7e0003259e9a6, 0xbc6208e426f622b7),
437    (0x3fe81fffedb4b2d2, 0xbc8402461ea5c92f),
438    (0x3fe860002dfafcc3, 0x3c6df7f4a2f29a1f),
439    (0x3fe89ffff78c6b50, 0xbc8e0453094995fd),
440    (0x3fe8e00039671566, 0xbc8a04f3bec77b45),
441    (0x3fe91fffe2bf1745, 0xbc77fa34400e203c),
442    (0x3fe95fffcc5c9fd1, 0xbc76ff8005a0695d),
443    (0x3fe9a0003bba4767, 0x3c70f8c4c4ec7e03),
444    (0x3fe9dfffe7b92da5, 0x3c8e7fd9478c4602),
445    (0x3fea1fffd72efdaf, 0xbc6a0c554dcdae7e),
446    (0x3fea5fffde04ff95, 0x3c867da98ce9b26b),
447    (0x3fea9fffca5e8d2b, 0xbc8284c9b54c13de),
448    (0x3feadfffddad03ea, 0x3c5812c8ea602e3c),
449    (0x3feb1ffff10d3d4d, 0xbc8efaddad27789c),
450    (0x3feb5fffce21165a, 0x3c53cb1719c61237),
451    (0x3feb9fffd950e674, 0x3c73f7d94194ce00),
452    (0x3febe000139ca8af, 0x3c750ac4215d9bc0),
453    (0x3fec20005b46df99, 0x3c6beea653e9c1c9),
454    (0x3fec600040b9f7ae, 0xbc7c079f274a70d6),
455    (0x3feca0006255fd8a, 0xbc7a0b4076e84c1f),
456    (0x3fecdfffd94c095d, 0x3c88f933f99ab5d7),
457    (0x3fed1ffff975d6cf, 0xbc582c08665fe1be),
458    (0x3fed5fffa2561c93, 0xbc7b04289bd295f3),
459    (0x3fed9fff9d228b0c, 0x3c870251340fa236),
460    (0x3fede00065bc7e16, 0xbc75011e16a4d80c),
461    (0x3fee200002f64791, 0x3c89802f09ef62e0),
462    (0x3fee600057d7a6d8, 0xbc7e0b75580cf7fa),
463    (0x3feea00027edc00c, 0xbc8c848309459811),
464    (0x3feee0006cf5cb7c, 0xbc8f8027951576f4),
465    (0x3fef2000782b7dcc, 0xbc8f81d97274538f),
466    (0x3fef6000260c450a, 0xbc4071002727ffdc),
467    (0x3fef9fffe88cd533, 0xbc581bdce1fda8b0),
468    (0x3fefdfffd50f8689, 0x3c87f91acb918e6e),
469    (0x3ff0200004292367, 0x3c9b7ff365324681),
470    (0x3ff05fffe3e3d668, 0x3c86fa08ddae957b),
471    (0x3ff0a0000a85a757, 0xbc57e2de80d3fb91),
472    (0x3ff0e0001a5f3fcc, 0xbc91823305c5f014),
473    (0x3ff11ffff8afbaf5, 0xbc8bfabb6680bac2),
474    (0x3ff15fffe54d91ad, 0xbc9d7f121737e7ef),
475    (0x3ff1a00011ac36e1, 0x3c9c000a0516f5ff),
476    (0x3ff1e00019c84248, 0xbc9082fbe4da5da0),
477    (0x3ff220000ffe5e6e, 0xbc88fdd04c9cfb43),
478    (0x3ff26000269fd891, 0x3c8cfe2a7994d182),
479    (0x3ff2a00029a6e6da, 0xbc700273715e8bc5),
480    (0x3ff2dfffe0293e39, 0x3c9b7c39dab2a6f9),
481    (0x3ff31ffff7dcf082, 0x3c7df1336edc5254),
482    (0x3ff35ffff05a8b60, 0xbc9e03564ccd31eb),
483    (0x3ff3a0002e0eaecc, 0x3c75f0e74bd3a477),
484    (0x3ff3e000043bb236, 0x3c9c7dcb149d8833),
485    (0x3ff4200002d187ff, 0x3c7e08afcf2d3d28),
486    (0x3ff460000d387cb1, 0x3c820837856599a6),
487    (0x3ff4a00004569f89, 0xbc89fa5c904fbcd2),
488    (0x3ff4e000043543f3, 0xbc781125ed175329),
489    (0x3ff51fffcc027f0f, 0x3c9883d8847754dc),
490    (0x3ff55ffffd87b36f, 0xbc8709e731d02807),
491    (0x3ff59ffff21df7ba, 0x3c87f79f68727b02),
492    (0x3ff5dfffebfc3481, 0xbc9180902e30e93e),
493];
494
495// ============================================================
496// log2f (single precision) data
497// ============================================================
498
499const LOG2F_TABLE_BITS: u32 = 4;
500const LOG2F_N: usize = 1 << LOG2F_TABLE_BITS;
501
502const LOG2F_TAB: [(f64, f64); 16] = [
503    (asdouble(0x3ff661ec79f8f3be), asdouble(0xbfdefec65b963019)),
504    (asdouble(0x3ff571ed4aaf883d), asdouble(0xbfdb0b6832d4fca4)),
505    (asdouble(0x3ff49539f0f010b0), asdouble(0xbfd7418b0a1fb77b)),
506    (asdouble(0x3ff3c995b0b80385), asdouble(0xbfd39de91a6dcf7b)),
507    (asdouble(0x3ff30d190c8864a5), asdouble(0xbfd01d9bf3f2b631)),
508    (asdouble(0x3ff25e227b0b8ea0), asdouble(0xbfc97c1d1b3b7af0)),
509    (asdouble(0x3ff1bb4a4a1a343f), asdouble(0xbfc2f9e393af3c9f)),
510    (asdouble(0x3ff12358f08ae5ba), asdouble(0xbfb960cbbf788d5c)),
511    (asdouble(0x3ff0953f419900a7), asdouble(0xbfaa6f9db6475fce)),
512    (asdouble(0x3ff0000000000000), asdouble(0x0000000000000000)),
513    (asdouble(0x3fee608cfd9a47ac), asdouble(0x3fb338ca9f24f53d)),
514    (asdouble(0x3feca4b31f026aa0), asdouble(0x3fc476a9543891ba)),
515    (asdouble(0x3feb2036576afce6), asdouble(0x3fce840b4ac4e4d2)),
516    (asdouble(0x3fe9c2d163a1aa2d), asdouble(0x3fd40645f0c6651c)),
517    (asdouble(0x3fe886e6037841ed), asdouble(0x3fd88e9c2c1b9ff8)),
518    (asdouble(0x3fe767dcf5534862), asdouble(0x3fdce0a44eb17bcc)),
519];
520
521const LOG2F_POLY: [f64; 4] = [
522    asdouble(0xbfd712b6f70a7e4d),
523    asdouble(0x3fdecabf496832e0),
524    asdouble(0xbfe715479ffae3de),
525    asdouble(0x3ff715475f35c8b8),
526];
527
528// ============================================================
529// log10 / log10f scalar constants
530// ============================================================
531
532const IVLN10HI: f64 = asdouble(0x3fdbcb7b15200000);
533const IVLN10LO: f64 = asdouble(0x3dbb9438ca9aadd5);
534const LOG10_2HI: f64 = asdouble(0x3fd34413509f6000);
535const LOG10_2LO: f64 = asdouble(0x3d59fef311f12b36);
536const LG1: f64 = asdouble(0x3fe5555555555593);
537const LG2: f64 = asdouble(0x3fd999999997fa04);
538const LG3: f64 = asdouble(0x3fd2492494229359);
539const LG4: f64 = asdouble(0x3fcc71c51d8e78af);
540const LG5: f64 = asdouble(0x3fc7466496cb03de);
541const LG6: f64 = asdouble(0x3fc39a09d078c69f);
542const LG7: f64 = asdouble(0x3fc2f112df3e5244);
543
544const IVLN10HI_F: f32 = asfloat(0x3ede6000);
545const IVLN10LO_F: f32 = asfloat(0xb804ead9);
546const LOG10_2HI_F: f32 = asfloat(0x3e9a2080);
547const LOG10_2LO_F: f32 = asfloat(0x355427db);
548const LG1_F: f32 = asfloat(0x3f2aaaaa);
549const LG2_F: f32 = asfloat(0x3eccce13);
550const LG3_F: f32 = asfloat(0x3e91e9ee);
551const LG4_F: f32 = asfloat(0x3e789e26);
552
553// ============================================================
554// Helper
555// ============================================================
556
557#[inline]
558fn log_top16(x: f64) -> u32 {
559    (asuint64(x) >> 48) as u32
560}
561
562// ============================================================
563// log (double precision)
564// ============================================================
565
566const LOG_LO: u64 = asuint64(asdouble(0x3fee000000000000));
567const LOG_HI: u64 = asuint64(asdouble(0x3ff1090000000000));
568
569#[no_mangle]
570pub extern "C" fn log(x: f64) -> f64 {
571    let ix = asuint64(x);
572    let top = log_top16(x);
573
574    if predict_false(ix.wrapping_sub(LOG_LO) < LOG_HI.wrapping_sub(LOG_LO)) {
575        if predict_false(ix == asuint64(1.0)) {
576            return 0.0;
577        }
578        let r = x - 1.0;
579        let r2 = r * r;
580        let r3 = r * r2;
581        let y = r3
582            * (LOG_POLY1[1]
583                + r * LOG_POLY1[2]
584                + r2 * LOG_POLY1[3]
585                + r3
586                    * (LOG_POLY1[4]
587                        + r * LOG_POLY1[5]
588                        + r2 * LOG_POLY1[6]
589                        + r3
590                            * (LOG_POLY1[7]
591                                + r * LOG_POLY1[8]
592                                + r2 * LOG_POLY1[9]
593                                + r3 * LOG_POLY1[10])));
594        let w = r * 134217728.0;
595        let rhi = r + w - w;
596        let rlo = r - rhi;
597        let w = rhi * rhi * LOG_POLY1[0];
598        let hi = r + w;
599        let lo = r - hi + w;
600        let lo = lo + LOG_POLY1[0] * rlo * (rhi + r);
601        let y = y + lo;
602        let y = y + hi;
603        return eval_as_double(y);
604    }
605
606    if predict_false(top.wrapping_sub(0x0010) >= 0x7ff0 - 0x0010) {
607        if ix * 2 == 0 {
608            return __math_divzero(1);
609        }
610        if ix == asuint64(f64::INFINITY) {
611            return x;
612        }
613        if (top & 0x8000) != 0 || (top & 0x7ff0) == 0x7ff0 {
614            return __math_invalid(x);
615        }
616        let ix = asuint64(x * 4503599627370496.0);
617        let ix = ix - (52u64 << 52);
618        return log_main(ix);
619    }
620
621    log_main(ix)
622}
623
624#[inline]
625fn log_main(ix: u64) -> f64 {
626    const OFF: u64 = 0x3fe6000000000000;
627
628    let tmp = ix - OFF;
629    let i = ((tmp >> (52 - LOG_TABLE_BITS)) % LOG_N as u64) as usize;
630    let k = (tmp as i64 >> 52) as i32;
631    let iz = ix - (tmp & (0xfffu64 << 52));
632    let invc = asdouble(LOG_TAB[i].0);
633    let logc = asdouble(LOG_TAB[i].1);
634    let z = asdouble(iz);
635
636    let r = (z - asdouble(LOG_TAB2[i].0) - asdouble(LOG_TAB2[i].1)) * invc;
637    let kd = k as f64;
638
639    let w = kd * LN2HI + logc;
640    let hi = w + r;
641    let lo = w - hi + r + kd * LN2LO;
642
643    let r2 = r * r;
644    let y = lo
645        + r2 * LOG_POLY[0]
646        + r * r2 * (LOG_POLY[1] + r * LOG_POLY[2] + r2 * (LOG_POLY[3] + r * LOG_POLY[4]))
647        + hi;
648    eval_as_double(y)
649}
650
651// ============================================================
652// logf (single precision)
653// ============================================================
654
655#[no_mangle]
656pub extern "C" fn logf(x: f32) -> f32 {
657    let ix = asuint(x);
658
659    if predict_false(ix == 0x3f800000) {
660        return 0.0;
661    }
662    if predict_false(ix.wrapping_sub(0x00800000) >= 0x7f800000 - 0x00800000) {
663        if ix * 2 == 0 {
664            return __math_divzerof(1);
665        }
666        if ix == 0x7f800000 {
667            return x;
668        }
669        if (ix & 0x80000000) != 0 || ix * 2 >= 0xff000000 {
670            return __math_invalidf(x);
671        }
672        let ix = asuint(x * 8388608.0f32);
673        let ix = ix - (23u32 << 23);
674        return logf_main(ix);
675    }
676
677    logf_main(ix)
678}
679
680#[inline]
681fn logf_main(ix: u32) -> f32 {
682    const OFF: u32 = 0x3f330000;
683
684    let tmp = ix - OFF;
685    let i = ((tmp >> (23 - LOGF_TABLE_BITS)) % LOGF_N as u32) as usize;
686    let k = (tmp as i32 >> 23) as i32;
687    let iz = ix - (tmp & 0xff800000);
688    let invc = LOGF_TAB[i].0;
689    let logc = LOGF_TAB[i].1;
690    let z = asfloat(iz) as f64;
691
692    let r = z * invc - 1.0;
693    let y0 = logc + k as f64 * LOGF_LN2;
694
695    let r2 = r * r;
696    let y = LOGF_POLY[1] * r + LOGF_POLY[2];
697    let y = LOGF_POLY[0] * r2 + y;
698    let y = y * r2 + (y0 + r);
699    eval_as_float(y as f32)
700}
701
702// ============================================================
703// log2 (double precision)
704// ============================================================
705
706const LOG2_LO: u64 = asuint64(asdouble(0x3feea4af00000000));
707const LOG2_HI: u64 = asuint64(asdouble(0x3ff0b55900000000));
708
709#[no_mangle]
710pub extern "C" fn log2(x: f64) -> f64 {
711    let ix = asuint64(x);
712    let top = log_top16(x);
713
714    if predict_false(ix.wrapping_sub(LOG2_LO) < LOG2_HI.wrapping_sub(LOG2_LO)) {
715        if predict_false(ix == asuint64(1.0)) {
716            return 0.0;
717        }
718        let r = x - 1.0;
719        let rhi = asdouble(asuint64(r) & (0xFFFFFFFFu64 << 32));
720        let rlo = r - rhi;
721        let hi = rhi * INVLN2HI;
722        let lo = rlo * INVLN2HI + r * INVLN2LO;
723        let r2 = r * r;
724        let r4 = r2 * r2;
725        let p = r2 * (LOG2_POLY1[0] + r * LOG2_POLY1[1]);
726        let y = hi + p;
727        let lo = lo + hi - y + p;
728        let lo = lo
729            + r4
730                * (LOG2_POLY1[2]
731                    + r * LOG2_POLY1[3]
732                    + r2 * (LOG2_POLY1[4] + r * LOG2_POLY1[5])
733                    + r4
734                        * (LOG2_POLY1[6]
735                            + r * LOG2_POLY1[7]
736                            + r2 * (LOG2_POLY1[8] + r * LOG2_POLY1[9])));
737        let y = y + lo;
738        return eval_as_double(y);
739    }
740
741    if predict_false(top.wrapping_sub(0x0010) >= 0x7ff0 - 0x0010) {
742        if ix * 2 == 0 {
743            return __math_divzero(1);
744        }
745        if ix == asuint64(f64::INFINITY) {
746            return x;
747        }
748        if (top & 0x8000) != 0 || (top & 0x7ff0) == 0x7ff0 {
749            return __math_invalid(x);
750        }
751        let ix = asuint64(x * 4503599627370496.0);
752        let ix = ix - (52u64 << 52);
753        return log2_main(ix);
754    }
755
756    log2_main(ix)
757}
758
759#[inline]
760fn log2_main(ix: u64) -> f64 {
761    const OFF: u64 = 0x3fe6000000000000;
762
763    let tmp = ix - OFF;
764    let i = ((tmp >> (52 - LOG2_TABLE_BITS)) % LOG2_N as u64) as usize;
765    let k = (tmp as i64 >> 52) as i32;
766    let iz = ix - (tmp & (0xfffu64 << 52));
767    let invc = asdouble(LOG2_TAB[i].0);
768    let logc = asdouble(LOG2_TAB[i].1);
769    let z = asdouble(iz);
770    let kd = k as f64;
771
772    let r = (z - asdouble(LOG2_TAB2[i].0) - asdouble(LOG2_TAB2[i].1)) * invc;
773    let rhi = asdouble(asuint64(r) & (0xFFFFFFFFu64 << 32));
774    let rlo = r - rhi;
775    let t1 = rhi * INVLN2HI;
776    let t2 = rlo * INVLN2HI + r * INVLN2LO;
777
778    let t3 = kd + logc;
779    let hi = t3 + t1;
780    let lo = t3 - hi + t1 + t2;
781
782    let r2 = r * r;
783    let r4 = r2 * r2;
784    let p = LOG2_POLY[0]
785        + r * LOG2_POLY[1]
786        + r2 * (LOG2_POLY[2] + r * LOG2_POLY[3])
787        + r4 * (LOG2_POLY[4] + r * LOG2_POLY[5]);
788    let y = lo + r2 * p + hi;
789    eval_as_double(y)
790}
791
792// ============================================================
793// log2f (single precision)
794// ============================================================
795
796#[no_mangle]
797pub extern "C" fn log2f(x: f32) -> f32 {
798    let ix = asuint(x);
799
800    if predict_false(ix == 0x3f800000) {
801        return 0.0;
802    }
803    if predict_false(ix.wrapping_sub(0x00800000) >= 0x7f800000 - 0x00800000) {
804        if ix * 2 == 0 {
805            return __math_divzerof(1);
806        }
807        if ix == 0x7f800000 {
808            return x;
809        }
810        if (ix & 0x80000000) != 0 || ix * 2 >= 0xff000000 {
811            return __math_invalidf(x);
812        }
813        let ix = asuint(x * 8388608.0f32);
814        let ix = ix - (23u32 << 23);
815        return log2f_main(ix);
816    }
817
818    log2f_main(ix)
819}
820
821#[inline]
822fn log2f_main(ix: u32) -> f32 {
823    const OFF: u32 = 0x3f330000;
824
825    let tmp = ix - OFF;
826    let i = ((tmp >> (23 - LOG2F_TABLE_BITS)) % LOG2F_N as u32) as usize;
827    let top = tmp & 0xff800000;
828    let iz = ix - top;
829    let k = (tmp as i32 >> 23) as i32;
830    let invc = LOG2F_TAB[i].0;
831    let logc = LOG2F_TAB[i].1;
832    let z = asfloat(iz) as f64;
833
834    let r = z * invc - 1.0;
835    let y0 = logc + k as f64;
836
837    let r2 = r * r;
838    let y = LOG2F_POLY[1] * r + LOG2F_POLY[2];
839    let y = LOG2F_POLY[0] * r2 + y;
840    let p = LOG2F_POLY[3] * r + y0;
841    let y = y * r2 + p;
842    eval_as_float(y as f32)
843}
844
845// ============================================================
846// log10 (double precision)
847// ============================================================
848
849#[no_mangle]
850pub extern "C" fn log10(x: f64) -> f64 {
851    let mut hx = (asuint64(x) >> 32) as u32;
852    let mut k: i32 = 0;
853    let mut x = x;
854
855    if hx < 0x00100000 || (hx >> 31) != 0 {
856        if asuint64(x) << 1 == 0 {
857            return -1.0 / (x * x);
858        }
859        if (hx >> 31) != 0 {
860            return (x - x) / 0.0;
861        }
862        k -= 54;
863        x *= 18014398509481984.0;
864        hx = (asuint64(x) >> 32) as u32;
865    } else if hx >= 0x7ff00000 {
866        return x;
867    } else if hx == 0x3ff00000 && (asuint64(x) << 32) == 0 {
868        return 0.0;
869    }
870
871    hx = hx.wrapping_add(0x3ff00000 - 0x3fe6a09e);
872    k += (hx >> 20) as i32 - 0x3ff;
873    hx = (hx & 0x000fffff) + 0x3fe6a09e;
874    x = asdouble(((hx as u64) << 32) | (asuint64(x) & 0xffffffff));
875
876    let f = x - 1.0;
877    let hfsq = 0.5 * f * f;
878    let s = f / (2.0 + f);
879    let z = s * s;
880    let w = z * z;
881    let t1 = w * (LG2 + w * (LG4 + w * LG6));
882    let t2 = z * (LG1 + w * (LG3 + w * (LG5 + w * LG7)));
883    let r = t2 + t1;
884
885    let hi = f - hfsq;
886    let hi = asdouble(asuint64(hi) & (0xFFFFFFFFu64 << 32));
887    let lo = f - hi - hfsq + s * (hfsq + r);
888
889    let val_hi = hi * IVLN10HI;
890    let dk = k as f64;
891    let y = dk * LOG10_2HI;
892    let val_lo = dk * LOG10_2LO + (lo + hi) * IVLN10LO + lo * IVLN10HI;
893
894    let w = y + val_hi;
895    let val_lo = val_lo + (y - w) + val_hi;
896    let val_hi = w;
897
898    val_lo + val_hi
899}
900
901// ============================================================
902// log10f (single precision)
903// ============================================================
904
905#[no_mangle]
906pub extern "C" fn log10f(x: f32) -> f32 {
907    let mut ix = asuint(x);
908    let mut k: i32 = 0;
909    let mut x = x;
910
911    if ix < 0x00800000 || (ix >> 31) != 0 {
912        if ix << 1 == 0 {
913            return -1.0 / (x * x);
914        }
915        if (ix >> 31) != 0 {
916            return (x - x) / 0.0f32;
917        }
918        k -= 25;
919        x *= 33554432.0f32;
920        ix = asuint(x);
921    } else if ix >= 0x7f800000 {
922        return x;
923    } else if ix == 0x3f800000 {
924        return 0.0;
925    }
926
927    ix = ix.wrapping_add(0x3f800000 - 0x3f3504f3);
928    k += (ix >> 23) as i32 - 0x7f;
929    ix = (ix & 0x007fffff) + 0x3f3504f3;
930    x = asfloat(ix);
931
932    let f = x - 1.0f32;
933    let s = f / (2.0f32 + f);
934    let z = s * s;
935    let w = z * z;
936    let t1 = w * (LG2_F + w * LG4_F);
937    let t2 = z * (LG1_F + w * LG3_F);
938    let r = t2 + t1;
939    let hfsq = 0.5f32 * f * f;
940
941    let hi = f - hfsq;
942    let hi = asfloat(asuint(hi) & 0xfffff000);
943    let lo = f - hi - hfsq + s * (hfsq + r);
944    let dk = k as f32;
945
946    dk * LOG10_2LO_F + (lo + hi) * IVLN10LO_F + lo * IVLN10HI_F + hi * IVLN10HI_F + dk * LOG10_2HI_F
947}
c/math_log.rs

c/
math_log.rs
