Rev 1657 | Only display areas with differences | Regard whitespace | Details | Blame | Last modification | View Log | RSS feed
Rev 1657 | Rev 1740 | ||
---|---|---|---|
1 | /* |
1 | /* |
2 | * Copyright (C) 2005 Josef Cejka |
2 | * Copyright (C) 2005 Josef Cejka |
3 | * All rights reserved. |
3 | * All rights reserved. |
4 | * |
4 | * |
5 | * Redistribution and use in source and binary forms, with or without |
5 | * Redistribution and use in source and binary forms, with or without |
6 | * modification, are permitted provided that the following conditions |
6 | * modification, are permitted provided that the following conditions |
7 | * are met: |
7 | * are met: |
8 | * |
8 | * |
9 | * - Redistributions of source code must retain the above copyright |
9 | * - Redistributions of source code must retain the above copyright |
10 | * notice, this list of conditions and the following disclaimer. |
10 | * notice, this list of conditions and the following disclaimer. |
11 | * - Redistributions in binary form must reproduce the above copyright |
11 | * - Redistributions in binary form must reproduce the above copyright |
12 | * notice, this list of conditions and the following disclaimer in the |
12 | * notice, this list of conditions and the following disclaimer in the |
13 | * documentation and/or other materials provided with the distribution. |
13 | * documentation and/or other materials provided with the distribution. |
14 | * - The name of the author may not be used to endorse or promote products |
14 | * - The name of the author may not be used to endorse or promote products |
15 | * derived from this software without specific prior written permission. |
15 | * derived from this software without specific prior written permission. |
16 | * |
16 | * |
17 | * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR |
17 | * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR |
18 | * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES |
18 | * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES |
19 | * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. |
19 | * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. |
20 | * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, |
20 | * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, |
21 | * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT |
21 | * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT |
22 | * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, |
22 | * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, |
23 | * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY |
23 | * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY |
24 | * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT |
24 | * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT |
25 | * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF |
25 | * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF |
26 | * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
26 | * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
27 | */ |
27 | */ |
28 | 28 | ||
29 | /** @addtogroup softfloat |
29 | /** @addtogroup softfloat |
30 | * @{ |
30 | * @{ |
31 | */ |
31 | */ |
32 | /** @file |
32 | /** @file |
33 | */ |
33 | */ |
34 | 34 | ||
35 | #include<sftypes.h> |
35 | #include<sftypes.h> |
36 | #include<add.h> |
36 | #include<add.h> |
37 | #include<div.h> |
37 | #include<div.h> |
38 | #include<comparison.h> |
38 | #include<comparison.h> |
39 | #include<mul.h> |
39 | #include<mul.h> |
40 | #include<common.h> |
40 | #include<common.h> |
41 | 41 | ||
42 | 42 | ||
43 | float32 divFloat32(float32 a, float32 b) |
43 | float32 divFloat32(float32 a, float32 b) |
44 | { |
44 | { |
45 | float32 result; |
45 | float32 result; |
46 | int32_t aexp, bexp, cexp; |
46 | int32_t aexp, bexp, cexp; |
47 | uint64_t afrac, bfrac, cfrac; |
47 | uint64_t afrac, bfrac, cfrac; |
48 | 48 | ||
49 | result.parts.sign = a.parts.sign ^ b.parts.sign; |
49 | result.parts.sign = a.parts.sign ^ b.parts.sign; |
50 | 50 | ||
51 | if (isFloat32NaN(a)) { |
51 | if (isFloat32NaN(a)) { |
52 | if (isFloat32SigNaN(a)) { |
52 | if (isFloat32SigNaN(a)) { |
53 | /*FIXME: SigNaN*/ |
53 | /*FIXME: SigNaN*/ |
54 | } |
54 | } |
55 | /*NaN*/ |
55 | /*NaN*/ |
56 | return a; |
56 | return a; |
57 | } |
57 | } |
58 | 58 | ||
59 | if (isFloat32NaN(b)) { |
59 | if (isFloat32NaN(b)) { |
60 | if (isFloat32SigNaN(b)) { |
60 | if (isFloat32SigNaN(b)) { |
61 | /*FIXME: SigNaN*/ |
61 | /*FIXME: SigNaN*/ |
62 | } |
62 | } |
63 | /*NaN*/ |
63 | /*NaN*/ |
64 | return b; |
64 | return b; |
65 | } |
65 | } |
66 | 66 | ||
67 | if (isFloat32Infinity(a)) { |
67 | if (isFloat32Infinity(a)) { |
68 | if (isFloat32Infinity(b)) { |
68 | if (isFloat32Infinity(b)) { |
69 | /*FIXME: inf / inf */ |
69 | /*FIXME: inf / inf */ |
70 | result.binary = FLOAT32_NAN; |
70 | result.binary = FLOAT32_NAN; |
71 | return result; |
71 | return result; |
72 | } |
72 | } |
73 | /* inf / num */ |
73 | /* inf / num */ |
74 | result.parts.exp = a.parts.exp; |
74 | result.parts.exp = a.parts.exp; |
75 | result.parts.fraction = a.parts.fraction; |
75 | result.parts.fraction = a.parts.fraction; |
76 | return result; |
76 | return result; |
77 | } |
77 | } |
78 | 78 | ||
79 | if (isFloat32Infinity(b)) { |
79 | if (isFloat32Infinity(b)) { |
80 | if (isFloat32Zero(a)) { |
80 | if (isFloat32Zero(a)) { |
81 | /* FIXME 0 / inf */ |
81 | /* FIXME 0 / inf */ |
82 | result.parts.exp = 0; |
82 | result.parts.exp = 0; |
83 | result.parts.fraction = 0; |
83 | result.parts.fraction = 0; |
84 | return result; |
84 | return result; |
85 | } |
85 | } |
86 | /* FIXME: num / inf*/ |
86 | /* FIXME: num / inf*/ |
87 | result.parts.exp = 0; |
87 | result.parts.exp = 0; |
88 | result.parts.fraction = 0; |
88 | result.parts.fraction = 0; |
89 | return result; |
89 | return result; |
90 | } |
90 | } |
91 | 91 | ||
92 | if (isFloat32Zero(b)) { |
92 | if (isFloat32Zero(b)) { |
93 | if (isFloat32Zero(a)) { |
93 | if (isFloat32Zero(a)) { |
94 | /*FIXME: 0 / 0*/ |
94 | /*FIXME: 0 / 0*/ |
95 | result.binary = FLOAT32_NAN; |
95 | result.binary = FLOAT32_NAN; |
96 | return result; |
96 | return result; |
97 | } |
97 | } |
98 | /* FIXME: division by zero */ |
98 | /* FIXME: division by zero */ |
99 | result.parts.exp = 0; |
99 | result.parts.exp = 0; |
100 | result.parts.fraction = 0; |
100 | result.parts.fraction = 0; |
101 | return result; |
101 | return result; |
102 | } |
102 | } |
103 | 103 | ||
104 | 104 | ||
105 | afrac = a.parts.fraction; |
105 | afrac = a.parts.fraction; |
106 | aexp = a.parts.exp; |
106 | aexp = a.parts.exp; |
107 | bfrac = b.parts.fraction; |
107 | bfrac = b.parts.fraction; |
108 | bexp = b.parts.exp; |
108 | bexp = b.parts.exp; |
109 | 109 | ||
110 | /* denormalized numbers */ |
110 | /* denormalized numbers */ |
111 | if (aexp == 0) { |
111 | if (aexp == 0) { |
112 | if (afrac == 0) { |
112 | if (afrac == 0) { |
113 | result.parts.exp = 0; |
113 | result.parts.exp = 0; |
114 | result.parts.fraction = 0; |
114 | result.parts.fraction = 0; |
115 | return result; |
115 | return result; |
116 | } |
116 | } |
117 | /* normalize it*/ |
117 | /* normalize it*/ |
118 | 118 | ||
119 | afrac <<= 1; |
119 | afrac <<= 1; |
120 | /* afrac is nonzero => it must stop */ |
120 | /* afrac is nonzero => it must stop */ |
121 | while (! (afrac & FLOAT32_HIDDEN_BIT_MASK) ) { |
121 | while (! (afrac & FLOAT32_HIDDEN_BIT_MASK) ) { |
122 | afrac <<= 1; |
122 | afrac <<= 1; |
123 | aexp--; |
123 | aexp--; |
124 | } |
124 | } |
125 | } |
125 | } |
126 | 126 | ||
127 | if (bexp == 0) { |
127 | if (bexp == 0) { |
128 | bfrac <<= 1; |
128 | bfrac <<= 1; |
129 | /* bfrac is nonzero => it must stop */ |
129 | /* bfrac is nonzero => it must stop */ |
130 | while (! (bfrac & FLOAT32_HIDDEN_BIT_MASK) ) { |
130 | while (! (bfrac & FLOAT32_HIDDEN_BIT_MASK) ) { |
131 | bfrac <<= 1; |
131 | bfrac <<= 1; |
132 | bexp--; |
132 | bexp--; |
133 | } |
133 | } |
134 | } |
134 | } |
135 | 135 | ||
136 | afrac = (afrac | FLOAT32_HIDDEN_BIT_MASK ) << (32 - FLOAT32_FRACTION_SIZE - 1 ); |
136 | afrac = (afrac | FLOAT32_HIDDEN_BIT_MASK ) << (32 - FLOAT32_FRACTION_SIZE - 1 ); |
137 | bfrac = (bfrac | FLOAT32_HIDDEN_BIT_MASK ) << (32 - FLOAT32_FRACTION_SIZE ); |
137 | bfrac = (bfrac | FLOAT32_HIDDEN_BIT_MASK ) << (32 - FLOAT32_FRACTION_SIZE ); |
138 | 138 | ||
139 | if ( bfrac <= (afrac << 1) ) { |
139 | if ( bfrac <= (afrac << 1) ) { |
140 | afrac >>= 1; |
140 | afrac >>= 1; |
141 | aexp++; |
141 | aexp++; |
142 | } |
142 | } |
143 | 143 | ||
144 | cexp = aexp - bexp + FLOAT32_BIAS - 2; |
144 | cexp = aexp - bexp + FLOAT32_BIAS - 2; |
145 | 145 | ||
146 | cfrac = (afrac << 32) / bfrac; |
146 | cfrac = (afrac << 32) / bfrac; |
147 | if (( cfrac & 0x3F ) == 0) { |
147 | if (( cfrac & 0x3F ) == 0) { |
148 | cfrac |= ( bfrac * cfrac != afrac << 32 ); |
148 | cfrac |= ( bfrac * cfrac != afrac << 32 ); |
149 | } |
149 | } |
150 | 150 | ||
151 | /* pack and round */ |
151 | /* pack and round */ |
152 | 152 | ||
153 | /* find first nonzero digit and shift result and detect possibly underflow */ |
153 | /* find first nonzero digit and shift result and detect possibly underflow */ |
154 | while ((cexp > 0) && (cfrac) && (!(cfrac & (FLOAT32_HIDDEN_BIT_MASK << 7 )))) { |
154 | while ((cexp > 0) && (cfrac) && (!(cfrac & (FLOAT32_HIDDEN_BIT_MASK << 7 )))) { |
155 | cexp--; |
155 | cexp--; |
156 | cfrac <<= 1; |
156 | cfrac <<= 1; |
157 | /* TODO: fix underflow */ |
157 | /* TODO: fix underflow */ |
158 | }; |
158 | }; |
159 | 159 | ||
160 | cfrac += (0x1 << 6); /* FIXME: 7 is not sure*/ |
160 | cfrac += (0x1 << 6); /* FIXME: 7 is not sure*/ |
161 | 161 | ||
162 | if (cfrac & (FLOAT32_HIDDEN_BIT_MASK << 7)) { |
162 | if (cfrac & (FLOAT32_HIDDEN_BIT_MASK << 7)) { |
163 | ++cexp; |
163 | ++cexp; |
164 | cfrac >>= 1; |
164 | cfrac >>= 1; |
165 | } |
165 | } |
166 | 166 | ||
167 | /* check overflow */ |
167 | /* check overflow */ |
168 | if (cexp >= FLOAT32_MAX_EXPONENT ) { |
168 | if (cexp >= FLOAT32_MAX_EXPONENT ) { |
169 | /* FIXME: overflow, return infinity */ |
169 | /* FIXME: overflow, return infinity */ |
170 | result.parts.exp = FLOAT32_MAX_EXPONENT; |
170 | result.parts.exp = FLOAT32_MAX_EXPONENT; |
171 | result.parts.fraction = 0; |
171 | result.parts.fraction = 0; |
172 | return result; |
172 | return result; |
173 | } |
173 | } |
174 | 174 | ||
175 | if (cexp < 0) { |
175 | if (cexp < 0) { |
176 | /* FIXME: underflow */ |
176 | /* FIXME: underflow */ |
177 | result.parts.exp = 0; |
177 | result.parts.exp = 0; |
178 | if ((cexp + FLOAT32_FRACTION_SIZE) < 0) { |
178 | if ((cexp + FLOAT32_FRACTION_SIZE) < 0) { |
179 | result.parts.fraction = 0; |
179 | result.parts.fraction = 0; |
180 | return result; |
180 | return result; |
181 | } |
181 | } |
182 | cfrac >>= 1; |
182 | cfrac >>= 1; |
183 | while (cexp < 0) { |
183 | while (cexp < 0) { |
184 | cexp ++; |
184 | cexp ++; |
185 | cfrac >>= 1; |
185 | cfrac >>= 1; |
186 | } |
186 | } |
187 | 187 | ||
188 | } else { |
188 | } else { |
189 | result.parts.exp = (uint32_t)cexp; |
189 | result.parts.exp = (uint32_t)cexp; |
190 | } |
190 | } |
191 | 191 | ||
192 | result.parts.fraction = ((cfrac >> 6) & (~FLOAT32_HIDDEN_BIT_MASK)); |
192 | result.parts.fraction = ((cfrac >> 6) & (~FLOAT32_HIDDEN_BIT_MASK)); |
193 | 193 | ||
194 | return result; |
194 | return result; |
195 | } |
195 | } |
196 | 196 | ||
197 | float64 divFloat64(float64 a, float64 b) |
197 | float64 divFloat64(float64 a, float64 b) |
198 | { |
198 | { |
199 | float64 result; |
199 | float64 result; |
200 | int64_t aexp, bexp, cexp; |
200 | int64_t aexp, bexp, cexp; |
201 | uint64_t afrac, bfrac, cfrac; |
201 | uint64_t afrac, bfrac, cfrac; |
202 | uint64_t remlo, remhi; |
202 | uint64_t remlo, remhi; |
203 | 203 | ||
204 | result.parts.sign = a.parts.sign ^ b.parts.sign; |
204 | result.parts.sign = a.parts.sign ^ b.parts.sign; |
205 | 205 | ||
206 | if (isFloat64NaN(a)) { |
206 | if (isFloat64NaN(a)) { |
207 | 207 | ||
208 | if (isFloat64SigNaN(b)) { |
208 | if (isFloat64SigNaN(b)) { |
209 | /*FIXME: SigNaN*/ |
209 | /*FIXME: SigNaN*/ |
210 | return b; |
210 | return b; |
211 | } |
211 | } |
212 | 212 | ||
213 | if (isFloat64SigNaN(a)) { |
213 | if (isFloat64SigNaN(a)) { |
214 | /*FIXME: SigNaN*/ |
214 | /*FIXME: SigNaN*/ |
215 | } |
215 | } |
216 | /*NaN*/ |
216 | /*NaN*/ |
217 | return a; |
217 | return a; |
218 | } |
218 | } |
219 | 219 | ||
220 | if (isFloat64NaN(b)) { |
220 | if (isFloat64NaN(b)) { |
221 | if (isFloat64SigNaN(b)) { |
221 | if (isFloat64SigNaN(b)) { |
222 | /*FIXME: SigNaN*/ |
222 | /*FIXME: SigNaN*/ |
223 | } |
223 | } |
224 | /*NaN*/ |
224 | /*NaN*/ |
225 | return b; |
225 | return b; |
226 | } |
226 | } |
227 | 227 | ||
228 | if (isFloat64Infinity(a)) { |
228 | if (isFloat64Infinity(a)) { |
229 | if (isFloat64Infinity(b) || isFloat64Zero(b)) { |
229 | if (isFloat64Infinity(b) || isFloat64Zero(b)) { |
230 | /*FIXME: inf / inf */ |
230 | /*FIXME: inf / inf */ |
231 | result.binary = FLOAT64_NAN; |
231 | result.binary = FLOAT64_NAN; |
232 | return result; |
232 | return result; |
233 | } |
233 | } |
234 | /* inf / num */ |
234 | /* inf / num */ |
235 | result.parts.exp = a.parts.exp; |
235 | result.parts.exp = a.parts.exp; |
236 | result.parts.fraction = a.parts.fraction; |
236 | result.parts.fraction = a.parts.fraction; |
237 | return result; |
237 | return result; |
238 | } |
238 | } |
239 | 239 | ||
240 | if (isFloat64Infinity(b)) { |
240 | if (isFloat64Infinity(b)) { |
241 | if (isFloat64Zero(a)) { |
241 | if (isFloat64Zero(a)) { |
242 | /* FIXME 0 / inf */ |
242 | /* FIXME 0 / inf */ |
243 | result.parts.exp = 0; |
243 | result.parts.exp = 0; |
244 | result.parts.fraction = 0; |
244 | result.parts.fraction = 0; |
245 | return result; |
245 | return result; |
246 | } |
246 | } |
247 | /* FIXME: num / inf*/ |
247 | /* FIXME: num / inf*/ |
248 | result.parts.exp = 0; |
248 | result.parts.exp = 0; |
249 | result.parts.fraction = 0; |
249 | result.parts.fraction = 0; |
250 | return result; |
250 | return result; |
251 | } |
251 | } |
252 | 252 | ||
253 | if (isFloat64Zero(b)) { |
253 | if (isFloat64Zero(b)) { |
254 | if (isFloat64Zero(a)) { |
254 | if (isFloat64Zero(a)) { |
255 | /*FIXME: 0 / 0*/ |
255 | /*FIXME: 0 / 0*/ |
256 | result.binary = FLOAT64_NAN; |
256 | result.binary = FLOAT64_NAN; |
257 | return result; |
257 | return result; |
258 | } |
258 | } |
259 | /* FIXME: division by zero */ |
259 | /* FIXME: division by zero */ |
260 | result.parts.exp = 0; |
260 | result.parts.exp = 0; |
261 | result.parts.fraction = 0; |
261 | result.parts.fraction = 0; |
262 | return result; |
262 | return result; |
263 | } |
263 | } |
264 | 264 | ||
265 | 265 | ||
266 | afrac = a.parts.fraction; |
266 | afrac = a.parts.fraction; |
267 | aexp = a.parts.exp; |
267 | aexp = a.parts.exp; |
268 | bfrac = b.parts.fraction; |
268 | bfrac = b.parts.fraction; |
269 | bexp = b.parts.exp; |
269 | bexp = b.parts.exp; |
270 | 270 | ||
271 | /* denormalized numbers */ |
271 | /* denormalized numbers */ |
272 | if (aexp == 0) { |
272 | if (aexp == 0) { |
273 | if (afrac == 0) { |
273 | if (afrac == 0) { |
274 | result.parts.exp = 0; |
274 | result.parts.exp = 0; |
275 | result.parts.fraction = 0; |
275 | result.parts.fraction = 0; |
276 | return result; |
276 | return result; |
277 | } |
277 | } |
278 | /* normalize it*/ |
278 | /* normalize it*/ |
279 | 279 | ||
280 | aexp++; |
280 | aexp++; |
281 | /* afrac is nonzero => it must stop */ |
281 | /* afrac is nonzero => it must stop */ |
282 | while (! (afrac & FLOAT64_HIDDEN_BIT_MASK) ) { |
282 | while (! (afrac & FLOAT64_HIDDEN_BIT_MASK) ) { |
283 | afrac <<= 1; |
283 | afrac <<= 1; |
284 | aexp--; |
284 | aexp--; |
285 | } |
285 | } |
286 | } |
286 | } |
287 | 287 | ||
288 | if (bexp == 0) { |
288 | if (bexp == 0) { |
289 | bexp++; |
289 | bexp++; |
290 | /* bfrac is nonzero => it must stop */ |
290 | /* bfrac is nonzero => it must stop */ |
291 | while (! (bfrac & FLOAT64_HIDDEN_BIT_MASK) ) { |
291 | while (! (bfrac & FLOAT64_HIDDEN_BIT_MASK) ) { |
292 | bfrac <<= 1; |
292 | bfrac <<= 1; |
293 | bexp--; |
293 | bexp--; |
294 | } |
294 | } |
295 | } |
295 | } |
296 | 296 | ||
297 | afrac = (afrac | FLOAT64_HIDDEN_BIT_MASK ) << (64 - FLOAT64_FRACTION_SIZE - 2 ); |
297 | afrac = (afrac | FLOAT64_HIDDEN_BIT_MASK ) << (64 - FLOAT64_FRACTION_SIZE - 2 ); |
298 | bfrac = (bfrac | FLOAT64_HIDDEN_BIT_MASK ) << (64 - FLOAT64_FRACTION_SIZE - 1); |
298 | bfrac = (bfrac | FLOAT64_HIDDEN_BIT_MASK ) << (64 - FLOAT64_FRACTION_SIZE - 1); |
299 | 299 | ||
300 | if ( bfrac <= (afrac << 1) ) { |
300 | if ( bfrac <= (afrac << 1) ) { |
301 | afrac >>= 1; |
301 | afrac >>= 1; |
302 | aexp++; |
302 | aexp++; |
303 | } |
303 | } |
304 | 304 | ||
305 | cexp = aexp - bexp + FLOAT64_BIAS - 2; |
305 | cexp = aexp - bexp + FLOAT64_BIAS - 2; |
306 | 306 | ||
307 | cfrac = divFloat64estim(afrac, bfrac); |
307 | cfrac = divFloat64estim(afrac, bfrac); |
308 | 308 | ||
309 | if (( cfrac & 0x1FF ) <= 2) { /*FIXME:?? */ |
309 | if (( cfrac & 0x1FF ) <= 2) { /*FIXME:?? */ |
310 | mul64integers( bfrac, cfrac, &remlo, &remhi); |
310 | mul64integers( bfrac, cfrac, &remlo, &remhi); |
311 | /* (__u128)afrac << 64 - ( ((__u128)remhi<<64) + (__u128)remlo )*/ |
311 | /* (__u128)afrac << 64 - ( ((__u128)remhi<<64) + (__u128)remlo )*/ |
312 | remhi = afrac - remhi - ( remlo > 0); |
312 | remhi = afrac - remhi - ( remlo > 0); |
313 | remlo = - remlo; |
313 | remlo = - remlo; |
314 | 314 | ||
315 | while ((int64_t) remhi < 0) { |
315 | while ((int64_t) remhi < 0) { |
316 | cfrac--; |
316 | cfrac--; |
317 | remlo += bfrac; |
317 | remlo += bfrac; |
318 | remhi += ( remlo < bfrac ); |
318 | remhi += ( remlo < bfrac ); |
319 | } |
319 | } |
320 | cfrac |= ( remlo != 0 ); |
320 | cfrac |= ( remlo != 0 ); |
321 | } |
321 | } |
322 | 322 | ||
323 | /* round and shift */ |
323 | /* round and shift */ |
324 | result = finishFloat64(cexp, cfrac, result.parts.sign); |
324 | result = finishFloat64(cexp, cfrac, result.parts.sign); |
325 | return result; |
325 | return result; |
326 | 326 | ||
327 | } |
327 | } |
328 | 328 | ||
329 | uint64_t divFloat64estim(uint64_t a, uint64_t b) |
329 | uint64_t divFloat64estim(uint64_t a, uint64_t b) |
330 | { |
330 | { |
331 | uint64_t bhi; |
331 | uint64_t bhi; |
332 | uint64_t remhi, remlo; |
332 | uint64_t remhi, remlo; |
333 | uint64_t result; |
333 | uint64_t result; |
334 | 334 | ||
335 | if ( b <= a ) { |
335 | if ( b <= a ) { |
336 | return 0xFFFFFFFFFFFFFFFFull; |
336 | return 0xFFFFFFFFFFFFFFFFull; |
337 | } |
337 | } |
338 | 338 | ||
339 | bhi = b >> 32; |
339 | bhi = b >> 32; |
340 | result = ((bhi << 32) <= a) ?( 0xFFFFFFFFull << 32) : ( a / bhi) << 32; |
340 | result = ((bhi << 32) <= a) ?( 0xFFFFFFFFull << 32) : ( a / bhi) << 32; |
341 | mul64integers(b, result, &remlo, &remhi); |
341 | mul64integers(b, result, &remlo, &remhi); |
342 | 342 | ||
343 | remhi = a - remhi - (remlo > 0); |
343 | remhi = a - remhi - (remlo > 0); |
344 | remlo = - remlo; |
344 | remlo = - remlo; |
345 | 345 | ||
346 | b <<= 32; |
346 | b <<= 32; |
347 | while ( (int64_t) remhi < 0 ) { |
347 | while ( (int64_t) remhi < 0 ) { |
348 | result -= 0x1ll << 32; |
348 | result -= 0x1ll << 32; |
349 | remlo += b; |
349 | remlo += b; |
350 | remhi += bhi + ( remlo < b ); |
350 | remhi += bhi + ( remlo < b ); |
351 | } |
351 | } |
352 | remhi = (remhi << 32) | (remlo >> 32); |
352 | remhi = (remhi << 32) | (remlo >> 32); |
353 | if (( bhi << 32) <= remhi) { |
353 | if (( bhi << 32) <= remhi) { |
354 | result |= 0xFFFFFFFF; |
354 | result |= 0xFFFFFFFF; |
355 | } else { |
355 | } else { |
356 | result |= remhi / bhi; |
356 | result |= remhi / bhi; |
357 | } |
357 | } |
358 | 358 | ||
359 | 359 | ||
360 | return result; |
360 | return result; |
361 | } |
361 | } |
362 | 362 | ||
363 | - | ||
364 | /** @} |
363 | /** @} |
365 | */ |
364 | */ |
366 | - | ||
367 | 365 |