]>
Commit | Line | Data |
---|---|---|
2aa62f2b | 1 | /* $NetBSD: misc.c,v 1.3.12.1 2008/04/08 21:10:55 jdc Exp $ */\r |
2 | \r | |
3 | /****************************************************************\r | |
4 | \r | |
5 | The author of this software is David M. Gay.\r | |
6 | \r | |
7 | Copyright (C) 1998, 1999 by Lucent Technologies\r | |
8 | All Rights Reserved\r | |
9 | \r | |
10 | Permission to use, copy, modify, and distribute this software and\r | |
11 | its documentation for any purpose and without fee is hereby\r | |
12 | granted, provided that the above copyright notice appear in all\r | |
13 | copies and that both that the copyright notice and this\r | |
14 | permission notice and warranty disclaimer appear in supporting\r | |
15 | documentation, and that the name of Lucent or any of its entities\r | |
16 | not be used in advertising or publicity pertaining to\r | |
17 | distribution of the software without specific, written prior\r | |
18 | permission.\r | |
19 | \r | |
20 | LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,\r | |
21 | INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.\r | |
22 | IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY\r | |
23 | SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES\r | |
24 | WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER\r | |
25 | IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,\r | |
26 | ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF\r | |
27 | THIS SOFTWARE.\r | |
28 | \r | |
29 | ****************************************************************/\r | |
30 | \r | |
31 | /* Please send bug reports to David M. Gay (dmg at acm dot org,\r | |
32 | * with " at " changed at "@" and " dot " changed to "."). */\r | |
33 | #include <LibConfig.h>\r | |
34 | \r | |
35 | #include "gdtoaimp.h"\r | |
36 | \r | |
37 | #if defined(_MSC_VER)\r | |
38 | // Disable warnings about assignment within conditional expressions.\r | |
39 | #pragma warning ( disable : 4706 )\r | |
40 | #endif\r | |
41 | \r | |
42 | static Bigint *freelist[Kmax+1];\r | |
43 | #ifndef Omit_Private_Memory\r | |
44 | #ifndef PRIVATE_MEM\r | |
45 | #define PRIVATE_MEM 2304\r | |
46 | #endif\r | |
47 | #define PRIVATE_mem ((PRIVATE_MEM+sizeof(double)-1)/sizeof(double))\r | |
48 | static double private_mem[PRIVATE_mem], *pmem_next = private_mem;\r | |
49 | #endif\r | |
50 | \r | |
51 | Bigint *\r | |
52 | Balloc\r | |
53 | #ifdef KR_headers\r | |
54 | (k) int k;\r | |
55 | #else\r | |
56 | (int k)\r | |
57 | #endif\r | |
58 | {\r | |
59 | int x;\r | |
60 | Bigint *rv;\r | |
61 | #ifndef Omit_Private_Memory\r | |
62 | unsigned int len;\r | |
63 | #endif\r | |
64 | \r | |
65 | ACQUIRE_DTOA_LOCK(0);\r | |
66 | if ( (rv = freelist[k]) !=0) {\r | |
67 | freelist[k] = rv->next;\r | |
68 | }\r | |
69 | else {\r | |
70 | x = 1 << k;\r | |
71 | #ifdef Omit_Private_Memory\r | |
72 | rv = (Bigint *)MALLOC(sizeof(Bigint) + (x-1)*sizeof(ULong));\r | |
73 | #else\r | |
74 | len = (sizeof(Bigint) + (x-1)*sizeof(ULong) + sizeof(double) - 1)\r | |
75 | /sizeof(double);\r | |
76 | if (pmem_next - private_mem + len <= PRIVATE_mem) {\r | |
77 | rv = (Bigint*)(void *)pmem_next;\r | |
78 | pmem_next += len;\r | |
79 | }\r | |
80 | else\r | |
81 | rv = (Bigint*)MALLOC(len*sizeof(double));\r | |
82 | #endif\r | |
83 | if (rv == NULL)\r | |
84 | return NULL;\r | |
85 | rv->k = k;\r | |
86 | rv->maxwds = x;\r | |
87 | }\r | |
88 | FREE_DTOA_LOCK(0);\r | |
89 | rv->sign = rv->wds = 0;\r | |
90 | return rv;\r | |
91 | }\r | |
92 | \r | |
93 | void\r | |
94 | Bfree\r | |
95 | #ifdef KR_headers\r | |
96 | (v) Bigint *v;\r | |
97 | #else\r | |
98 | (Bigint *v)\r | |
99 | #endif\r | |
100 | {\r | |
101 | if (v) {\r | |
102 | ACQUIRE_DTOA_LOCK(0);\r | |
103 | v->next = freelist[v->k];\r | |
104 | freelist[v->k] = v;\r | |
105 | FREE_DTOA_LOCK(0);\r | |
106 | }\r | |
107 | }\r | |
108 | \r | |
109 | int\r | |
110 | lo0bits\r | |
111 | #ifdef KR_headers\r | |
112 | (y) ULong *y;\r | |
113 | #else\r | |
114 | (ULong *y)\r | |
115 | #endif\r | |
116 | {\r | |
117 | int k;\r | |
118 | ULong x = *y;\r | |
119 | \r | |
120 | if (x & 7) {\r | |
121 | if (x & 1)\r | |
122 | return 0;\r | |
123 | if (x & 2) {\r | |
124 | *y = x >> 1;\r | |
125 | return 1;\r | |
126 | }\r | |
127 | *y = x >> 2;\r | |
128 | return 2;\r | |
129 | }\r | |
130 | k = 0;\r | |
131 | if (!(x & 0xffff)) {\r | |
132 | k = 16;\r | |
133 | x >>= 16;\r | |
134 | }\r | |
135 | if (!(x & 0xff)) {\r | |
136 | k += 8;\r | |
137 | x >>= 8;\r | |
138 | }\r | |
139 | if (!(x & 0xf)) {\r | |
140 | k += 4;\r | |
141 | x >>= 4;\r | |
142 | }\r | |
143 | if (!(x & 0x3)) {\r | |
144 | k += 2;\r | |
145 | x >>= 2;\r | |
146 | }\r | |
147 | if (!(x & 1)) {\r | |
148 | k++;\r | |
149 | x >>= 1;\r | |
150 | if (!x)\r | |
151 | return 32;\r | |
152 | }\r | |
153 | *y = x;\r | |
154 | return k;\r | |
155 | }\r | |
156 | \r | |
157 | Bigint *\r | |
158 | multadd\r | |
159 | #ifdef KR_headers\r | |
160 | (b, m, a) Bigint *b; int m, a;\r | |
161 | #else\r | |
162 | (Bigint *b, int m, int a) /* multiply by m and add a */\r | |
163 | #endif\r | |
164 | {\r | |
165 | int i, wds;\r | |
166 | #ifdef ULLong\r | |
167 | ULong *x;\r | |
168 | ULLong carry, y;\r | |
169 | #else\r | |
170 | ULong carry, *x, y;\r | |
171 | #ifdef Pack_32\r | |
172 | ULong xi, z;\r | |
173 | #endif\r | |
174 | #endif\r | |
175 | Bigint *b1;\r | |
176 | \r | |
177 | wds = b->wds;\r | |
178 | x = b->x;\r | |
179 | i = 0;\r | |
180 | carry = a;\r | |
181 | do {\r | |
182 | #ifdef ULLong\r | |
183 | y = *x * (ULLong)m + carry;\r | |
184 | carry = y >> 32;\r | |
185 | /* LINTED conversion */\r | |
186 | *x++ = (uint32_t)(y & 0xffffffffUL);\r | |
187 | #else\r | |
188 | #ifdef Pack_32\r | |
189 | xi = *x;\r | |
190 | y = (xi & 0xffff) * m + carry;\r | |
191 | z = (xi >> 16) * m + (y >> 16);\r | |
192 | carry = z >> 16;\r | |
193 | *x++ = (z << 16) + (y & 0xffff);\r | |
194 | #else\r | |
195 | y = *x * m + carry;\r | |
196 | carry = y >> 16;\r | |
197 | *x++ = y & 0xffff;\r | |
198 | #endif\r | |
199 | #endif\r | |
200 | }\r | |
201 | while(++i < wds);\r | |
202 | if (carry) {\r | |
203 | if (wds >= b->maxwds) {\r | |
204 | b1 = Balloc(b->k+1);\r | |
205 | if (b1 == NULL) {\r | |
206 | Bfree(b);\r | |
207 | return NULL;\r | |
208 | }\r | |
209 | Bcopy(b1, b);\r | |
210 | Bfree(b);\r | |
211 | b = b1;\r | |
212 | }\r | |
213 | /* LINTED conversion */\r | |
214 | b->x[wds++] = (uint32_t)carry;\r | |
215 | b->wds = wds;\r | |
216 | }\r | |
217 | return b;\r | |
218 | }\r | |
219 | \r | |
220 | int\r | |
221 | hi0bits_D2A\r | |
222 | #ifdef KR_headers\r | |
223 | (x) ULong x;\r | |
224 | #else\r | |
225 | (ULong x)\r | |
226 | #endif\r | |
227 | {\r | |
228 | int k = 0;\r | |
229 | \r | |
230 | if (!(x & 0xffff0000)) {\r | |
231 | k = 16;\r | |
232 | x <<= 16;\r | |
233 | }\r | |
234 | if (!(x & 0xff000000)) {\r | |
235 | k += 8;\r | |
236 | x <<= 8;\r | |
237 | }\r | |
238 | if (!(x & 0xf0000000)) {\r | |
239 | k += 4;\r | |
240 | x <<= 4;\r | |
241 | }\r | |
242 | if (!(x & 0xc0000000)) {\r | |
243 | k += 2;\r | |
244 | x <<= 2;\r | |
245 | }\r | |
246 | if (!(x & 0x80000000)) {\r | |
247 | k++;\r | |
248 | if (!(x & 0x40000000))\r | |
249 | return 32;\r | |
250 | }\r | |
251 | return k;\r | |
252 | }\r | |
253 | \r | |
254 | Bigint *\r | |
255 | i2b\r | |
256 | #ifdef KR_headers\r | |
257 | (i) int i;\r | |
258 | #else\r | |
259 | (int i)\r | |
260 | #endif\r | |
261 | {\r | |
262 | Bigint *b;\r | |
263 | \r | |
264 | b = Balloc(1);\r | |
265 | if (b == NULL)\r | |
266 | return NULL;\r | |
267 | b->x[0] = i;\r | |
268 | b->wds = 1;\r | |
269 | return b;\r | |
270 | }\r | |
271 | \r | |
272 | Bigint *\r | |
273 | mult\r | |
274 | #ifdef KR_headers\r | |
275 | (a, b) Bigint *a, *b;\r | |
276 | #else\r | |
277 | (Bigint *a, Bigint *b)\r | |
278 | #endif\r | |
279 | {\r | |
280 | Bigint *c;\r | |
281 | int k, wa, wb, wc;\r | |
282 | ULong *x, *xa, *xae, *xb, *xbe, *xc, *xc0;\r | |
283 | ULong y;\r | |
284 | #ifdef ULLong\r | |
285 | ULLong carry, z;\r | |
286 | #else\r | |
287 | ULong carry, z;\r | |
288 | #ifdef Pack_32\r | |
289 | ULong z2;\r | |
290 | #endif\r | |
291 | #endif\r | |
292 | \r | |
293 | if (a->wds < b->wds) {\r | |
294 | c = a;\r | |
295 | a = b;\r | |
296 | b = c;\r | |
297 | }\r | |
298 | k = a->k;\r | |
299 | wa = a->wds;\r | |
300 | wb = b->wds;\r | |
301 | wc = wa + wb;\r | |
302 | if (wc > a->maxwds)\r | |
303 | k++;\r | |
304 | c = Balloc(k);\r | |
305 | if (c == NULL)\r | |
306 | return NULL;\r | |
307 | for(x = c->x, xa = x + wc; x < xa; x++)\r | |
308 | *x = 0;\r | |
309 | xa = a->x;\r | |
310 | xae = xa + wa;\r | |
311 | xb = b->x;\r | |
312 | xbe = xb + wb;\r | |
313 | xc0 = c->x;\r | |
314 | #ifdef ULLong\r | |
315 | for(; xb < xbe; xc0++) {\r | |
316 | if ( (y = *xb++) !=0) {\r | |
317 | x = xa;\r | |
318 | xc = xc0;\r | |
319 | carry = 0;\r | |
320 | do {\r | |
321 | z = *x++ * (ULLong)y + *xc + carry;\r | |
322 | carry = z >> 32;\r | |
323 | /* LINTED conversion */\r | |
324 | *xc++ = (uint32_t)(z & 0xffffffffUL);\r | |
325 | }\r | |
326 | while(x < xae);\r | |
327 | /* LINTED conversion */\r | |
328 | *xc = (uint32_t)carry;\r | |
329 | }\r | |
330 | }\r | |
331 | #else\r | |
332 | #ifdef Pack_32\r | |
333 | for(; xb < xbe; xb++, xc0++) {\r | |
334 | if ( (y = *xb & 0xffff) !=0) {\r | |
335 | x = xa;\r | |
336 | xc = xc0;\r | |
337 | carry = 0;\r | |
338 | do {\r | |
339 | z = (*x & 0xffff) * y + (*xc & 0xffff) + carry;\r | |
340 | carry = z >> 16;\r | |
341 | z2 = (*x++ >> 16) * y + (*xc >> 16) + carry;\r | |
342 | carry = z2 >> 16;\r | |
343 | Storeinc(xc, z2, z);\r | |
344 | }\r | |
345 | while(x < xae);\r | |
346 | *xc = carry;\r | |
347 | }\r | |
348 | if ( (y = *xb >> 16) !=0) {\r | |
349 | x = xa;\r | |
350 | xc = xc0;\r | |
351 | carry = 0;\r | |
352 | z2 = *xc;\r | |
353 | do {\r | |
354 | z = (*x & 0xffff) * y + (*xc >> 16) + carry;\r | |
355 | carry = z >> 16;\r | |
356 | Storeinc(xc, z, z2);\r | |
357 | z2 = (*x++ >> 16) * y + (*xc & 0xffff) + carry;\r | |
358 | carry = z2 >> 16;\r | |
359 | }\r | |
360 | while(x < xae);\r | |
361 | *xc = z2;\r | |
362 | }\r | |
363 | }\r | |
364 | #else\r | |
365 | for(; xb < xbe; xc0++) {\r | |
366 | if ( (y = *xb++) !=0) {\r | |
367 | x = xa;\r | |
368 | xc = xc0;\r | |
369 | carry = 0;\r | |
370 | do {\r | |
371 | z = *x++ * y + *xc + carry;\r | |
372 | carry = z >> 16;\r | |
373 | *xc++ = z & 0xffff;\r | |
374 | }\r | |
375 | while(x < xae);\r | |
376 | *xc = carry;\r | |
377 | }\r | |
378 | }\r | |
379 | #endif\r | |
380 | #endif\r | |
381 | for(xc0 = c->x, xc = xc0 + wc; wc > 0 && !*--xc; --wc) ;\r | |
382 | c->wds = wc;\r | |
383 | return c;\r | |
384 | }\r | |
385 | \r | |
386 | static Bigint *p5s;\r | |
387 | \r | |
388 | Bigint *\r | |
389 | pow5mult\r | |
390 | #ifdef KR_headers\r | |
391 | (b, k) Bigint *b; int k;\r | |
392 | #else\r | |
393 | (Bigint *b, int k)\r | |
394 | #endif\r | |
395 | {\r | |
396 | Bigint *b1, *p5, *p51;\r | |
397 | int i;\r | |
398 | static CONST int p05[3] = { 5, 25, 125 };\r | |
399 | \r | |
400 | if ( (i = k & 3) !=0) {\r | |
401 | b = multadd(b, p05[i-1], 0);\r | |
402 | if (b == NULL)\r | |
403 | return NULL;\r | |
404 | }\r | |
405 | \r | |
406 | if ((k = (unsigned int)k >> 2) == 0)\r | |
407 | return b;\r | |
408 | if ((p5 = p5s) == 0) {\r | |
409 | /* first time */\r | |
410 | #ifdef MULTIPLE_THREADS\r | |
411 | ACQUIRE_DTOA_LOCK(1);\r | |
412 | if (!(p5 = p5s)) {\r | |
413 | p5 = p5s = i2b(625);\r | |
414 | if (p5 == NULL)\r | |
415 | return NULL;\r | |
416 | p5->next = 0;\r | |
417 | }\r | |
418 | FREE_DTOA_LOCK(1);\r | |
419 | #else\r | |
420 | p5 = p5s = i2b(625);\r | |
421 | if (p5 == NULL)\r | |
422 | return NULL;\r | |
423 | p5->next = 0;\r | |
424 | #endif\r | |
425 | }\r | |
426 | for(;;) {\r | |
427 | if (k & 1) {\r | |
428 | b1 = mult(b, p5);\r | |
429 | if (b1 == NULL)\r | |
430 | return NULL;\r | |
431 | b = b1;\r | |
432 | }\r | |
433 | if ((k = (unsigned int)k >> 1) == 0)\r | |
434 | break;\r | |
435 | if ((p51 = p5->next) == 0) {\r | |
436 | #ifdef MULTIPLE_THREADS\r | |
437 | ACQUIRE_DTOA_LOCK(1);\r | |
438 | if (!(p51 = p5->next)) {\r | |
439 | p51 = p5->next = mult(p5,p5);\r | |
440 | if (p51 == NULL)\r | |
441 | return NULL;\r | |
442 | p51->next = 0;\r | |
443 | }\r | |
444 | FREE_DTOA_LOCK(1);\r | |
445 | #else\r | |
446 | p51 = p5->next = mult(p5,p5);\r | |
447 | if (p51 == NULL)\r | |
448 | return NULL;\r | |
449 | p51->next = 0;\r | |
450 | #endif\r | |
451 | }\r | |
452 | p5 = p51;\r | |
453 | }\r | |
454 | return b;\r | |
455 | }\r | |
456 | \r | |
457 | Bigint *\r | |
458 | lshift\r | |
459 | #ifdef KR_headers\r | |
460 | (b, k) Bigint *b; int k;\r | |
461 | #else\r | |
462 | (Bigint *b, int k)\r | |
463 | #endif\r | |
464 | {\r | |
465 | int i, k1, n, n1;\r | |
466 | Bigint *b1;\r | |
467 | ULong *x, *x1, *xe, z;\r | |
468 | \r | |
469 | n = (unsigned int)k >> kshift;\r | |
470 | k1 = b->k;\r | |
471 | n1 = n + b->wds + 1;\r | |
472 | for(i = b->maxwds; n1 > i; i <<= 1)\r | |
473 | k1++;\r | |
474 | b1 = Balloc(k1);\r | |
475 | if (b1 == NULL)\r | |
476 | return NULL;\r | |
477 | x1 = b1->x;\r | |
478 | for(i = 0; i < n; i++)\r | |
479 | *x1++ = 0;\r | |
480 | x = b->x;\r | |
481 | xe = x + b->wds;\r | |
482 | if (k &= kmask) {\r | |
483 | #ifdef Pack_32\r | |
484 | k1 = 32 - k;\r | |
485 | z = 0;\r | |
486 | do {\r | |
487 | *x1++ = *x << k | z;\r | |
488 | z = *x++ >> k1;\r | |
489 | }\r | |
490 | while(x < xe);\r | |
491 | if ((*x1 = z) !=0)\r | |
492 | ++n1;\r | |
493 | #else\r | |
494 | k1 = 16 - k;\r | |
495 | z = 0;\r | |
496 | do {\r | |
497 | *x1++ = *x << k & 0xffff | z;\r | |
498 | z = *x++ >> k1;\r | |
499 | }\r | |
500 | while(x < xe);\r | |
501 | if (*x1 = z)\r | |
502 | ++n1;\r | |
503 | #endif\r | |
504 | }\r | |
505 | else do\r | |
506 | *x1++ = *x++;\r | |
507 | while(x < xe);\r | |
508 | b1->wds = n1 - 1;\r | |
509 | Bfree(b);\r | |
510 | return b1;\r | |
511 | }\r | |
512 | \r | |
513 | int\r | |
514 | cmp\r | |
515 | #ifdef KR_headers\r | |
516 | (a, b) Bigint *a, *b;\r | |
517 | #else\r | |
518 | (Bigint *a, Bigint *b)\r | |
519 | #endif\r | |
520 | {\r | |
521 | ULong *xa, *xa0, *xb, *xb0;\r | |
522 | int i, j;\r | |
523 | \r | |
524 | i = a->wds;\r | |
525 | j = b->wds;\r | |
526 | #ifdef DEBUG\r | |
527 | if (i > 1 && !a->x[i-1])\r | |
528 | Bug("cmp called with a->x[a->wds-1] == 0");\r | |
529 | if (j > 1 && !b->x[j-1])\r | |
530 | Bug("cmp called with b->x[b->wds-1] == 0");\r | |
531 | #endif\r | |
532 | if (i -= j)\r | |
533 | return i;\r | |
534 | xa0 = a->x;\r | |
535 | xa = xa0 + j;\r | |
536 | xb0 = b->x;\r | |
537 | xb = xb0 + j;\r | |
538 | for(;;) {\r | |
539 | if (*--xa != *--xb)\r | |
540 | return *xa < *xb ? -1 : 1;\r | |
541 | if (xa <= xa0)\r | |
542 | break;\r | |
543 | }\r | |
544 | return 0;\r | |
545 | }\r | |
546 | \r | |
547 | Bigint *\r | |
548 | diff\r | |
549 | #ifdef KR_headers\r | |
550 | (a, b) Bigint *a, *b;\r | |
551 | #else\r | |
552 | (Bigint *a, Bigint *b)\r | |
553 | #endif\r | |
554 | {\r | |
555 | Bigint *c;\r | |
556 | int i, wa, wb;\r | |
557 | ULong *xa, *xae, *xb, *xbe, *xc;\r | |
558 | #ifdef ULLong\r | |
559 | ULLong borrow, y;\r | |
560 | #else\r | |
561 | ULong borrow, y;\r | |
562 | #ifdef Pack_32\r | |
563 | ULong z;\r | |
564 | #endif\r | |
565 | #endif\r | |
566 | \r | |
567 | i = cmp(a,b);\r | |
568 | if (!i) {\r | |
569 | c = Balloc(0);\r | |
570 | if (c == NULL)\r | |
571 | return NULL;\r | |
572 | c->wds = 1;\r | |
573 | c->x[0] = 0;\r | |
574 | return c;\r | |
575 | }\r | |
576 | if (i < 0) {\r | |
577 | c = a;\r | |
578 | a = b;\r | |
579 | b = c;\r | |
580 | i = 1;\r | |
581 | }\r | |
582 | else\r | |
583 | i = 0;\r | |
584 | c = Balloc(a->k);\r | |
585 | if (c == NULL)\r | |
586 | return NULL;\r | |
587 | c->sign = i;\r | |
588 | wa = a->wds;\r | |
589 | xa = a->x;\r | |
590 | xae = xa + wa;\r | |
591 | wb = b->wds;\r | |
592 | xb = b->x;\r | |
593 | xbe = xb + wb;\r | |
594 | xc = c->x;\r | |
595 | borrow = 0;\r | |
596 | #ifdef ULLong\r | |
597 | do {\r | |
598 | y = (ULLong)*xa++ - *xb++ - borrow;\r | |
599 | borrow = y >> 32 & 1UL;\r | |
600 | /* LINTED conversion */\r | |
601 | *xc++ = (uint32_t)(y & 0xffffffffUL);\r | |
602 | }\r | |
603 | while(xb < xbe);\r | |
604 | while(xa < xae) {\r | |
605 | y = *xa++ - borrow;\r | |
606 | borrow = y >> 32 & 1UL;\r | |
607 | /* LINTED conversion */\r | |
608 | *xc++ = (uint32_t)(y & 0xffffffffUL);\r | |
609 | }\r | |
610 | #else\r | |
611 | #ifdef Pack_32\r | |
612 | do {\r | |
613 | y = (*xa & 0xffff) - (*xb & 0xffff) - borrow;\r | |
614 | borrow = (y & 0x10000) >> 16;\r | |
615 | z = (*xa++ >> 16) - (*xb++ >> 16) - borrow;\r | |
616 | borrow = (z & 0x10000) >> 16;\r | |
617 | Storeinc(xc, z, y);\r | |
618 | }\r | |
619 | while(xb < xbe);\r | |
620 | while(xa < xae) {\r | |
621 | y = (*xa & 0xffff) - borrow;\r | |
622 | borrow = (y & 0x10000) >> 16;\r | |
623 | z = (*xa++ >> 16) - borrow;\r | |
624 | borrow = (z & 0x10000) >> 16;\r | |
625 | Storeinc(xc, z, y);\r | |
626 | }\r | |
627 | #else\r | |
628 | do {\r | |
629 | y = *xa++ - *xb++ - borrow;\r | |
630 | borrow = (y & 0x10000) >> 16;\r | |
631 | *xc++ = y & 0xffff;\r | |
632 | }\r | |
633 | while(xb < xbe);\r | |
634 | while(xa < xae) {\r | |
635 | y = *xa++ - borrow;\r | |
636 | borrow = (y & 0x10000) >> 16;\r | |
637 | *xc++ = y & 0xffff;\r | |
638 | }\r | |
639 | #endif\r | |
640 | #endif\r | |
641 | while(!*--xc)\r | |
642 | wa--;\r | |
643 | c->wds = wa;\r | |
644 | return c;\r | |
645 | }\r | |
646 | \r | |
647 | double\r | |
648 | b2d\r | |
649 | #ifdef KR_headers\r | |
650 | (a, e) Bigint *a; int *e;\r | |
651 | #else\r | |
652 | (Bigint *a, int *e)\r | |
653 | #endif\r | |
654 | {\r | |
655 | ULong *xa, *xa0, w, y, z;\r | |
656 | int k;\r | |
657 | double d;\r | |
658 | #ifdef VAX\r | |
659 | ULong d0, d1;\r | |
660 | #else\r | |
661 | #define d0 word0(d)\r | |
662 | #define d1 word1(d)\r | |
663 | #endif\r | |
664 | \r | |
665 | xa0 = a->x;\r | |
666 | xa = xa0 + a->wds;\r | |
667 | y = *--xa;\r | |
668 | #ifdef DEBUG\r | |
669 | if (!y) Bug("zero y in b2d");\r | |
670 | #endif\r | |
671 | k = hi0bits(y);\r | |
672 | *e = 32 - k;\r | |
673 | #ifdef Pack_32\r | |
674 | if (k < Ebits) {\r | |
675 | d0 = (UINT32)(Exp_1 | y >> (Ebits - k));\r | |
676 | w = xa > xa0 ? *--xa : 0;\r | |
677 | d1 = (UINT32)(y << ((32-Ebits) + k) | w >> (Ebits - k));\r | |
678 | goto ret_d;\r | |
679 | }\r | |
680 | z = xa > xa0 ? *--xa : 0;\r | |
681 | if (k -= Ebits) {\r | |
682 | d0 = (UINT32)(Exp_1 | y << k | z >> (32 - k));\r | |
683 | y = xa > xa0 ? *--xa : 0;\r | |
684 | d1 = (UINT32)(z << k | y >> (32 - k));\r | |
685 | }\r | |
686 | else {\r | |
687 | d0 = (UINT32)(Exp_1 | y);\r | |
688 | d1 = (UINT32)z;\r | |
689 | }\r | |
690 | #else\r | |
691 | if (k < Ebits + 16) {\r | |
692 | z = xa > xa0 ? *--xa : 0;\r | |
693 | d0 = Exp_1 | y << k - Ebits | z >> Ebits + 16 - k;\r | |
694 | w = xa > xa0 ? *--xa : 0;\r | |
695 | y = xa > xa0 ? *--xa : 0;\r | |
696 | d1 = z << k + 16 - Ebits | w << k - Ebits | y >> 16 + Ebits - k;\r | |
697 | goto ret_d;\r | |
698 | }\r | |
699 | z = xa > xa0 ? *--xa : 0;\r | |
700 | w = xa > xa0 ? *--xa : 0;\r | |
701 | k -= Ebits + 16;\r | |
702 | d0 = Exp_1 | y << k + 16 | z << k | w >> 16 - k;\r | |
703 | y = xa > xa0 ? *--xa : 0;\r | |
704 | d1 = w << k + 16 | y << k;\r | |
705 | #endif\r | |
706 | ret_d:\r | |
707 | #ifdef VAX\r | |
708 | word0(d) = d0 >> 16 | d0 << 16;\r | |
709 | word1(d) = d1 >> 16 | d1 << 16;\r | |
710 | #endif\r | |
711 | return dval(d);\r | |
712 | }\r | |
713 | #undef d0\r | |
714 | #undef d1\r | |
715 | \r | |
716 | Bigint *\r | |
717 | d2b\r | |
718 | #ifdef KR_headers\r | |
719 | (d, e, bits) double d; int *e, *bits;\r | |
720 | #else\r | |
721 | (double d, int *e, int *bits)\r | |
722 | #endif\r | |
723 | {\r | |
724 | Bigint *b;\r | |
725 | #ifndef Sudden_Underflow\r | |
726 | int i;\r | |
727 | #endif\r | |
728 | int de, k;\r | |
729 | ULong *x, y, z;\r | |
730 | #ifdef VAX\r | |
731 | ULong d0, d1;\r | |
732 | d0 = word0(d) >> 16 | word0(d) << 16;\r | |
733 | d1 = word1(d) >> 16 | word1(d) << 16;\r | |
734 | #else\r | |
735 | #define d0 word0(d)\r | |
736 | #define d1 word1(d)\r | |
737 | #endif\r | |
738 | \r | |
739 | #ifdef Pack_32\r | |
740 | b = Balloc(1);\r | |
741 | #else\r | |
742 | b = Balloc(2);\r | |
743 | #endif\r | |
744 | if (b == NULL)\r | |
745 | return NULL;\r | |
746 | x = b->x;\r | |
747 | \r | |
748 | z = d0 & Frac_mask;\r | |
749 | d0 &= 0x7fffffff; /* clear sign bit, which we ignore */\r | |
750 | #ifdef Sudden_Underflow\r | |
751 | de = (int)(d0 >> Exp_shift);\r | |
752 | #ifndef IBM\r | |
753 | z |= Exp_msk11;\r | |
754 | #endif\r | |
755 | #else\r | |
756 | if ( (de = (int)(d0 >> Exp_shift)) !=0)\r | |
757 | z |= Exp_msk1;\r | |
758 | #endif\r | |
759 | #ifdef Pack_32\r | |
760 | if ( (y = d1) !=0) {\r | |
761 | if ( (k = lo0bits(&y)) !=0) {\r | |
762 | x[0] = y | z << (32 - k);\r | |
763 | z >>= k;\r | |
764 | }\r | |
765 | else\r | |
766 | x[0] = y;\r | |
767 | #ifndef Sudden_Underflow\r | |
768 | i =\r | |
769 | #endif\r | |
770 | b->wds = (x[1] = z) !=0 ? 2 : 1;\r | |
771 | }\r | |
772 | else {\r | |
773 | #ifdef DEBUG\r | |
774 | if (!z)\r | |
775 | Bug("Zero passed to d2b");\r | |
776 | #endif\r | |
777 | k = lo0bits(&z);\r | |
778 | x[0] = z;\r | |
779 | #ifndef Sudden_Underflow\r | |
780 | i =\r | |
781 | #endif\r | |
782 | b->wds = 1;\r | |
783 | k += 32;\r | |
784 | }\r | |
785 | #else\r | |
786 | if ( (y = d1) !=0) {\r | |
787 | if ( (k = lo0bits(&y)) !=0)\r | |
788 | if (k >= 16) {\r | |
789 | x[0] = y | z << 32 - k & 0xffff;\r | |
790 | x[1] = z >> k - 16 & 0xffff;\r | |
791 | x[2] = z >> k;\r | |
792 | i = 2;\r | |
793 | }\r | |
794 | else {\r | |
795 | x[0] = y & 0xffff;\r | |
796 | x[1] = y >> 16 | z << 16 - k & 0xffff;\r | |
797 | x[2] = z >> k & 0xffff;\r | |
798 | x[3] = z >> k+16;\r | |
799 | i = 3;\r | |
800 | }\r | |
801 | else {\r | |
802 | x[0] = y & 0xffff;\r | |
803 | x[1] = y >> 16;\r | |
804 | x[2] = z & 0xffff;\r | |
805 | x[3] = z >> 16;\r | |
806 | i = 3;\r | |
807 | }\r | |
808 | }\r | |
809 | else {\r | |
810 | #ifdef DEBUG\r | |
811 | if (!z)\r | |
812 | Bug("Zero passed to d2b");\r | |
813 | #endif\r | |
814 | k = lo0bits(&z);\r | |
815 | if (k >= 16) {\r | |
816 | x[0] = z;\r | |
817 | i = 0;\r | |
818 | }\r | |
819 | else {\r | |
820 | x[0] = z & 0xffff;\r | |
821 | x[1] = z >> 16;\r | |
822 | i = 1;\r | |
823 | }\r | |
824 | k += 32;\r | |
825 | }\r | |
826 | while(!x[i])\r | |
827 | --i;\r | |
828 | b->wds = i + 1;\r | |
829 | #endif\r | |
830 | #ifndef Sudden_Underflow\r | |
831 | if (de) {\r | |
832 | #endif\r | |
833 | #ifdef IBM\r | |
834 | *e = (de - Bias - (P-1) << 2) + k;\r | |
835 | *bits = 4*P + 8 - k - hi0bits(word0(d) & Frac_mask);\r | |
836 | #else\r | |
837 | *e = de - Bias - (P-1) + k;\r | |
838 | *bits = P - k;\r | |
839 | #endif\r | |
840 | #ifndef Sudden_Underflow\r | |
841 | }\r | |
842 | else {\r | |
843 | *e = de - Bias - (P-1) + 1 + k;\r | |
844 | #ifdef Pack_32\r | |
845 | *bits = 32*i - hi0bits(x[i-1]);\r | |
846 | #else\r | |
847 | *bits = (i+2)*16 - hi0bits(x[i]);\r | |
848 | #endif\r | |
849 | }\r | |
850 | #endif\r | |
851 | return b;\r | |
852 | }\r | |
853 | #undef d0\r | |
854 | #undef d1\r | |
855 | \r | |
856 | CONST double\r | |
857 | #ifdef IEEE_Arith\r | |
858 | bigtens[] = { 1e16, 1e32, 1e64, 1e128, 1e256 };\r | |
859 | CONST double tinytens[] = { 1e-16, 1e-32, 1e-64, 1e-128, 1e-256\r | |
860 | };\r | |
861 | #else\r | |
862 | #ifdef IBM\r | |
863 | bigtens[] = { 1e16, 1e32, 1e64 };\r | |
864 | CONST double tinytens[] = { 1e-16, 1e-32, 1e-64 };\r | |
865 | #else\r | |
866 | bigtens[] = { 1e16, 1e32 };\r | |
867 | CONST double tinytens[] = { 1e-16, 1e-32 };\r | |
868 | #endif\r | |
869 | #endif\r | |
870 | \r | |
871 | CONST double\r | |
872 | tens[] = {\r | |
873 | 1e0, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9,\r | |
874 | 1e10, 1e11, 1e12, 1e13, 1e14, 1e15, 1e16, 1e17, 1e18, 1e19,\r | |
875 | 1e20, 1e21, 1e22\r | |
876 | #ifdef VAX\r | |
877 | , 1e23, 1e24\r | |
878 | #endif\r | |
879 | };\r | |
880 | \r | |
881 | char *\r | |
882 | #ifdef KR_headers\r | |
883 | strcp_D2A(a, b) char *a; char *b;\r | |
884 | #else\r | |
885 | strcp_D2A(char *a, CONST char *b)\r | |
886 | #endif\r | |
887 | {\r | |
888 | while((*a = *b++))\r | |
889 | a++;\r | |
890 | return a;\r | |
891 | }\r | |
892 | \r | |
893 | #ifdef NO_STRING_H\r | |
894 | \r | |
895 | Char *\r | |
896 | #ifdef KR_headers\r | |
897 | memcpy_D2A(a, b, len) Char *a; Char *b; size_t len;\r | |
898 | #else\r | |
899 | memcpy_D2A(void *a1, void *b1, size_t len)\r | |
900 | #endif\r | |
901 | {\r | |
902 | char *a = (char*)a1, *ae = a + len;\r | |
903 | char *b = (char*)b1, *a0 = a;\r | |
904 | while(a < ae)\r | |
905 | *a++ = *b++;\r | |
906 | return a0;\r | |
907 | }\r | |
908 | \r | |
909 | #endif /* NO_STRING_H */\r |