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