contrib/gdtoa/misc.c

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