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