Home | History | Annotate | Line # | Download | only in gdtoa
      1 /* $NetBSD: dmisc.c,v 1.4 2008/03/21 23:13:48 christos Exp $ */
      2 
      3 /****************************************************************
      4 
      5 The author of this software is David M. Gay.
      6 
      7 Copyright (C) 1998 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 #ifndef MULTIPLE_THREADS
     37  char *dtoa_result;
     38 #endif
     39 
     40  char *
     41 #ifdef KR_headers
     42 rv_alloc(i) size_t i;
     43 #else
     44 rv_alloc(size_t i)
     45 #endif
     46 {
     47 	size_t j;
     48 	int k, *r;
     49 
     50 	j = sizeof(ULong);
     51 	for(k = 0;
     52 		sizeof(Bigint) - sizeof(ULong) - sizeof(int) + j <= i;
     53 		j <<= 1)
     54 			k++;
     55 	r = (int*)(void*)Balloc(k);
     56 	if (r == NULL)
     57 		return NULL;
     58 	*r = k;
     59 	return
     60 #ifndef MULTIPLE_THREADS
     61 	dtoa_result =
     62 #endif
     63 		(char *)(void *)(r+1);
     64 	}
     65 
     66  char *
     67 #ifdef KR_headers
     68 nrv_alloc(s, rve, n) CONST char *s; char **rve; size_t n;
     69 #else
     70 nrv_alloc(CONST char *s, char **rve, size_t n)
     71 #endif
     72 {
     73 	char *rv, *t;
     74 
     75 	t = rv = rv_alloc(n);
     76 	if (t == NULL)
     77 		return NULL;
     78 	while((*t = *s++) !=0)
     79 		t++;
     80 	if (rve)
     81 		*rve = t;
     82 	return rv;
     83 	}
     84 
     85 /* freedtoa(s) must be used to free values s returned by dtoa
     86  * when MULTIPLE_THREADS is #defined.  It should be used in all cases,
     87  * but for consistency with earlier versions of dtoa, it is optional
     88  * when MULTIPLE_THREADS is not defined.
     89  */
     90 
     91  void
     92 #ifdef KR_headers
     93 freedtoa(s) char *s;
     94 #else
     95 freedtoa(char *s)
     96 #endif
     97 {
     98 	Bigint *b = (Bigint *)(void *)((int *)(void *)s - 1);
     99 	b->maxwds = 1 << (b->k = *(int*)(void*)b);
    100 	Bfree(b);
    101 #ifndef MULTIPLE_THREADS
    102 	if (s == dtoa_result)
    103 		dtoa_result = 0;
    104 #endif
    105 	}
    106 
    107  int
    108 quorem
    109 #ifdef KR_headers
    110 	(b, S) Bigint *b, *S;
    111 #else
    112 	(Bigint *b, Bigint *S)
    113 #endif
    114 {
    115 	int n;
    116 	ULong *bx, *bxe, q, *sx, *sxe;
    117 #ifdef ULLong
    118 	ULLong borrow, carry, y, ys;
    119 #else
    120 	ULong borrow, carry, y, ys;
    121 #ifdef Pack_32
    122 	ULong si, z, zs;
    123 #endif
    124 #endif
    125 
    126 	n = S->wds;
    127 #ifdef DEBUG
    128 	/*debug*/ if (b->wds > n)
    129 	/*debug*/	Bug("oversize b in quorem");
    130 #endif
    131 	if (b->wds < n)
    132 		return 0;
    133 	sx = S->x;
    134 	sxe = sx + --n;
    135 	bx = b->x;
    136 	bxe = bx + n;
    137 	q = *bxe / (*sxe + 1);	/* ensure q <= true quotient */
    138 #ifdef DEBUG
    139 	/*debug*/ if (q > 9)
    140 	/*debug*/	Bug("oversized quotient in quorem");
    141 #endif
    142 	if (q) {
    143 		borrow = 0;
    144 		carry = 0;
    145 		do {
    146 #ifdef ULLong
    147 			ys = *sx++ * (ULLong)q + carry;
    148 			carry = ys >> 32;
    149 			/* LINTED conversion */
    150 			y = *bx - (ys & 0xffffffffUL) - borrow;
    151 			borrow = y >> 32 & 1UL;
    152 			/* LINTED conversion */
    153 			*bx++ = y & 0xffffffffUL;
    154 #else
    155 #ifdef Pack_32
    156 			si = *sx++;
    157 			ys = (si & 0xffff) * q + carry;
    158 			zs = (si >> 16) * q + (ys >> 16);
    159 			carry = zs >> 16;
    160 			y = (*bx & 0xffff) - (ys & 0xffff) - borrow;
    161 			borrow = (y & 0x10000) >> 16;
    162 			z = (*bx >> 16) - (zs & 0xffff) - borrow;
    163 			borrow = (z & 0x10000) >> 16;
    164 			Storeinc(bx, z, y);
    165 #else
    166 			ys = *sx++ * q + carry;
    167 			carry = ys >> 16;
    168 			y = *bx - (ys & 0xffff) - borrow;
    169 			borrow = (y & 0x10000) >> 16;
    170 			*bx++ = y & 0xffff;
    171 #endif
    172 #endif
    173 			}
    174 			while(sx <= sxe);
    175 		if (!*bxe) {
    176 			bx = b->x;
    177 			while(--bxe > bx && !*bxe)
    178 				--n;
    179 			b->wds = n;
    180 			}
    181 		}
    182 	if (cmp(b, S) >= 0) {
    183 		q++;
    184 		borrow = 0;
    185 		carry = 0;
    186 		bx = b->x;
    187 		sx = S->x;
    188 		do {
    189 #ifdef ULLong
    190 			ys = *sx++ + carry;
    191 			carry = ys >> 32;
    192 			/* LINTED conversion */
    193 			y = *bx - (ys & 0xffffffffUL) - borrow;
    194 			borrow = y >> 32 & 1UL;
    195 			/* LINTED conversion */
    196 			*bx++ = y & 0xffffffffUL;
    197 #else
    198 #ifdef Pack_32
    199 			si = *sx++;
    200 			ys = (si & 0xffff) + carry;
    201 			zs = (si >> 16) + (ys >> 16);
    202 			carry = zs >> 16;
    203 			y = (*bx & 0xffff) - (ys & 0xffff) - borrow;
    204 			borrow = (y & 0x10000) >> 16;
    205 			z = (*bx >> 16) - (zs & 0xffff) - borrow;
    206 			borrow = (z & 0x10000) >> 16;
    207 			Storeinc(bx, z, y);
    208 #else
    209 			ys = *sx++ + carry;
    210 			carry = ys >> 16;
    211 			y = *bx - (ys & 0xffff) - borrow;
    212 			borrow = (y & 0x10000) >> 16;
    213 			*bx++ = y & 0xffff;
    214 #endif
    215 #endif
    216 			}
    217 			while(sx <= sxe);
    218 		bx = b->x;
    219 		bxe = bx + n;
    220 		if (!*bxe) {
    221 			while(--bxe > bx && !*bxe)
    222 				--n;
    223 			b->wds = n;
    224 			}
    225 		}
    226 	return q;
    227 	}
    228