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