/usr.bin/bc/bc.library

  1/*	$FreeBSD$							*/
  2/*      $OpenBSD: bc.library,v 1.3 2007/02/03 21:15:06 otto Exp $	*/
  3
  4/*
  5 * Copyright (C) Caldera International Inc.  2001-2002.
  6 * All rights reserved.
  7 *
  8 * Redistribution and use in source and binary forms, with or without
  9 * modification, are permitted provided that the following conditions
 10 * are met:
 11 * 1. Redistributions of source code and documentation must retain the above
 12 *    copyright notice, this list of conditions and the following disclaimer.
 13 * 2. Redistributions in binary form must reproduce the above copyright
 14 *    notice, this list of conditions and the following disclaimer in the
 15 *    documentation and/or other materials provided with the distribution.
 16 * 3. All advertising materials mentioning features or use of this software
 17 *    must display the following acknowledgement:
 18 *      This product includes software developed or owned by Caldera
 19 *      International, Inc.
 20 * 4. Neither the name of Caldera International, Inc. nor the names of other
 21 *    contributors may be used to endorse or promote products derived from
 22 *    this software without specific prior written permission.
 23 *
 24 * USE OF THE SOFTWARE PROVIDED FOR UNDER THIS LICENSE BY CALDERA
 25 * INTERNATIONAL, INC. AND CONTRIBUTORS ``AS IS'' AND ANY EXPRESS OR
 26 * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
 27 * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
 28 * IN NO EVENT SHALL CALDERA INTERNATIONAL, INC. BE LIABLE FOR ANY DIRECT,
 29 * INDIRECT INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
 30 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 31 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
 32 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
 33 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING
 34 * IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 35 * POSSIBILITY OF SUCH DAMAGE.
 36 */
 37
 38/*
 39 *	@(#)bc.library	5.1 (Berkeley) 4/17/91
 40 */
 41
 42scale = 20
 43define e(x) {
 44	auto a, b, c, d, e, g, t, w, y, r
 45
 46	r = ibase
 47	ibase = A
 48	t = scale
 49	scale = 0
 50	if (x > 0) scale = (0.435*x)/1
 51	scale = scale + t + length(scale + t) + 1
 52
 53	w = 0
 54	if (x < 0) {
 55		x = -x
 56		w = 1
 57	}
 58	y = 0
 59	while (x > 2) {
 60		x = x/2
 61		y = y + 1
 62	}
 63
 64	a = 1
 65	b = 1
 66	c = b
 67	d = 1
 68	e = 1
 69	for (a = 1; 1 == 1; a++) {
 70		b = b*x
 71		c = c*a + b
 72		d = d*a
 73		g = c/d
 74		if (g == e) {
 75			g = g/1
 76			while (y--) {
 77				g = g*g
 78			}
 79			scale = t
 80			ibase = r
 81			if (w == 1) return (1/g)
 82			return (g/1)
 83		}
 84		e = g
 85	}
 86}
 87
 88define l(x) {
 89	auto a, b, c, d, e, f, g, u, s, t, r
 90	r = ibase
 91	ibase = A
 92	if (x <= 0) {
 93		a = (1 - 10^scale)
 94		ibase = r
 95		return (a)
 96	}
 97	t = scale
 98
 99	f = 1
100	if (x < 1) {
101		s = scale(x)
102	} else {
103		s = length(x) - scale(x)
104	}
105	scale = 0
106	a = (2.31*s)/1 /* estimated integer part of the answer */
107	s = t + length(a) + 2 /* estimated length of the answer */
108	while (x > 2) {
109		scale=0
110		scale = (length(x) + scale(x))/2 + 1
111		if (scale < s) scale = s
112		x = sqrt(x)
113		f = f*2
114	}
115	while (x < .5) {
116		scale = 0
117		scale = scale(x)/2 + 1
118		if (scale < s) scale = s
119		x = sqrt(x)
120		f = f*2
121	}
122
123	scale = 0
124	scale = t + length(f) + length((1.05*(t+length(f))/1)) + 1
125	u = (x - 1)/(x + 1)
126	s = u*u
127	scale = t + 2
128	b = 2*f
129	c = b
130	d = 1
131	e = 1
132	for (a = 3; 1 == 1 ; a = a + 2) {
133		b = b*s
134		c = c*a + d*b
135		d = d*a
136		g = c/d
137		if (g == e) {
138			scale = t
139			ibase = r
140			return (u*c/d)
141		}
142		e = g
143	}
144}
145
146define s(x) {
147	auto a, b, c, s, t, y, p, n, i, r
148	r = ibase
149	ibase = A
150	t = scale
151	y = x/.7853
152	s = t + length(y) - scale(y)
153	if (s < t) s = t
154	scale = s
155	p = a(1)
156
157	scale = 0
158	if (x >= 0) n = (x/(2*p) + 1)/2
159	if (x < 0) n = (x/(2*p) - 1)/2
160	x = x - 4*n*p
161	if (n % 2 != 0) x = -x
162
163	scale = t + length(1.2*t) - scale(1.2*t)
164	y = -x*x
165	a = x
166	b = 1
167	s = x
168	for (i =3 ; 1 == 1; i = i + 2) {
169		a = a*y
170		b = b*i*(i - 1)
171		c = a/b
172		if (c == 0) {
173			scale = t
174			ibase = r
175			return (s/1)
176		}
177		s = s + c
178	}
179}
180
181define c(x) {
182	auto t, r
183	r = ibase
184	ibase = A
185	t = scale
186	scale = scale + 1
187	x = s(x + 2*a(1))
188	scale = t
189	ibase = r
190	return (x/1)
191}
192
193define a(x) {
194	auto a, b, c, d, e, f, g, s, t, r
195	if (x == 0) return(0)
196
197	r = ibase
198	ibase = A
199	if (x == 1) {
200		if (scale < 52) {
201			 a = .7853981633974483096156608458198757210492923498437764/1
202			 ibase = r
203			 return (a)
204		}
205	}
206	t = scale
207	f = 1
208	while (x > .5) {
209		scale = scale + 1
210		x = -(1 - sqrt(1. + x*x))/x
211		f = f*2
212	}
213	while (x < -.5) {
214		scale = scale + 1
215		x = -(1 - sqrt(1. + x*x))/x
216		f = f*2
217	}
218	s = -x*x
219	b = f
220	c = f
221	d = 1
222	e = 1
223	for (a = 3; 1 == 1; a = a + 2) {
224		b = b*s
225		c = c*a + d*b
226		d = d*a
227		g = c/d
228		if (g == e) {
229			ibase = r
230			scale = t
231			return (x*c/d)
232		}
233		e = g
234	}
235}
236
237define j(n,x) {
238	auto a, b, c, d, e, g, i, s, k, t, r
239
240	r = ibase
241	ibase = A
242	t = scale
243	k = 1.36*x + 1.16*t - n
244	k = length(k) - scale(k)
245	if (k > 0) scale = scale + k
246
247	s = -x*x/4
248	if (n < 0) {
249		n = -n
250		x = -x
251	}
252	a = 1
253	c = 1
254	for (i = 1; i <= n; i++) {
255		a = a*x
256		c = c*2*i
257	}
258	b = a
259	d = 1
260	e = 1
261	for (i = 1; 1; i++) {
262		a = a*s
263		b = b*i*(n + i) + a
264		c = c*i*(n + i)
265		g = b/c
266		if (g == e) {
267			ibase = r
268			scale = t
269			return (g/1)
270		}
271		e = g
272	}
273}
274/* vim: set filetype=bc shiftwidth=8 noexpandtab: */