Annotation of src/usr.bin/bc/bc.library, Revision 1.3
1.3 ! otto 1: /* $OpenBSD: bc.library,v 1.2 2007/01/29 11:02:53 otto Exp $ */
1.1 otto 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: * @(#)bc.library 5.1 (Berkeley) 4/17/91
39: */
40:
41: scale = 20
1.2 otto 42: define e(x) {
1.3 ! otto 43: auto a, b, c, d, e, g, t, w, y, r
1.1 otto 44:
1.3 ! otto 45: r = ibase
! 46: ibase = A
1.1 otto 47: t = scale
48: scale = t + .434*x + 1
49:
50: w = 0
1.2 otto 51: if (x < 0) {
1.1 otto 52: x = -x
53: w = 1
54: }
55: y = 0
1.2 otto 56: while (x > 2) {
1.1 otto 57: x = x/2
58: y = y + 1
59: }
60:
1.2 otto 61: a = 1
62: b = 1
63: c = b
64: d = 1
65: e = 1
66: for (a = 1; 1 == 1; a++) {
67: b = b*x
68: c = c*a + b
69: d = d*a
1.1 otto 70: g = c/d
1.2 otto 71: if (g == e) {
1.1 otto 72: g = g/1
1.2 otto 73: while (y--) {
1.1 otto 74: g = g*g
75: }
76: scale = t
1.3 ! otto 77: ibase = r
! 78: if (w == 1) return (1/g)
! 79: return (g/1)
1.1 otto 80: }
1.2 otto 81: e = g
1.1 otto 82: }
83: }
84:
1.2 otto 85: define l(x) {
1.3 ! otto 86: auto a, b, c, d, e, f, g, u, s, t, r
! 87: r = ibase
! 88: ibase = A
! 89: if (x <= 0) {
! 90: a = (1 - 10^scale)
! 91: ibase = r
! 92: return (a)
! 93: }
1.1 otto 94: t = scale
95:
1.2 otto 96: f = 1
1.1 otto 97: scale = scale + scale(x) - length(x) + 1
1.2 otto 98: s = scale
99: while (x > 2) {
100: s = s + (length(x) - scale(x))/2 + 1
101: if (s > 0) scale = s
1.1 otto 102: x = sqrt(x)
1.2 otto 103: f = f*2
1.1 otto 104: }
1.2 otto 105: while (x < .5) {
106: s = s + (length(x) - scale(x))/2 + 1
107: if (s > 0) scale = s
1.1 otto 108: x = sqrt(x)
1.2 otto 109: f = f*2
1.1 otto 110: }
111:
112: scale = t + length(f) - scale(f) + 1
1.2 otto 113: u = (x - 1)/(x + 1)
1.1 otto 114:
115: scale = scale + 1.1*length(t) - 1.1*scale(t)
116: s = u*u
117: b = 2*f
118: c = b
119: d = 1
120: e = 1
1.2 otto 121: for (a = 3; 1 == 1 ; a = a + 2) {
122: b = b*s
123: c = c*a + d*b
124: d = d*a
125: g = c/d
126: if (g == e) {
1.1 otto 127: scale = t
1.3 ! otto 128: ibase = r
1.2 otto 129: return (u*c/d)
1.1 otto 130: }
1.2 otto 131: e = g
1.1 otto 132: }
133: }
134:
1.2 otto 135: define s(x) {
1.3 ! otto 136: auto a, b, c, s, t, y, p, n, i, r
! 137: r = ibase
! 138: ibase = A
1.1 otto 139: t = scale
140: y = x/.7853
141: s = t + length(y) - scale(y)
1.2 otto 142: if (s < t) s = t
1.1 otto 143: scale = s
144: p = a(1)
145:
146: scale = 0
1.2 otto 147: if (x >= 0) n = (x/(2*p) + 1)/2
148: if (x < 0) n = (x/(2*p) - 1)/2
1.1 otto 149: x = x - 4*n*p
1.2 otto 150: if (n % 2 != 0) x = -x
1.1 otto 151:
152: scale = t + length(1.2*t) - scale(1.2*t)
153: y = -x*x
154: a = x
155: b = 1
156: s = x
1.2 otto 157: for (i =3 ; 1 == 1; i = i + 2) {
1.1 otto 158: a = a*y
1.2 otto 159: b = b*i*(i - 1)
1.1 otto 160: c = a/b
1.2 otto 161: if (c == 0) {
162: scale = t
1.3 ! otto 163: ibase = r
1.2 otto 164: return (s/1)
165: }
166: s = s + c
1.1 otto 167: }
168: }
169:
1.2 otto 170: define c(x) {
1.3 ! otto 171: auto t, r
! 172: r = ibase
! 173: ibase = A
1.1 otto 174: t = scale
1.2 otto 175: scale = scale + 1
176: x = s(x + 2*a(1))
1.1 otto 177: scale = t
1.3 ! otto 178: ibase = r
1.2 otto 179: return (x/1)
1.1 otto 180: }
181:
1.2 otto 182: define a(x) {
1.3 ! otto 183: auto a, b, c, d, e, f, g, s, t, r
1.2 otto 184: if (x == 0) return(0)
1.3 ! otto 185:
! 186: r = ibase
! 187: ibase = A
1.2 otto 188: if (x == 1) {
189: if (scale < 52) {
1.3 ! otto 190: a = .7853981633974483096156608458198757210492923498437764/1
! 191: ibase = r
! 192: return (a)
1.1 otto 193: }
194: }
195: t = scale
1.2 otto 196: f = 1
197: while (x > .5) {
1.1 otto 198: scale = scale + 1
1.2 otto 199: x = -(1 - sqrt(1. + x*x))/x
200: f = f*2
1.1 otto 201: }
1.2 otto 202: while (x < -.5) {
1.1 otto 203: scale = scale + 1
1.2 otto 204: x = -(1 - sqrt(1. + x*x))/x
205: f = f*2
1.1 otto 206: }
207: s = -x*x
208: b = f
209: c = f
210: d = 1
211: e = 1
1.2 otto 212: for (a = 3; 1 == 1; a = a + 2) {
213: b = b*s
214: c = c*a + d*b
215: d = d*a
216: g = c/d
217: if (g == e) {
1.3 ! otto 218: ibase = r
1.1 otto 219: scale = t
1.2 otto 220: return (x*c/d)
1.1 otto 221: }
1.2 otto 222: e = g
1.1 otto 223: }
224: }
225:
1.2 otto 226: define j(n,x) {
1.3 ! otto 227: auto a, b, c, d, e, g, i, s, k, t, r
1.1 otto 228:
1.3 ! otto 229: r = ibase
! 230: ibase = A
1.1 otto 231: t = scale
232: k = 1.36*x + 1.16*t - n
233: k = length(k) - scale(k)
1.2 otto 234: if (k > 0) scale = scale + k
1.1 otto 235:
1.2 otto 236: s = -x*x/4
237: if (n < 0) {
238: n = -n
239: x = -x
240: }
241: a = 1
242: c = 1
243: for (i = 1; i <= n; i++) {
244: a = a*x
1.1 otto 245: c = c*2*i
246: }
1.2 otto 247: b = a
248: d = 1
249: e = 1
250: for (i = 1; 1; i++) {
251: a = a*s
252: b = b*i*(n + i) + a
253: c = c*i*(n + i)
254: g = b/c
255: if (g == e) {
1.3 ! otto 256: ibase = r
1.1 otto 257: scale = t
1.2 otto 258: return (g/1)
1.1 otto 259: }
1.2 otto 260: e = g
1.1 otto 261: }
262: }