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