1
2 /*
3 * Copyright (C) Dmitry Volyntsev
4 * Copyright (C) NGINX, Inc.
5 *
6 * An internal diy_fp implementation.
7 * For details, see Loitsch, Florian. "Printing floating-point numbers quickly
8 * and accurately with integers." ACM Sigplan Notices 45.6 (2010): 233-243.
9 */
10
11 #ifndef _NJS_DIYFP_H_INCLUDED_
12 #define _NJS_DIYFP_H_INCLUDED_
13
14
15 typedef struct {
16 uint64_t significand;
17 int exp;
18 } njs_diyfp_t;
19
20
21 typedef union {
22 double d;
23 uint64_t u64;
24 } njs_diyfp_conv_t;
25
26
27 #define njs_diyfp(_s, _e) (njs_diyfp_t) \
28 { .significand = (_s), .exp = (_e) }
29 #define njs_uint64(h, l) (((uint64_t) (h) << 32) + (l))
30
31
32 #define NJS_DBL_SIGNIFICAND_SIZE 52
33 #define NJS_DBL_EXPONENT_OFFSET ((int64_t) 0x3ff)
34 #define NJS_DBL_EXPONENT_BIAS (NJS_DBL_EXPONENT_OFFSET \
35 + NJS_DBL_SIGNIFICAND_SIZE)
36 #define NJS_DBL_EXPONENT_MIN (-NJS_DBL_EXPONENT_BIAS)
37 #define NJS_DBL_EXPONENT_MAX (0x7ff - NJS_DBL_EXPONENT_BIAS)
38 #define NJS_DBL_EXPONENT_DENORMAL (-NJS_DBL_EXPONENT_BIAS + 1)
39
40 #define NJS_DBL_SIGNIFICAND_MASK njs_uint64(0x000FFFFF, 0xFFFFFFFF)
41 #define NJS_DBL_HIDDEN_BIT njs_uint64(0x00100000, 0x00000000)
42 #define NJS_DBL_EXPONENT_MASK njs_uint64(0x7FF00000, 0x00000000)
43 #define NJS_DBL_SIGN_MASK njs_uint64(0x80000000, 0x00000000)
44
45 #define NJS_DIYFP_SIGNIFICAND_SIZE 64
46
47 #define NJS_SIGNIFICAND_SIZE 53
48 #define NJS_SIGNIFICAND_SHIFT (NJS_DIYFP_SIGNIFICAND_SIZE \
49 - NJS_DBL_SIGNIFICAND_SIZE)
50
51 #define NJS_DECIMAL_EXPONENT_OFF 348
52 #define NJS_DECIMAL_EXPONENT_MIN (-348)
53 #define NJS_DECIMAL_EXPONENT_MAX 340
54 #define NJS_DECIMAL_EXPONENT_DIST 8
55
56
57 njs_inline njs_diyfp_t
njs_d2diyfp(double d)58 njs_d2diyfp(double d)
59 {
60 int biased_exp;
61 uint64_t significand;
62 njs_diyfp_t r;
63
64 union {
65 double d;
66 uint64_t u64;
67 } u;
68
69 u.d = d;
70
71 biased_exp = (u.u64 & NJS_DBL_EXPONENT_MASK) >> NJS_DBL_SIGNIFICAND_SIZE;
72 significand = u.u64 & NJS_DBL_SIGNIFICAND_MASK;
73
74 if (biased_exp != 0) {
75 r.significand = significand + NJS_DBL_HIDDEN_BIT;
76 r.exp = biased_exp - NJS_DBL_EXPONENT_BIAS;
77
78 } else {
79 r.significand = significand;
80 r.exp = NJS_DBL_EXPONENT_MIN + 1;
81 }
82
83 return r;
84 }
85
86
87 njs_inline double
njs_diyfp2d(njs_diyfp_t v)88 njs_diyfp2d(njs_diyfp_t v)
89 {
90 int exp;
91 uint64_t significand, biased_exp;
92 njs_diyfp_conv_t conv;
93
94 exp = v.exp;
95 significand = v.significand;
96
97 while (significand > NJS_DBL_HIDDEN_BIT + NJS_DBL_SIGNIFICAND_MASK) {
98 significand >>= 1;
99 exp++;
100 }
101
102 if (exp >= NJS_DBL_EXPONENT_MAX) {
103 return INFINITY;
104 }
105
106 if (exp < NJS_DBL_EXPONENT_DENORMAL) {
107 return 0.0;
108 }
109
110 while (exp > NJS_DBL_EXPONENT_DENORMAL
111 && (significand & NJS_DBL_HIDDEN_BIT) == 0)
112 {
113 significand <<= 1;
114 exp--;
115 }
116
117 if (exp == NJS_DBL_EXPONENT_DENORMAL
118 && (significand & NJS_DBL_HIDDEN_BIT) == 0)
119 {
120 biased_exp = 0;
121
122 } else {
123 biased_exp = (uint64_t) (exp + NJS_DBL_EXPONENT_BIAS);
124 }
125
126 conv.u64 = (significand & NJS_DBL_SIGNIFICAND_MASK)
127 | (biased_exp << NJS_DBL_SIGNIFICAND_SIZE);
128
129 return conv.d;
130 }
131
132
133 njs_inline njs_diyfp_t
njs_diyfp_shift_left(njs_diyfp_t v,unsigned shift)134 njs_diyfp_shift_left(njs_diyfp_t v, unsigned shift)
135 {
136 return njs_diyfp(v.significand << shift, v.exp - shift);
137 }
138
139
140 njs_inline njs_diyfp_t
njs_diyfp_shift_right(njs_diyfp_t v,unsigned shift)141 njs_diyfp_shift_right(njs_diyfp_t v, unsigned shift)
142 {
143 return njs_diyfp(v.significand >> shift, v.exp + shift);
144 }
145
146
147 njs_inline njs_diyfp_t
njs_diyfp_sub(njs_diyfp_t lhs,njs_diyfp_t rhs)148 njs_diyfp_sub(njs_diyfp_t lhs, njs_diyfp_t rhs)
149 {
150 return njs_diyfp(lhs.significand - rhs.significand, lhs.exp);
151 }
152
153
154 njs_inline njs_diyfp_t
njs_diyfp_mul(njs_diyfp_t lhs,njs_diyfp_t rhs)155 njs_diyfp_mul(njs_diyfp_t lhs, njs_diyfp_t rhs)
156 {
157 #if (NJS_HAVE_UNSIGNED_INT128)
158
159 uint64_t l, h;
160 njs_uint128_t u128;
161
162 u128 = (njs_uint128_t) (lhs.significand)
163 * (njs_uint128_t) (rhs.significand);
164
165 h = u128 >> 64;
166 l = (uint64_t) u128;
167
168 /* rounding. */
169
170 if (l & ((uint64_t) 1 << 63)) {
171 h++;
172 }
173
174 return njs_diyfp(h, lhs.exp + rhs.exp + 64);
175
176 #else
177
178 uint64_t a, b, c, d, ac, bc, ad, bd, tmp;
179
180 a = lhs.significand >> 32;
181 b = lhs.significand & 0xffffffff;
182 c = rhs.significand >> 32;
183 d = rhs.significand & 0xffffffff;
184
185 ac = a * c;
186 bc = b * c;
187 ad = a * d;
188 bd = b * d;
189
190 tmp = (bd >> 32) + (ad & 0xffffffff) + (bc & 0xffffffff);
191
192 /* mult_round. */
193
194 tmp += 1U << 31;
195
196 return njs_diyfp(ac + (ad >> 32) + (bc >> 32) + (tmp >> 32),
197 lhs.exp + rhs.exp + 64);
198
199 #endif
200 }
201
202
203 njs_inline njs_diyfp_t
njs_diyfp_normalize(njs_diyfp_t v)204 njs_diyfp_normalize(njs_diyfp_t v)
205 {
206 return njs_diyfp_shift_left(v, njs_leading_zeros64(v.significand));
207 }
208
209
210 njs_diyfp_t njs_cached_power_dec(int exp, int *dec_exp);
211 njs_diyfp_t njs_cached_power_bin(int exp, int *dec_exp);
212
213
214 #endif /* _NJS_DIYFP_H_INCLUDED_ */
215