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