master
  1/****************************************************************
  2
  3The author of this software is David M. Gay.
  4
  5Copyright (C) 1998 by Lucent Technologies
  6All Rights Reserved
  7
  8Permission to use, copy, modify, and distribute this software and
  9its documentation for any purpose and without fee is hereby
 10granted, provided that the above copyright notice appear in all
 11copies and that both that the copyright notice and this
 12permission notice and warranty disclaimer appear in supporting
 13documentation, and that the name of Lucent or any of its entities
 14not be used in advertising or publicity pertaining to
 15distribution of the software without specific, written prior
 16permission.
 17
 18LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
 19INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.
 20IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY
 21SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
 22WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
 23IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
 24ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
 25THIS 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#ifndef MULTIPLE_THREADS
 35char *dtoa_result;
 36#endif
 37
 38char *rv_alloc (int i)
 39{
 40	int j, k, *r;
 41
 42	j = sizeof(ULong);
 43	for(k = 0;
 44		(int) (sizeof(Bigint) - sizeof(ULong) - sizeof(int)) + j <= i;
 45		j <<= 1)
 46			k++;
 47	r = (int*)Balloc(k);
 48	*r = k;
 49	return
 50#ifndef MULTIPLE_THREADS
 51	dtoa_result =
 52#endif
 53		(char *)(r+1);
 54}
 55
 56char *nrv_alloc (char *s, char **rve, int n)
 57{
 58	char *rv, *t;
 59
 60	t = rv = rv_alloc(n);
 61	while((*t = *s++) !=0)
 62		t++;
 63	if (rve)
 64		*rve = t;
 65	return rv;
 66}
 67
 68/* freedtoa(s) must be used to free values s returned by dtoa
 69 * when MULTIPLE_THREADS is #defined.  It should be used in all cases,
 70 * but for consistency with earlier versions of dtoa, it is optional
 71 * when MULTIPLE_THREADS is not defined.
 72 */
 73
 74void __freedtoa (char *s)
 75{
 76	Bigint *b = (Bigint *)((int *)s - 1);
 77	b->maxwds = 1 << (b->k = *(int*)b);
 78	Bfree(b);
 79#ifndef MULTIPLE_THREADS
 80	if (s == dtoa_result)
 81		dtoa_result = 0;
 82#endif
 83}
 84
 85int quorem (Bigint *b, Bigint *S)
 86{
 87	int n;
 88	ULong *bx, *bxe, q, *sx, *sxe;
 89#ifdef ULLong
 90	ULLong borrow, carry, y, ys;
 91#else
 92	ULong borrow, carry, y, ys;
 93#ifdef Pack_32
 94	ULong si, z, zs;
 95#endif
 96#endif
 97
 98	n = S->wds;
 99#ifdef DEBUG
100	/*debug*/ if (b->wds > n)
101	/*debug*/	Bug("oversize b in quorem");
102#endif
103	if (b->wds < n)
104		return 0;
105	sx = S->x;
106	sxe = sx + --n;
107	bx = b->x;
108	bxe = bx + n;
109	q = *bxe / (*sxe + 1);	/* ensure q <= true quotient */
110#ifdef DEBUG
111	/*debug*/ if (q > 9)
112	/*debug*/	Bug("oversized quotient in quorem");
113#endif
114	if (q) {
115		borrow = 0;
116		carry = 0;
117		do {
118#ifdef ULLong
119			ys = *sx++ * (ULLong)q + carry;
120			carry = ys >> 32;
121			y = *bx - (ys & 0xffffffffUL) - borrow;
122			borrow = y >> 32 & 1UL;
123			*bx++ = y & 0xffffffffUL;
124#else
125#ifdef Pack_32
126			si = *sx++;
127			ys = (si & 0xffff) * q + carry;
128			zs = (si >> 16) * q + (ys >> 16);
129			carry = zs >> 16;
130			y = (*bx & 0xffff) - (ys & 0xffff) - borrow;
131			borrow = (y & 0x10000) >> 16;
132			z = (*bx >> 16) - (zs & 0xffff) - borrow;
133			borrow = (z & 0x10000) >> 16;
134			Storeinc(bx, z, y);
135#else
136			ys = *sx++ * q + carry;
137			carry = ys >> 16;
138			y = *bx - (ys & 0xffff) - borrow;
139			borrow = (y & 0x10000) >> 16;
140			*bx++ = y & 0xffff;
141#endif
142#endif
143		} while(sx <= sxe);
144
145		if (!*bxe) {
146			bx = b->x;
147			while(--bxe > bx && !*bxe)
148				--n;
149			b->wds = n;
150		}
151	}
152
153	if (cmp(b, S) >= 0) {
154		q++;
155		borrow = 0;
156		carry = 0;
157		bx = b->x;
158		sx = S->x;
159		do {
160#ifdef ULLong
161			ys = *sx++ + carry;
162			carry = ys >> 32;
163			y = *bx - (ys & 0xffffffffUL) - borrow;
164			borrow = y >> 32 & 1UL;
165			*bx++ = y & 0xffffffffUL;
166#else
167#ifdef Pack_32
168			si = *sx++;
169			ys = (si & 0xffff) + carry;
170			zs = (si >> 16) + (ys >> 16);
171			carry = zs >> 16;
172			y = (*bx & 0xffff) - (ys & 0xffff) - borrow;
173			borrow = (y & 0x10000) >> 16;
174			z = (*bx >> 16) - (zs & 0xffff) - borrow;
175			borrow = (z & 0x10000) >> 16;
176			Storeinc(bx, z, y);
177#else
178			ys = *sx++ + carry;
179			carry = ys >> 16;
180			y = *bx - (ys & 0xffff) - borrow;
181			borrow = (y & 0x10000) >> 16;
182			*bx++ = y & 0xffff;
183#endif
184#endif
185		} while(sx <= sxe);
186
187		bx = b->x;
188		bxe = bx + n;
189		if (!*bxe) {
190			while(--bxe > bx && !*bxe)
191				--n;
192			b->wds = n;
193		}
194	}
195	return q;
196}