lib/libc/gdtoa/misc.c

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