2  * Double-precision log2(x) function.
4  * Copyright (c) 2018-2024, Arm Limited.
5  * SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
22 #define OFF 0x3fe6000000000000
44 # define LO asuint64 (1.0 - 0x1.5b51p-5)  in log2()
45 # define HI asuint64 (1.0 + 0x1.6ab2p-5)  in log2()
47   if (unlikely (ix - LO < HI - LO))  in log2()
52 	return 0;  in log2()
53       r = x - 1.0;  in log2()
56       lo = r * InvLn2lo + fma (r, InvLn2hi, -hi);  in log2()
59       rhi = asdouble (asuint64 (r) & -1ULL << 32);  in log2()
60       rlo = r - rhi;  in log2()
64       r2 = r * r; /* rounding error: 0x1p-62.  */  in log2()
67       /* Worst-case error is less than 0.54 ULP (0.55 ULP without fma).  */  in log2()
68       p = r2 * (B[0] + r * B[1]);  in log2()
70       lo += hi - y + p;  in log2()
77   if (unlikely (top - 0x0010 >= 0x7ff0 - 0x0010))  in log2()
79       /* x < 0x1p-1022 or inf or nan.  */  in log2()
80       if (ix * 2 == 0)  in log2()
84       if ((top & 0x8000) || (top & 0x7ff0) == 0x7ff0)  in log2()
87       ix = asuint64 (x * 0x1p52);  in log2()
88       ix -= 52ULL << 52;  in log2()
94   tmp = ix - OFF;  in log2()
95   i = (tmp >> (52 - LOG2_TABLE_BITS)) % N;  in log2()
97   iz = ix - (tmp & 0xfffULL << 52);  in log2()
104   /* r ~= z/c - 1, |r| < 1/(2*N).  */  in log2()
106   /* rounding error: 0x1p-55/N.  */  in log2()
107   r = fma (z, invc, -1.0);  in log2()
109   t2 = r * InvLn2lo + fma (r, InvLn2hi, -t1);  in log2()
112   /* rounding error: 0x1p-55/N + 0x1p-65.  */  in log2()
113   r = (z - T2[i].chi - T2[i].clo) * invc;  in log2()
114   rhi = asdouble (asuint64 (r) & -1ULL << 32);  in log2()
115   rlo = r - rhi;  in log2()
123   lo = t3 - hi + t1 + t2;  in log2()
127   r2 = r * r; /* rounding error: 0x1p-54/N^2.  */  in log2()
130   /* Worst-case error if |y| > 0x1p-4: 0.547 ULP (0.550 ULP without fma).  in log2()
131      ~ 0.5 + 2/N/ln2 + abs-poly-error*0x1p56 ULP (+ 0.003 ULP without fma).  */  in log2()
132   p = A[0] + r * A[1] + r2 * (A[2] + r * A[3]) + r4 * (A[4] + r * A[5]);  in log2()
148 TEST_INTERVAL (log2, 0, 0xffff000000000000, 10000)
149 TEST_INTERVAL (log2, 0x1p-4, 0x1p4, 40000)
150 TEST_INTERVAL (log2, 0, inf, 40000)