[796] | 1 | ;"16-Feb-1999, 16:54:35
|
---|
| 2 | ;"Routine Save for all M[UMPS] Library Functions
|
---|
| 3 | ;
|
---|
| 4 | ;" Unless otherwise noted, the code below
|
---|
| 5 | ;" was approved in document X11/95-11
|
---|
| 6 | ;
|
---|
| 7 | ;" If corrections have been applied,
|
---|
| 8 | ;" first the original line appears,
|
---|
| 9 | ;" with three semicolons at the beginning of the line.
|
---|
| 10 | ;
|
---|
| 11 | ;" Then the source of the correction is acknowledged,
|
---|
| 12 | ;" then the corrected line appears, followed by a
|
---|
| 13 | ;" line containing three semicolons.
|
---|
| 14 | ;
|
---|
| 15 | ;"Downloaded from http://www.jacquardsystems.com/Examples/lib/mlibfunc.rs
|
---|
| 16 | ;"on 5/21/07
|
---|
| 17 | ABS(X) Quit $Translate(+X,"-")
|
---|
| 18 | ;===
|
---|
| 19 | ;
|
---|
| 20 | ;
|
---|
| 21 | ARCCOS(X) ;
|
---|
| 22 | ;;;" ;" Number ~~
|
---|
| 23 | ;" Winfried Gerum (8 June 1995)
|
---|
| 24 | ;" Comment: This version of the function is
|
---|
| 25 | ;" optimized for speed, not for precision.
|
---|
| 26 | ;" The 'precision' parameter is not supported,
|
---|
| 27 | ;" and the precision is at best 2 in 10**-8.
|
---|
| 28 | ;;;
|
---|
| 29 | ;
|
---|
| 30 | New A,N,R,SIGN,XX
|
---|
| 31 | If X<-1 Set $Ecode=",M28,"
|
---|
| 32 | If X>1 Set $Ecode=",M28,"
|
---|
| 33 | Set SIGN=1 Set:X<0 X=-X,SIGN=-1
|
---|
| 34 | Set A(0)=1.5707963050,A(1)=-0.2145988016,A(2)=0.0889789874
|
---|
| 35 | Set A(3)=-0.0501743046,A(4)=0.0308918810,A(5)=-0.0170881256
|
---|
| 36 | Set A(6)=0.0066700901,A(7)=-0.0012624911
|
---|
| 37 | Set R=A(0),XX=1 For N=1:1:7 Set XX=XX*X,R=A(N)*XX+R
|
---|
| 38 | ;
|
---|
| 39 | ;;;" Set R=$%SQRT^MATH(1-X)*R ;" Number ~~
|
---|
| 40 | ;" Winfried Gerum (8 June 1995)
|
---|
| 41 | Set R=$%SQRT^MATH(1-X,11)*R
|
---|
| 42 | ;;;
|
---|
| 43 | ;
|
---|
| 44 | Quit R*SIGN
|
---|
| 45 | ;===
|
---|
| 46 | ;
|
---|
| 47 | ;
|
---|
| 48 | ARCCOS(X,PREC) ;
|
---|
| 49 | ;
|
---|
| 50 | ;;;" New L,LIM,K,SIG,SIGS ;" Number ~~
|
---|
| 51 | ;" Winfried Gerum (8 June 1995)
|
---|
| 52 | New L,LIM,K,SIG,SIGS,VALUE
|
---|
| 53 | ;;;
|
---|
| 54 | ;
|
---|
| 55 | If X<-1 Set $Ecode=",M28,"
|
---|
| 56 | If X>1 Set $Ecode=",M28,"
|
---|
| 57 | Set PREC=$Get(PREC,11)
|
---|
| 58 | ;
|
---|
| 59 | ;;;" If $Translate(X,"-")=1 Set VALUE=0 Quit ;" Number ~~
|
---|
| 60 | ;" Winfried Gerum (8 June 1995)
|
---|
| 61 | ;" Eli Reidler (28 June 1996)
|
---|
| 62 | If $Translate(X,"-")=1 Quit 0
|
---|
| 63 | ;;;
|
---|
| 64 | ;
|
---|
| 65 | Set SIG=$Select(X<0:-1,1:1),VALUE=1-(X*X)
|
---|
| 66 | ;
|
---|
| 67 | ;;;" Set X=$%SQRT^MATH(VALUE) ;" Number ~~
|
---|
| 68 | ;" Winfried Gerum (8 June 1995)
|
---|
| 69 | Set X=$%SQRT^MATH(VALUE,PREC)
|
---|
| 70 | ;;;
|
---|
| 71 | ;
|
---|
| 72 | ;;;" If $Translate(X,"-")=1 Do Quit ;" Number ~~
|
---|
| 73 | ;" Winfried Gerum (8 June 1995)
|
---|
| 74 | ;" Eli Reidler (28 June 1996)
|
---|
| 75 | If $Translate(X,"-")=1 Do Quit VALUE
|
---|
| 76 | . ;;;
|
---|
| 77 | . ;
|
---|
| 78 | . Set VALUE=$%PI^MATH()/2*X
|
---|
| 79 | . Quit
|
---|
| 80 | ;
|
---|
| 81 | ;;;" If X>0.9 Do Quit ;" Number ~~
|
---|
| 82 | ;" Winfried Gerum (8 June 1995)
|
---|
| 83 | ;" Eli Reidler (28 June 1996)
|
---|
| 84 | If X>0.9 Do Quit VALUE
|
---|
| 85 | . ;;;
|
---|
| 86 | . ;
|
---|
| 87 | . Set SIGS=$Select(X<0:-1,1:1)
|
---|
| 88 | . Set VALUE=1/(1/X/X-1)
|
---|
| 89 | . ;
|
---|
| 90 | . ;;;" Set X=$%SQRT^MATH(VALUE) ;" Number ~~
|
---|
| 91 | . ;" Winfried Gerum (8 June 1995)
|
---|
| 92 | . Set X=$%SQRT^MATH(VALUE,PREC)
|
---|
| 93 | . ;;;
|
---|
| 94 | . ;
|
---|
| 95 | . ;
|
---|
| 96 | . ;;;" Set VALUE=$%ARCTAN^MATH(X,10)*SIGS ;" Number ~~
|
---|
| 97 | . ;" Winfried Gerum (8 June 1995)
|
---|
| 98 | . Set VALUE=$%ARCTAN^MATH(X,PREC)*SIGS
|
---|
| 99 | . ;;;
|
---|
| 100 | ;
|
---|
| 101 | . Quit
|
---|
| 102 | Set (VALUE,L)=X
|
---|
| 103 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 104 | For K=3:2 Do Quit:($Translate(L,"-")<LIM)
|
---|
| 105 | . Set L=L*X*X*(K-2)/(K-1)*(K-2)/K,VALUE=VALUE+L
|
---|
| 106 | . Quit
|
---|
| 107 | Quit $Select(SIG<0:$%PI^MATH()-VALUE,1:VALUE)
|
---|
| 108 | ;===
|
---|
| 109 | ;
|
---|
| 110 | ;
|
---|
| 111 | ARCCOSH(X,PREC) ;
|
---|
| 112 | If X<1 Set $Ecode=",M28,"
|
---|
| 113 | New SQ
|
---|
| 114 | ;
|
---|
| 115 | ;;;" ;" Number ~~
|
---|
| 116 | ;" Winfried Gerum (8 June 1995)
|
---|
| 117 | ;" Alan Frank (October 1995)
|
---|
| 118 | Set PREC=$Get(PREC,11)
|
---|
| 119 | ;;;
|
---|
| 120 | ;
|
---|
| 121 | Set SQ=$%SQRT^MATH(X*X-1,PREC)
|
---|
| 122 | Quit $%LOG^MATH(X+SQ,PREC)
|
---|
| 123 | ;===
|
---|
| 124 | ;
|
---|
| 125 | ;
|
---|
| 126 | ARCCOT(X,PREC) ;
|
---|
| 127 | Set PREC=$Get(PREC,11)
|
---|
| 128 | Set X=1/X
|
---|
| 129 | Quit $%ARCTAN^MATH(X,PREC)
|
---|
| 130 | ;===
|
---|
| 131 | ;
|
---|
| 132 | ;
|
---|
| 133 | ARCCOTH(X,PREC) ;
|
---|
| 134 | New L1,L2
|
---|
| 135 | ;
|
---|
| 136 | ;;;" ;" Number ~~
|
---|
| 137 | ;" Winfried Gerum (8 June 1995)
|
---|
| 138 | ;" Alan Frank (October 1995)
|
---|
| 139 | Set PREC=$Get(PREC,11)
|
---|
| 140 | ;;;
|
---|
| 141 | ;
|
---|
| 142 | Set L1=$%LOG^MATH(X+1,PREC)
|
---|
| 143 | Set L2=$%LOG^MATH(X-1,PREC)
|
---|
| 144 | Quit L1-L2/2
|
---|
| 145 | ;===
|
---|
| 146 | ;
|
---|
| 147 | ;
|
---|
| 148 | ARCCSC(X,PREC) ;
|
---|
| 149 | Set PREC=$Get(PREC,11)
|
---|
| 150 | Set X=1/X
|
---|
| 151 | Quit $%ARCSIN^MATH(X,PREC)
|
---|
| 152 | ;===
|
---|
| 153 | ;
|
---|
| 154 | ;
|
---|
| 155 | ARCSEC(X,PREC) ;
|
---|
| 156 | Set PREC=$Get(PREC,11)
|
---|
| 157 | Set X=1/X
|
---|
| 158 | Quit $%ARCCOS^MATH(X,PREC)
|
---|
| 159 | ;===
|
---|
| 160 | ;
|
---|
| 161 | ;
|
---|
| 162 | ARCSIN(X) ;
|
---|
| 163 | ;;;" ;" Number ~~
|
---|
| 164 | ;" Winfried Gerum (8 June 1995)
|
---|
| 165 | ;" Comment: This version of the function is
|
---|
| 166 | ;" optimized for speed, not for precision.
|
---|
| 167 | ;" The 'precision' parameter is not supported,
|
---|
| 168 | ;" and the precision is at best 2 in 10**-8.
|
---|
| 169 | ;;;
|
---|
| 170 | ;
|
---|
| 171 | New A,N,R,SIGN,XX
|
---|
| 172 | If X<-1 Set $Ecode=",M28,"
|
---|
| 173 | If X>1 Set $Ecode=",M28,"
|
---|
| 174 | Set SIGN=1 Set:X<0 X=-X,SIGN=-1
|
---|
| 175 | Set A(0)=1.5707963050,A(1)=-0.2145988016,A(2)=0.0889789874
|
---|
| 176 | Set A(3)=-0.0501743046,A(4)=0.0308918810,A(5)=-0.0170881256
|
---|
| 177 | Set A(6)=0.0066700901,A(7)=-0.0012624911
|
---|
| 178 | Set R=A(0),XX=1 For N=1:1:7 Set XX=XX*X,R=A(N)*XX+R
|
---|
| 179 | ;
|
---|
| 180 | ;;;" Set R=$%SQRT^MATH(1-X)*R ;" Number ~~
|
---|
| 181 | ;" Winfried Gerum (8 June 1995)
|
---|
| 182 | Set R=$%SQRT^MATH(1-X,11)*R
|
---|
| 183 | ;;;
|
---|
| 184 | ;
|
---|
| 185 | Set R=$%PI^MATH()/2-R
|
---|
| 186 | Quit R*SIGN
|
---|
| 187 | ;===
|
---|
| 188 | ;
|
---|
| 189 | ;
|
---|
| 190 | ARCSIN(X,PREC) ;
|
---|
| 191 | New L,LIM,K,SIGS,VALUE
|
---|
| 192 | Set PREC=$Get(PREC,11)
|
---|
| 193 | ;
|
---|
| 194 | ;;;" If $Translate(X,"-")=1 Do Quit ;" Number ~~
|
---|
| 195 | ;" Winfried Gerum (8 June 1995)
|
---|
| 196 | ;" Eli Reidler (28 June 1996)
|
---|
| 197 | If $Translate(X,"-")=1 Do Quit VALUE
|
---|
| 198 | . ;;;
|
---|
| 199 | . ;
|
---|
| 200 | . Set VALUE=$%PI^MATH()/2*X
|
---|
| 201 | . Quit
|
---|
| 202 | ;
|
---|
| 203 | ;;;" If X>0.99999 Do Quit ;" Number ~~
|
---|
| 204 | ;" Winfried Gerum (8 June 1995)
|
---|
| 205 | ;" Eli Reidler (28 June 1996)
|
---|
| 206 | If X>0.99999 Do Quit VALUE
|
---|
| 207 | . ;;;
|
---|
| 208 | . ;
|
---|
| 209 | . Set SIGS=$Select(X<0:-1,1:1)
|
---|
| 210 | . Set VALUE=1/(1/X/X-1)
|
---|
| 211 | . ;
|
---|
| 212 | . ;;;" Set X=$%SQRT^MATH(VALUE) ;" Number ~~
|
---|
| 213 | . ;" Winfried Gerum (8 June 1995)
|
---|
| 214 | . Set X=$%SQRT^MATH(VALUE,PREC)
|
---|
| 215 | . ;;;
|
---|
| 216 | . ;
|
---|
| 217 | . ;;;" Set VALUE=$%ARCTAN^MATH(X,10)*SIGS ;" Number ~~
|
---|
| 218 | . ;" Winfried Gerum (8 June 1995)
|
---|
| 219 | . Set VALUE=$%ARCTAN^MATH(X,PREC)*SIGS
|
---|
| 220 | . ;;;
|
---|
| 221 | . ;
|
---|
| 222 | . Quit
|
---|
| 223 | Set (VALUE,L)=X
|
---|
| 224 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 225 | For K=3:2 Do Quit:($Translate(L,"-")<LIM)
|
---|
| 226 | . Set L=L*X*X*(K-2)/(K-1)*(K-2)/K,VALUE=VALUE+L
|
---|
| 227 | . Quit
|
---|
| 228 | Quit VALUE
|
---|
| 229 | ;===
|
---|
| 230 | ;
|
---|
| 231 | ;
|
---|
| 232 | ARCSINH(X,PREC) ;
|
---|
| 233 | If X<1 Set $Ecode=",M28,"
|
---|
| 234 | New SQ
|
---|
| 235 | ;
|
---|
| 236 | ;;;" ;" Number ~~
|
---|
| 237 | ;" Winfried Gerum (8 June 1995)
|
---|
| 238 | ;" Alan Frank (October 1995)
|
---|
| 239 | Set PREC=$Get(PREC,11)
|
---|
| 240 | ;;;
|
---|
| 241 | ;
|
---|
| 242 | Set SQ=$%SQRT^MATH(X*X+1,PREC)
|
---|
| 243 | Quit $%LOG^MATH(X+SQ,PREC)
|
---|
| 244 | ;===
|
---|
| 245 | ;
|
---|
| 246 | ;
|
---|
| 247 | ARCTAN(X,PREC) ;
|
---|
| 248 | New FOLD,HI,L,LIM,LO,K,SIGN,SIGS,SIGT,VALUE
|
---|
| 249 | Set PREC=$Get(PREC,11)
|
---|
| 250 | Set LO=0.0000000001,HI=9999999999
|
---|
| 251 | Set SIGT=$Select(X<0:-1,1:1),X=$Translate(X,"-")
|
---|
| 252 | Set X=$Select(X<LO:LO,X>HI:HI,1:X)
|
---|
| 253 | ;
|
---|
| 254 | ;;;" Set FOLD=$Select(X'<1:0,1:1), ;" Number ~~
|
---|
| 255 | ;" Eli Reidler (28 June 1996)
|
---|
| 256 | Set FOLD=$Select(X'<1:0,1:1)
|
---|
| 257 | ;;;
|
---|
| 258 | ;
|
---|
| 259 | Set X=$Select(FOLD:1/X,1:X)
|
---|
| 260 | Set L=X,VALUE=$%PI^MATH()/2-(1/X),SIGN=1
|
---|
| 261 | ;
|
---|
| 262 | ;;;" If X<1.3 Do Quit ;" Number ~~
|
---|
| 263 | ;" Winfried Gerum (8 June 1995)
|
---|
| 264 | ;" Eli Reidler (28 June 1996)
|
---|
| 265 | If X<1.3 Do Quit VALUE
|
---|
| 266 | . ;;;
|
---|
| 267 | . ;
|
---|
| 268 | . Set X=$Select(FOLD:1/X,1:X),VALUE=1/((1/X/X)+1)
|
---|
| 269 | . ;
|
---|
| 270 | . ;;;" Set $%SQRT^MATH(VALUE) ;" Number ~~
|
---|
| 271 | . ;" Winfried Gerum (8 June 1995)
|
---|
| 272 | . ;" Eli Reidler (28 June 1996)
|
---|
| 273 | . Set X=$%SQRT^MATH(VALUE,PREC)
|
---|
| 274 | . ;;;
|
---|
| 275 | . ;
|
---|
| 276 | . If $Translate(X,"-")=1 Do Quit
|
---|
| 277 | . . Set VALUE=$%PI^MATH()/2*X
|
---|
| 278 | . . Quit
|
---|
| 279 | . If X>0.9 Do Quit
|
---|
| 280 | . . Set SIGS=$Select(X<0:-1,1:1)
|
---|
| 281 | . . Set VALUE=1/(1/X/X-1)
|
---|
| 282 | . . Set X=$%SQRT^MATH(VALUE)
|
---|
| 283 | . . Set VALUE=$$ARCTAN(X,10)
|
---|
| 284 | . . Set VALUE=VALUE*SIGS
|
---|
| 285 | . . Quit
|
---|
| 286 | . Set (VALUE,L)=X
|
---|
| 287 | . Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 288 | . For K=3:2 Do Quit:($Translate(L,"-")<LIM)
|
---|
| 289 | . . Set L=L*X*X*(K-2)/(K-1)*(K-2)/K,VALUE=VALUE+L
|
---|
| 290 | . . Quit
|
---|
| 291 | . Set VALUE=$Select(SIGT<1:-VALUE,1:VALUE)
|
---|
| 292 | . Quit
|
---|
| 293 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 294 | For K=3:2 Do Quit:($Translate(1/L,"-")<LIM)
|
---|
| 295 | . ;
|
---|
| 296 | . ;;;" Set L=L*X*X,VALUE=VALUE+(1/(K*L)*SIGN), ;" Number ~~
|
---|
| 297 | . ;" Eli Reidler (28 June 1996)
|
---|
| 298 | . Set L=L*X*X,VALUE=VALUE+(1/(K*L)*SIGN)
|
---|
| 299 | . ;;;
|
---|
| 300 | . ;
|
---|
| 301 | . Set SIGN=SIGN*-1
|
---|
| 302 | . Quit
|
---|
| 303 | Set VALUE=$Select(FOLD:$%PI^MATH()/2-VALUE,1:VALUE)
|
---|
| 304 | Set VALUE=$Select(SIGT<1:-VALUE,1:VALUE)
|
---|
| 305 | Quit VALUE
|
---|
| 306 | ;===
|
---|
| 307 | ;
|
---|
| 308 | ;
|
---|
| 309 | ARCTANH(X,PREC) ;
|
---|
| 310 | If X<-1 Set $Ecode=",M28,"
|
---|
| 311 | If X>1 Set $Ecode=",M28,"
|
---|
| 312 | ;
|
---|
| 313 | ;;;" ;" Number ~~
|
---|
| 314 | ;" Winfried Gerum (8 June 1995)
|
---|
| 315 | ;" Alan Frank (October 1995)
|
---|
| 316 | Set PREC=$Get(PREC,11)
|
---|
| 317 | ;;;
|
---|
| 318 | ;
|
---|
| 319 | Quit $%LOG^MATH(1+X/(1-X),PREC)/2
|
---|
| 320 | ;===
|
---|
| 321 | ;
|
---|
| 322 | ;
|
---|
| 323 | CABS(Z) ;
|
---|
| 324 | New ZRE,ZIM
|
---|
| 325 | Set ZRE=+Z,ZIM=+$Piece(Z,"%",2)
|
---|
| 326 | Quit $%SQRT^MATH(ZRE*ZRE+(ZIM*ZIM))
|
---|
| 327 | ;===
|
---|
| 328 | ;
|
---|
| 329 | ;
|
---|
| 330 | CADD(X,Y) ;
|
---|
| 331 | New XRE,XIM,YRE,YIM
|
---|
| 332 | Set XRE=+X,XIM=+$Piece(X,"%",2)
|
---|
| 333 | Set YRE=+Y,YIM=+$Piece(Y,"%",2)
|
---|
| 334 | Quit XRE+YRE_"%"_(XIM+YIM)
|
---|
| 335 | ;===
|
---|
| 336 | ;
|
---|
| 337 | ;
|
---|
| 338 | CCOS(Z,PREC) ;
|
---|
| 339 | New E1,E2,IA
|
---|
| 340 | ;
|
---|
| 341 | ;;;" ;" Number ~~
|
---|
| 342 | ;" Alan Frank (October 1995)
|
---|
| 343 | Set PREC=$Get(PREC,11)
|
---|
| 344 | ;;;
|
---|
| 345 | ;
|
---|
| 346 | Set IA=$%CMUL^MATH(Z,"0%1")
|
---|
| 347 | Set E1=$%CEXP^MATH(IA,PREC)
|
---|
| 348 | Set IA=-IA_"%"_(-$Piece(IA,"%",2))
|
---|
| 349 | Set E2=$%CEXP^MATH(IA,PREC)
|
---|
| 350 | Set IA=$%CADD^MATH(E1,E2)
|
---|
| 351 | Quit $%CMUL^MATH(IA,"0.5%0")
|
---|
| 352 | ;===
|
---|
| 353 | ;
|
---|
| 354 | ;
|
---|
| 355 | CDIV(X,Y) ;
|
---|
| 356 | New D,IM,RE,XIM,XRE,YIM,YRE
|
---|
| 357 | Set XRE=+X,XIM=+$Piece(X,"%",2)
|
---|
| 358 | Set YRE=+Y,YIM=+$Piece(Y,"%",2)
|
---|
| 359 | Set D=YRE*YRE+(YIM*YIM)
|
---|
| 360 | Set RE=XRE*YRE+(XIM*YIM)/D
|
---|
| 361 | Set IM=XIM*YRE-(XRE*YIM)/D
|
---|
| 362 | Quit RE_"%"_IM
|
---|
| 363 | ;===
|
---|
| 364 | ;
|
---|
| 365 | ;
|
---|
| 366 | CEXP(Z,PREC) ;
|
---|
| 367 | New R,ZIM,ZRE
|
---|
| 368 | ;
|
---|
| 369 | ;;;" ;" Number ~~
|
---|
| 370 | ;" Alan Frank (October 1995)
|
---|
| 371 | Set PREC=$Get(PREC,11)
|
---|
| 372 | ;;;
|
---|
| 373 | ;
|
---|
| 374 | Set ZRE=+Z,ZIM=+$Piece(Z,"%",2)
|
---|
| 375 | Set R=$%EXP^MATH(ZRE,PREC)
|
---|
| 376 | Quit R*$%COS^MATH(ZIM,PREC)_"%"_(R*$%SIN^MATH(ZIM,PREC))
|
---|
| 377 | ;===
|
---|
| 378 | ;
|
---|
| 379 | ;
|
---|
| 380 | CLOG(Z,PREC) ;
|
---|
| 381 | New ABS,ARG,ZIM,ZRE
|
---|
| 382 | ;
|
---|
| 383 | ;;;" ;" Number ~~
|
---|
| 384 | ;" Alan Frank (October 1995)
|
---|
| 385 | Set PREC=$Get(PREC,11)
|
---|
| 386 | ;;;
|
---|
| 387 | ;
|
---|
| 388 | Set ABS=$%CABS^MATH(Z)
|
---|
| 389 | Set ZRE=+Z,ZIM=+$Piece(Z,"%",2)
|
---|
| 390 | ;
|
---|
| 391 | ;;;" Set ARG=$%ARCTAN^MATH(ZIM,ZRE,PREC) ;" Number ~~
|
---|
| 392 | ;" Alan Frank (October 1995)
|
---|
| 393 | Set ARG=$%ARCTAN^MATH(ZIM/ZRE,PREC)
|
---|
| 394 | ;;;
|
---|
| 395 | ;
|
---|
| 396 | Quit $%LOG^MATH(ABS,PREC)_"%"_ARG
|
---|
| 397 | ;===
|
---|
| 398 | ;
|
---|
| 399 | ;
|
---|
| 400 | CMUL(X,Y) ;
|
---|
| 401 | New XIM,XRE,YIM,YRE
|
---|
| 402 | Set XRE=+X,XIM=+$Piece(X,"%",2)
|
---|
| 403 | Set YRE=+Y,YIM=+$Piece(Y,"%",2)
|
---|
| 404 | Quit XRE*YRE-(XIM*YIM)_"%"_(XRE*YIM+(XIM*YRE))
|
---|
| 405 | ;===
|
---|
| 406 | ;
|
---|
| 407 | ;
|
---|
| 408 | COMPLEX(X) Quit +X_"%0"
|
---|
| 409 | ;===
|
---|
| 410 | ;
|
---|
| 411 | ;
|
---|
| 412 | CONJUG(Z) ;
|
---|
| 413 | New ZIM,ZRE
|
---|
| 414 | Set ZRE=+Z,ZIM=+$Piece(Z,"%",2)
|
---|
| 415 | Quit ZRE_"%"_(-ZIM)
|
---|
| 416 | ;===
|
---|
| 417 | ;
|
---|
| 418 | ;
|
---|
| 419 | COS(X,PREC) ;
|
---|
| 420 | New L,LIM,K,SIGN,VALUE
|
---|
| 421 | ;
|
---|
| 422 | ;;;" Set:X[":" X=$%DMSDEC^MATH(X,12) ;" Number ~~
|
---|
| 423 | ;" Winfried Gerum (8 June 1995)
|
---|
| 424 | ;" Comment: The official description does not mention than
|
---|
| 425 | ;" the function may also be called with the first
|
---|
| 426 | ;" parameter in degrees, minutes and seconds.
|
---|
| 427 | Set:X[":" X=$%DMSDEC^MATH(X)
|
---|
| 428 | ;;;
|
---|
| 429 | ;
|
---|
| 430 | Set PREC=$Get(PREC,11)
|
---|
| 431 | Set X=X#(2*$%PI^MATH())
|
---|
| 432 | Set (VALUE,L)=1,SIGN=-1
|
---|
| 433 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 434 | For K=2:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 435 | . Set L=L*X*X/(K-1*K),VALUE=VALUE+(SIGN*L)
|
---|
| 436 | . Quit
|
---|
| 437 | Quit VALUE
|
---|
| 438 | ;===
|
---|
| 439 | ;
|
---|
| 440 | ;
|
---|
| 441 | COS(X) ;
|
---|
| 442 | ;;;" ;" Number ~~
|
---|
| 443 | ;" Winfried Gerum (8 June 1995)
|
---|
| 444 | ;" Comment: This version of the function is
|
---|
| 445 | ;" optimized for speed, not for precision.
|
---|
| 446 | ;" The 'precision' parameter is not supported,
|
---|
| 447 | ;" and the precision is at best 1 in 10**-9.
|
---|
| 448 | ;" Note that this function does not accept its
|
---|
| 449 | ;" parameter in degrees, minutes and seconds.
|
---|
| 450 | ;;;
|
---|
| 451 | ;
|
---|
| 452 | New A,N,PI,R,SIGN,XX
|
---|
| 453 | ;
|
---|
| 454 | ;" This approximation only works for 0 <= x <= pi/2
|
---|
| 455 | ;" so reduce angle to correct quadrant.
|
---|
| 456 | ;
|
---|
| 457 | Set PI=$%PI^MATH(),X=X#(PI*2),SIGN=1
|
---|
| 458 | Set:X>PI X=2*PI-X
|
---|
| 459 | Set:X*2>PI X=PI-X,SIGN=-1
|
---|
| 460 | ;
|
---|
| 461 | Set XX=X*X,A(1)=-0.4999999963,A(2)=0.0416666418
|
---|
| 462 | Set A(3)=-0.0013888397,A(4)=0.0000247609,A(5)=-0.0000002605
|
---|
| 463 | Set (X,R)=1 For N=1:1:5 Set X=X*XX,R=A(N)*X+R
|
---|
| 464 | Quit R*SIGN
|
---|
| 465 | ;===
|
---|
| 466 | ;
|
---|
| 467 | ;
|
---|
| 468 | COSH(X,PREC) ;
|
---|
| 469 | ;
|
---|
| 470 | ;;;" New F,I,P,R,T,XX ;" Number ~~
|
---|
| 471 | ;" Winfried Gerum (8 June 1995)
|
---|
| 472 | New E,F,I,P,R,T,XX
|
---|
| 473 | ;;;
|
---|
| 474 | ;
|
---|
| 475 | Set PREC=$Get(PREC,11)+1
|
---|
| 476 | Set @("E=1E-"_PREC)
|
---|
| 477 | Set XX=X*X,F=1,(P,R,T)=1,I=1
|
---|
| 478 | For Set T=T*XX,F=I+1*I*F,R=T/F+R,P=P-R/R,I=I+2 If -E<P,P<E Quit
|
---|
| 479 | Quit R
|
---|
| 480 | ;===
|
---|
| 481 | ;
|
---|
| 482 | ;
|
---|
| 483 | COT(X,PREC) ;
|
---|
| 484 | New C,L,LIM,K,SIGN,VALUE
|
---|
| 485 | ;
|
---|
| 486 | ;;;" Set:X[":" X=$%DMSDEC^MATH(X,12) ;" Number ~~
|
---|
| 487 | ;" Winfried Gerum (8 June 1995)
|
---|
| 488 | ;" Comment: The official description does not mention than
|
---|
| 489 | ;" the function may also be called with the first
|
---|
| 490 | ;" parameter in degrees, minutes and seconds.
|
---|
| 491 | Set:X[":" X=$%DMSDEC^MATH(X)
|
---|
| 492 | ;;;
|
---|
| 493 | ;
|
---|
| 494 | Set PREC=$Get(PREC,11)
|
---|
| 495 | Set (VALUE,L)=1,SIGN=-1
|
---|
| 496 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 497 | For K=2:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 498 | . Set L=L*X*X/(K-1*K),VALUE=VALUE+(SIGN*L)
|
---|
| 499 | . Quit
|
---|
| 500 | Set C=VALUE
|
---|
| 501 | Set X=X#(2*$%PI^MATH())
|
---|
| 502 | Set (VALUE,L)=X,SIGN=-1
|
---|
| 503 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 504 | For K=3:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 505 | . Set L=L/(K-1)*X/K*X,VALUE=VALUE+(SIGN*L)
|
---|
| 506 | . Quit
|
---|
| 507 | If 'VALUE Quit "INFINITE"
|
---|
| 508 | Quit VALUE=C/VALUE
|
---|
| 509 | ;===
|
---|
| 510 | ;
|
---|
| 511 | ;
|
---|
| 512 | COTH(X,PREC) ;
|
---|
| 513 | New SINH
|
---|
| 514 | If 'X Quit "INFINITE"
|
---|
| 515 | ;
|
---|
| 516 | ;;;" ;" Number ~~
|
---|
| 517 | ;" Winfried Gerum (8 June 1995)
|
---|
| 518 | ;" Alan Frank (October 1995)
|
---|
| 519 | Set PREC=$Get(PREC,11)
|
---|
| 520 | ;;;
|
---|
| 521 | ;
|
---|
| 522 | Set SINH=$%SINH^MATH(X,PREC)
|
---|
| 523 | If 'SINH Quit "INFINITE"
|
---|
| 524 | Quit $%COSH^MATH(X,PREC)/SINH
|
---|
| 525 | ;===
|
---|
| 526 | ;
|
---|
| 527 | ;
|
---|
| 528 | CPOWER(Z,N,PREC) ;
|
---|
| 529 | New AR,NIM,NRE,PHI,PI,R,RHO,TH,ZIM,ZRE
|
---|
| 530 | ;
|
---|
| 531 | ;;;" ;" Number ~~
|
---|
| 532 | ;" Alan Frank (October 1995)
|
---|
| 533 | Set PREC=$Get(PREC,11)
|
---|
| 534 | ;;;
|
---|
| 535 | ;
|
---|
| 536 | Set ZRE=+Z,ZIM=+$Piece(Z,"%",2)
|
---|
| 537 | Set NRE=+N,NIM=+$Piece(N,"%",2)
|
---|
| 538 | If 'ZRE,'ZIM,'NRE,'NIM Set $Ecode=",M28,"
|
---|
| 539 | ;
|
---|
| 540 | ;;;" If 'ZRE,'ZIM Quit "0%0% ;" Number ~~
|
---|
| 541 | ;" Eli Reidler (28 June 1996)
|
---|
| 542 | If 'ZRE,'ZIM Quit "0%0"
|
---|
| 543 | ;;;
|
---|
| 544 | ;
|
---|
| 545 | Set PI=$%PI^MATH()
|
---|
| 546 | ;
|
---|
| 547 | ;;;" Set R=$%SQRT^MATH(ZRE*ZRE+(ZIM*ZIM,PREC)) ;" Number ~~
|
---|
| 548 | ;" Winfried Gerum (8 June 1995)
|
---|
| 549 | ;" Eli Reidler (28 June 1996)
|
---|
| 550 | Set R=$%SQRT^MATH(ZRE*ZRE+(ZIM*ZIM),PREC)
|
---|
| 551 | ;;;
|
---|
| 552 | ;
|
---|
| 553 | ;
|
---|
| 554 | ;;;" If ZRE Set TH=$%ARCTAN^MATH(ZIM,ZRE,PREC) ;" Number ~~
|
---|
| 555 | ;" Alan Frank (October 1995)
|
---|
| 556 | If ZRE Set TH=$%ARCTAN^MATH(ZIM/ZRE,PREC)
|
---|
| 557 | ;;;
|
---|
| 558 | ;
|
---|
| 559 | ;;;" Else Set TH=$SELECT(ZRE>0:PI/2,1:-PI/2) ;" Number ~~
|
---|
| 560 | ;" Winfried Gerum (8 June 1995)
|
---|
| 561 | Else Set TH=$SELECT(ZIM>0:PI/2,1:-PI/2)
|
---|
| 562 | ;;;
|
---|
| 563 | ;
|
---|
| 564 | Set RHO=$%LOG^MATH(R,PREC)
|
---|
| 565 | Set AR=$%EXP^MATH(RHO*NRE-(TH*NIM),PREC)
|
---|
| 566 | Set PHI=RHO*NIM+(NRE*TH)
|
---|
| 567 | Quit AR*$%COS^MATH(PHI,PREC)_"%"_(AR*$%SIN^MATH(PHI,PREC))
|
---|
| 568 | ;===
|
---|
| 569 | ;
|
---|
| 570 | ;
|
---|
| 571 | CSC(X,PREC) ;
|
---|
| 572 | New L,LIM,K,SIGN,VALUE
|
---|
| 573 | ;
|
---|
| 574 | ;;;" Set:X[":" X=$%DMSDEC^MATH(X,12) ;" Number ~~
|
---|
| 575 | ;" Winfried Gerum (8 June 1995)
|
---|
| 576 | ;" Comment: The official description does not mention than
|
---|
| 577 | ;" the function may also be called with the first
|
---|
| 578 | ;" parameter in degrees, minutes and seconds.
|
---|
| 579 | Set:X[":" X=$%DMSDEC^MATH(X)
|
---|
| 580 | ;;;
|
---|
| 581 | ;
|
---|
| 582 | ;;;" Set PREC=$Select($Data(PREC)#2:PREC,1:10) ;" Number ~~
|
---|
| 583 | ;" Winfried Gerum (8 June 1995)
|
---|
| 584 | Set PREC=$Get(PREC,11)
|
---|
| 585 | ;;;
|
---|
| 586 | ;
|
---|
| 587 | Set X=X#(2*$%PI^MATH())
|
---|
| 588 | Set (VALUE,L)=X,SIGN=-1
|
---|
| 589 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 590 | For K=3:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 591 | . Set L=L/(K-1)*X/K*X,VALUE=VALUE+(SIGN*L)
|
---|
| 592 | . Quit
|
---|
| 593 | If 'VALUE Quit "INFINITE"
|
---|
| 594 | Quit 1/VALUE
|
---|
| 595 | ;===
|
---|
| 596 | ;
|
---|
| 597 | ;
|
---|
| 598 | ;
|
---|
| 599 | CSCH(X,PREC) ;;;Quit 1/$%SINH^MATH(X,PREC) ;" Number ~~
|
---|
| 600 | ;" Winfried Gerum (8 June 1995)
|
---|
| 601 | ;" Alan Frank (October 1995)
|
---|
| 602 | Quit 1/$%SINH^MATH(X,$Get(PREC,11))
|
---|
| 603 | ;;;
|
---|
| 604 | ;
|
---|
| 605 | ;===
|
---|
| 606 | ;
|
---|
| 607 | ;
|
---|
| 608 | CSIN(Z,PREC) ;
|
---|
| 609 | New IA,E1,E2
|
---|
| 610 | ;
|
---|
| 611 | ;;;" ;" Number ~~
|
---|
| 612 | ;" Alan Frank (October 1995)
|
---|
| 613 | Set PREC=$Get(PREC,11)
|
---|
| 614 | ;;;
|
---|
| 615 | ;
|
---|
| 616 | Set IA=$%CMUL^MATH(Z,"0%1")
|
---|
| 617 | Set E1=$%CEXP^MATH(IA,PREC)
|
---|
| 618 | Set IA=-IA_"%"_(-$Piece(IA,"%",2))
|
---|
| 619 | Set E2=$%CEXP^MATH(IA,PREC)
|
---|
| 620 | Set IA=$%CSUB^MATH(E1,E2)
|
---|
| 621 | Set IA=$%CMUL^MATH(IA,"0.5%0")
|
---|
| 622 | Quit $%CMUL^MATH("0%-1",IA)
|
---|
| 623 | ;===
|
---|
| 624 | ;
|
---|
| 625 | ;
|
---|
| 626 | CSUB(X,Y) ;
|
---|
| 627 | New XIM,XRE,YIM,YRE
|
---|
| 628 | Set XRE=+X,XIM=+$Piece(X,"%",2)
|
---|
| 629 | Set YRE=+Y,YIM=+$Piece(Y,"%",2)
|
---|
| 630 | Quit XRE-YRE_"%"_(XIM-YIM)
|
---|
| 631 | ;===
|
---|
| 632 | ;
|
---|
| 633 | ;
|
---|
| 634 | DECDMS(X,PREC) ;
|
---|
| 635 | Set PREC=$Get(PREC,5)
|
---|
| 636 | Set X=X#360*3600
|
---|
| 637 | Set X=+$Justify(X,0,$Select((PREC-$Length(X\1))'<0:PREC-$Length(X\1),1:0))
|
---|
| 638 | Quit X\3600_":"_(X\60#60)_":"_(X#60)
|
---|
| 639 | ;===
|
---|
| 640 | ;
|
---|
| 641 | ;
|
---|
| 642 | DEGRAD(X) Quit X*3.14159265358979/180
|
---|
| 643 | ;===
|
---|
| 644 | ;
|
---|
| 645 | ;
|
---|
| 646 | DMSDEC(X) ;
|
---|
| 647 | Quit $Piece(X,":")+($Piece(X,":",2)/60)+($Piece(X,":",3)/3600)
|
---|
| 648 | ;===
|
---|
| 649 | ;
|
---|
| 650 | ;
|
---|
| 651 | E() Quit 2.71828182845905
|
---|
| 652 | ;===
|
---|
| 653 | ;
|
---|
| 654 | ;
|
---|
| 655 | EXP(X,PREC) ;
|
---|
| 656 | New L,LIM,K,VALUE
|
---|
| 657 | Set PREC=$Get(PREC,11)
|
---|
| 658 | Set L=X,VALUE=X+1
|
---|
| 659 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 660 | For K=2:1 Set L=L*X/K,VALUE=VALUE+L Quit:($Translate(L,"-")<LIM)
|
---|
| 661 | Quit VALUE
|
---|
| 662 | ;===
|
---|
| 663 | ;
|
---|
| 664 | ;
|
---|
| 665 | LOG(X,PREC) ;
|
---|
| 666 | New L,LIM,M,N,K,VALUE
|
---|
| 667 | If X'>0 Set $Ecode=",M28,"
|
---|
| 668 | Set PREC=$Get(PREC,11)
|
---|
| 669 | Set M=1
|
---|
| 670 | ;
|
---|
| 671 | ;;;" If X>0 For N=0:1 Quit:(X/M)<10 Set M=M*10 ;" Number ~~
|
---|
| 672 | ;" Winfried Gerum (8 June 1995)
|
---|
| 673 | For N=0:1 Quit:(X/M)<10 Set M=M*10
|
---|
| 674 | ;;;
|
---|
| 675 | ;
|
---|
| 676 | If X<1 For N=0:-1 Quit:(X/M)>0.1 Set M=M*0.1
|
---|
| 677 | Set X=X/M
|
---|
| 678 | Set X=(X-1)/(X+1),(VALUE,L)=X
|
---|
| 679 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 680 | For K=3:2 Set L=L*X*X,M=L/K,VALUE=M+VALUE Set:M<0 M=-M Quit:M<LIM
|
---|
| 681 | Set VALUE=VALUE*2+(N*2.30258509298749)
|
---|
| 682 | Quit VALUE
|
---|
| 683 | ;===
|
---|
| 684 | ;
|
---|
| 685 | ;
|
---|
| 686 | LOG10(X,PREC) ;
|
---|
| 687 | New L,LIM,M,N,K,VALUE
|
---|
| 688 | If X'>0 Set $Ecode=",M28,"
|
---|
| 689 | Set PREC=$Get(PREC,11)
|
---|
| 690 | Set M=1
|
---|
| 691 | ;
|
---|
| 692 | ;;;" If X>0 For N=0:1 Quit:(X/M)<10 Set M=M*10 ;" Number ~~
|
---|
| 693 | ;" Winfried Gerum (8 June 1995)
|
---|
| 694 | For N=0:1 Quit:(X/M)<10 Set M=M*10
|
---|
| 695 | ;;;
|
---|
| 696 | ;
|
---|
| 697 | If X<1 For N=0:-1 Quit:(X/M)>0.1 Set M=M*0.1
|
---|
| 698 | Set X=X/M
|
---|
| 699 | Set X=(X-1)/(X+1),(VALUE,L)=X
|
---|
| 700 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 701 | For K=3:2 Set L=L*X*X,M=L/K,VALUE=M+VALUE Set:M<0 M=-M Quit:M<LIM
|
---|
| 702 | Set VALUE=VALUE*2+(N*2.30258509298749)
|
---|
| 703 | Quit VALUE/2.30258509298749
|
---|
| 704 | ;===
|
---|
| 705 | ;
|
---|
| 706 | ;
|
---|
| 707 | MTXADD(A,B,R,ROWS,COLS) ;
|
---|
| 708 | ;" Add A[ROWS,COLS] to B[ROWS,COLS],
|
---|
| 709 | ;" result goes to R[ROWS,COLS]
|
---|
| 710 | IF $DATA(A)<10 QUIT 0
|
---|
| 711 | IF $DATA(B)<10 QUIT 0
|
---|
| 712 | IF $GET(ROWS)<1 QUIT 0
|
---|
| 713 | IF $GET(COLS)<1 QUIT 0
|
---|
| 714 | ;
|
---|
| 715 | NEW ROW,COL,ANY
|
---|
| 716 | FOR ROW=1:1:ROWS FOR COL=1:1:COLS DO
|
---|
| 717 | . KVALUE R(ROW,COL) SET ANY=0
|
---|
| 718 | . SET:$DATA(A(ROW,COL))#2 ANY=1
|
---|
| 719 | . SET:$DATA(B(ROW,COL))#2 ANY=1
|
---|
| 720 | . SET:ANY R(ROW,COL)=$GET(A(ROW,COL))+$GET(B(ROW,COL))
|
---|
| 721 | . QUIT
|
---|
| 722 | QUIT 1
|
---|
| 723 | ;===
|
---|
| 724 | ;
|
---|
| 725 | ;
|
---|
| 726 | MTXCOF(A,I,K,N) ;
|
---|
| 727 | ;" Compute cofactor for element [i,k]
|
---|
| 728 | ;" in matrix A[N,N]
|
---|
| 729 | NEW T,R,C,RR,CC
|
---|
| 730 | SET CC=0 FOR C=1:1:N DO:C'=K
|
---|
| 731 | . SET CC=CC+1,RR=0
|
---|
| 732 | . FOR R=1:1:N SET:R'=I RR=RR+1,T(RR,CC)=$GET(A(R,C))
|
---|
| 733 | . QUIT
|
---|
| 734 | QUIT $%MTXDET^MATH(.T,N-1)
|
---|
| 735 | ;===
|
---|
| 736 | ;
|
---|
| 737 | ;
|
---|
| 738 | MTXCOPY(A,R,ROWS,COLS) ;
|
---|
| 739 | ;" Copy A[ROWS,COLS] to R[ROWS,COLS]
|
---|
| 740 | IF $DATA(A)<10 QUIT 0
|
---|
| 741 | IF $GET(ROWS)<1 QUIT 0
|
---|
| 742 | IF $GET(COLS)<1 QUIT 0
|
---|
| 743 | ;
|
---|
| 744 | NEW ROW,COL
|
---|
| 745 | FOR ROW=1:1:ROWS FOR COL=1:1:COLS DO
|
---|
| 746 | . KVALUE R(ROW,COL)
|
---|
| 747 | . SET:$DATA(A(ROW,COL))#2 R(ROW,COL)=A(ROW,COL)
|
---|
| 748 | . QUIT
|
---|
| 749 | QUIT 1
|
---|
| 750 | ;===
|
---|
| 751 | ;
|
---|
| 752 | ;
|
---|
| 753 | MTXDET(A,N) ;
|
---|
| 754 | ;" Compute determinant of matrix A[N,N]
|
---|
| 755 | IF $DATA(A)<10 QUIT ""
|
---|
| 756 | IF $GET(N)<1 QUIT ""
|
---|
| 757 | ;
|
---|
| 758 | ;" First the simple cases
|
---|
| 759 | ;
|
---|
| 760 | IF N=1 QUIT $GET(A(1,1))
|
---|
| 761 | IF N=2 QUIT $GET(A(1,1))*$GET(A(2,2))-($GET(A(1,2))*$GET(A(2,1)))
|
---|
| 762 | ;
|
---|
| 763 | NEW DET,I,SIGN
|
---|
| 764 | ;
|
---|
| 765 | ;" Det A = sum (k=1:n) element (i,k) * cofactor [i,k]
|
---|
| 766 | ;
|
---|
| 767 | SET DET=0,SIGN=1
|
---|
| 768 | FOR I=1:1:N DO
|
---|
| 769 | . SET DET=$GET(A(1,I))*$%MTXCOF^MATH(.A,1,I,N)*SIGN+DET
|
---|
| 770 | . SET SIGN=-SIGN
|
---|
| 771 | . QUIT
|
---|
| 772 | QUIT DET
|
---|
| 773 | ;===
|
---|
| 774 | ;
|
---|
| 775 | ;
|
---|
| 776 | MTXEQU(A,B,R,N,M) ;
|
---|
| 777 | ;" Solve matrix equation A [M,M] * R [M,N] = B [M,N]
|
---|
| 778 | IF $GET(M)<1 QUIT ""
|
---|
| 779 | IF $GET(N)<1 QUIT ""
|
---|
| 780 | ;;;IF '$%MTXDET^MATH(.A) QUIT 0
|
---|
| 781 | ;" Ed de Moel, 29 August 1999
|
---|
| 782 | IF '$%MTXDET^MATH(.A,M) QUIT 0
|
---|
| 783 | ;;;
|
---|
| 784 | ;
|
---|
| 785 | NEW I,I1,J,J1,J2,K,L,T,T1,T2,TEMP,X
|
---|
| 786 | ;
|
---|
| 787 | SET X=$%MTXCOPY^MATH(.A,.T,N,N)
|
---|
| 788 | SET X=$%MTXCOPY^MATH(.B,.R,N,M)
|
---|
| 789 | ;
|
---|
| 790 | ;" Reduction of matrix A
|
---|
| 791 | ;" Steps of reduction are counted by index K
|
---|
| 792 | ;
|
---|
| 793 | FOR K=1:1:N-1 DO
|
---|
| 794 | . ;
|
---|
| 795 | . ;" Search for largest coefficient of T
|
---|
| 796 | . ;" (denoted by TEMP)
|
---|
| 797 | . ;" in first column of reduced system
|
---|
| 798 | . ;
|
---|
| 799 | . SET TEMP=0,J2=K
|
---|
| 800 | . FOR J1=K:1:N DO
|
---|
| 801 | . . QUIT:$TRANSLATE($GET(T(J1,K)),"-")>$TRANSLATE(TEMP,"-")
|
---|
| 802 | . . SET TEMP=T(J1,K),J2=J1
|
---|
| 803 | . . QUIT
|
---|
| 804 | . ;
|
---|
| 805 | . ;" Exchange row number K with row number J2,
|
---|
| 806 | . ;" if necessary
|
---|
| 807 | . ;
|
---|
| 808 | . DO:J2'=K
|
---|
| 809 | . . ;
|
---|
| 810 | . . FOR J=K:1:N DO
|
---|
| 811 | . . . SET T1=$GET(T(K,J)),T2=$GET(T(J2,J))
|
---|
| 812 | . . . KILL T(K,J),T(J2,J)
|
---|
| 813 | . . . IF T1'="" SET T(J2,J)=T1
|
---|
| 814 | . . . IF T2'="" SET T(K,J)=T2
|
---|
| 815 | . . . QUIT
|
---|
| 816 | . . FOR J=1:1:M DO
|
---|
| 817 | . . . SET T1=$GET(R(K,J)),T2=$GET(R(J2,J))
|
---|
| 818 | . . . KILL R(K,J),R(J2,J)
|
---|
| 819 | . . . IF T1'="" SET R(J2,J)=T1
|
---|
| 820 | . . . IF T2'="" SET R(K,J)=T2
|
---|
| 821 | . . . QUIT
|
---|
| 822 | . . QUIT
|
---|
| 823 | . ;
|
---|
| 824 | . ;" Actual reduction
|
---|
| 825 | . ;
|
---|
| 826 | . FOR I=K+1:1:N DO
|
---|
| 827 | . . FOR J=K+1:1:N DO
|
---|
| 828 | . . . QUIT:'$GET(T(K,K))
|
---|
| 829 | . . . SET T(I,J)=-$GET(T(K,J))*$GET(T(I,K))/T(K,K)+$GET(T(I,J))
|
---|
| 830 | . . . QUIT
|
---|
| 831 | . . FOR J=1:1:M DO
|
---|
| 832 | . . . QUIT:'$GET(T(K,K))
|
---|
| 833 | . . . SET R(I,J)=-$GET(R(K,J))*$GET(T(I,K))/T(K,K)+$GET(R(I,J))
|
---|
| 834 | . . . QUIT
|
---|
| 835 | . . QUIT
|
---|
| 836 | . QUIT
|
---|
| 837 | ;
|
---|
| 838 | ;" Backsubstitution
|
---|
| 839 | ;
|
---|
| 840 | FOR J=1:1:M DO
|
---|
| 841 | . IF $GET(T(N,N)) SET R(N,J)=$GET(R(N,J))/T(N,N)
|
---|
| 842 | . IF N-1>0 FOR I1=1:1:N-1 DO
|
---|
| 843 | . . SET I=N-I1
|
---|
| 844 | . . FOR L=I+1:1:N DO
|
---|
| 845 | . . . SET R(I,J)=-$GET(T(I,L))*$GET(R(L,J))+$GET(R(I,J))
|
---|
| 846 | . . . QUIT
|
---|
| 847 | . . IF $GET(T(I,I)) SET R(I,J)=$GET(R(I,J))/$GET(T(I,I))
|
---|
| 848 | . . QUIT
|
---|
| 849 | . QUIT
|
---|
| 850 | ;;;QUIT $%MTXDET^MATH(.R)
|
---|
| 851 | ;" Ed de Moel, 29 Aug 1999
|
---|
| 852 | QUIT $SELECT(M=N:$%MTXDET^MATH(.R,M),1:1)
|
---|
| 853 | ;;;
|
---|
| 854 | ;===
|
---|
| 855 | ;
|
---|
| 856 | MTXINV(A,R,N) ;
|
---|
| 857 | ;" Invert A[N,N], result goes to R[N,N]
|
---|
| 858 | IF $DATA(A)<10 QUIT 0
|
---|
| 859 | IF $GET(N)<1 QUIT 0
|
---|
| 860 | ;
|
---|
| 861 | NEW T,X
|
---|
| 862 | SET X=$%MTXUNIT^MATH(.T,N)
|
---|
| 863 | QUIT $%MTXEQU^MATH(.A,.T,.R,N,N)
|
---|
| 864 | ;===
|
---|
| 865 | ;
|
---|
| 866 | ;
|
---|
| 867 | MTXMUL(A,B,R,M,L,N) ;
|
---|
| 868 | ;" Multiply A[M,L] by B[L,N], result goes to R[M,N]
|
---|
| 869 | IF $DATA(A)<10 QUIT 0
|
---|
| 870 | IF $DATA(B)<10 QUIT 0
|
---|
| 871 | IF $GET(L)<1 QUIT 0
|
---|
| 872 | IF $GET(M)<1 QUIT 0
|
---|
| 873 | IF $GET(N)<1 QUIT 0
|
---|
| 874 | ;
|
---|
| 875 | NEW I,J,K,SUM,ANY
|
---|
| 876 | FOR I=1:1:M FOR J=1:1:N DO
|
---|
| 877 | . SET (SUM,ANY)=0
|
---|
| 878 | . KVALUE R(I,J)
|
---|
| 879 | . FOR K=1:1:L DO
|
---|
| 880 | . . SET:$DATA(A(I,K))#2 ANY=1
|
---|
| 881 | . . SET:$DATA(B(K,J))#2 ANY=1
|
---|
| 882 | . . SET SUM=$GET(A(I,K))*$GET(B(K,J))+SUM
|
---|
| 883 | . . QUIT
|
---|
| 884 | . SET:ANY R(I,J)=SUM
|
---|
| 885 | . QUIT
|
---|
| 886 | QUIT 1
|
---|
| 887 | ;===
|
---|
| 888 | ;
|
---|
| 889 | ;
|
---|
| 890 | MTXSCA(A,R,ROWS,COLS,S) ;
|
---|
| 891 | ;" Multiply A[ROWS,COLS] with the scalar S,
|
---|
| 892 | ;" result goes to R[ROWS,COLS]
|
---|
| 893 | IF $DATA(A)<10 QUIT 0
|
---|
| 894 | IF $GET(ROWS)<1 QUIT 0
|
---|
| 895 | IF $GET(COLS)<1 QUIT 0
|
---|
| 896 | IF '($DATA(S)#2) QUIT 0
|
---|
| 897 | ;
|
---|
| 898 | NEW ROW,COL
|
---|
| 899 | FOR ROW=1:1:ROWS FOR COL=1:1:COLS DO
|
---|
| 900 | . KVALUE R(ROW,COL)
|
---|
| 901 | . SET:$DATA(A(ROW,COL))#2 R(ROW,COL)=A(ROW,COL)*S
|
---|
| 902 | . QUIT
|
---|
| 903 | QUIT 1
|
---|
| 904 | ;===
|
---|
| 905 | ;
|
---|
| 906 | ;
|
---|
| 907 | MTXSUB(A,B,R,ROWS,COLS) ;
|
---|
| 908 | ;" Subtract B[ROWS,COLS] from A[ROWS,COLS],
|
---|
| 909 | ;" result goes to R[ROWS,COLS]
|
---|
| 910 | IF $DATA(A)<10 QUIT 0
|
---|
| 911 | IF $DATA(B)<10 QUIT 0
|
---|
| 912 | IF $GET(ROWS)<1 QUIT 0
|
---|
| 913 | IF $GET(COLS)<1 QUIT 0
|
---|
| 914 | ;
|
---|
| 915 | NEW ROW,COL,ANY
|
---|
| 916 | FOR ROW=1:1:ROWS FOR COL=1:1:COLS DO
|
---|
| 917 | . KVALUE R(ROW,COL) SET ANY=0
|
---|
| 918 | . SET:$DATA(A(ROW,COL))#2 ANY=1
|
---|
| 919 | . SET:$DATA(B(ROW,COL))#2 ANY=1
|
---|
| 920 | . ;
|
---|
| 921 | . ;;;" SET:ANY R(ROW,COL)=$GET(A(ROW,COL)-$GET(B(ROW,COL)) ;" Number ~~
|
---|
| 922 | . ;" Eli Reidler (28 June 1996)
|
---|
| 923 | . SET:ANY R(ROW,COL)=$GET(A(ROW,COL))-$GET(B(ROW,COL))
|
---|
| 924 | . ;;;
|
---|
| 925 | . ;
|
---|
| 926 | . QUIT
|
---|
| 927 | QUIT 1
|
---|
| 928 | ;===
|
---|
| 929 | ;
|
---|
| 930 | ;
|
---|
| 931 | MTXTRP(A,R,M,N) ;
|
---|
| 932 | ;" Transpose A[M,N], result goes to R[N,M]
|
---|
| 933 | IF $DATA(A)<10 QUIT 0
|
---|
| 934 | IF $GET(M)<1 QUIT 0
|
---|
| 935 | IF $GET(N)<1 QUIT 0
|
---|
| 936 | ;
|
---|
| 937 | NEW I,J,K,D1,V1,D2,V2
|
---|
| 938 | FOR I=1:1:M+N-1 FOR J=1:1:I+1\2 DO
|
---|
| 939 | . SET K=I-J+1
|
---|
| 940 | . IF K=J DO QUIT
|
---|
| 941 | . . SET V1=$GET(A(J,J)),D1=$DATA(A(J,J))#2
|
---|
| 942 | . . IF J'>N,J'>M KVALUE R(J,J) SET:D1 R(J,J)=V1
|
---|
| 943 | . . QUIT
|
---|
| 944 | . ;
|
---|
| 945 | . SET V1=$GET(A(K,J)),D1=$DATA(A(K,J))#2
|
---|
| 946 | . SET V2=$GET(A(J,K)),D2=$DATA(A(J,K))#2
|
---|
| 947 | . IF K'>M,J'>N KVALUE R(K,J) SET:D2 R(K,J)=V2
|
---|
| 948 | . IF J'>M,K'>N KVALUE R(J,K) SET:D1 R(J,K)=V1
|
---|
| 949 | . QUIT
|
---|
| 950 | QUIT 1
|
---|
| 951 | ;===
|
---|
| 952 | ;
|
---|
| 953 | ;
|
---|
| 954 | MTXUNIT(R,N,SPARSE) ;
|
---|
| 955 | ;" Create a unit matrix R[N,N]
|
---|
| 956 | IF $GET(N)<1 QUIT 0
|
---|
| 957 | ;
|
---|
| 958 | NEW ROW,COL
|
---|
| 959 | FOR ROW=1:1:N FOR COL=1:1:N DO
|
---|
| 960 | . KVALUE R(ROW,COL)
|
---|
| 961 | . IF $GET(SPARSE) QUIT:ROW'=COL
|
---|
| 962 | . SET R(ROW,COL)=$SELECT(ROW=COL:1,1:0)
|
---|
| 963 | . QUIT
|
---|
| 964 | QUIT 1
|
---|
| 965 | ;===
|
---|
| 966 | ;
|
---|
| 967 | ;
|
---|
| 968 | PI() Quit 3.14159265358979
|
---|
| 969 | ;===
|
---|
| 970 | ;
|
---|
| 971 | ;
|
---|
| 972 | PRODUCE(IN,SPEC,MAX) ;
|
---|
| 973 | NEW VALUE,AGAIN,P1,P2,I,COUNT
|
---|
| 974 | SET VALUE=IN,COUNT=0
|
---|
| 975 | FOR DO QUIT:'AGAIN
|
---|
| 976 | . SET AGAIN=0
|
---|
| 977 | . SET I=""
|
---|
| 978 | . FOR SET I=$ORDER(SPEC(I)) QUIT:I="" DO QUIT:COUNT<0
|
---|
| 979 | . . QUIT:$GET(SPEC(I,1))=""
|
---|
| 980 | . . QUIT:'($DATA(SPEC(I,2))#2)
|
---|
| 981 | . . FOR QUIT:VALUE'[SPEC(I,1) DO QUIT:COUNT<0
|
---|
| 982 | . . . SET P1=$PIECE(VALUE,SPEC(I,1),1)
|
---|
| 983 | . . . SET P2=$PIECE(VALUE,SPEC(I,1),2,$LENGTH(VALUE))
|
---|
| 984 | . . . SET VALUE=P1_SPEC(I,2)_P2,AGAIN=1
|
---|
| 985 | . . . SET COUNT=COUNT+1
|
---|
| 986 | . . . IF $DATA(MAX),COUNT>MAX SET COUNT=-1,AGAIN=0
|
---|
| 987 | . . . QUIT
|
---|
| 988 | . . QUIT
|
---|
| 989 | . QUIT
|
---|
| 990 | QUIT VALUE
|
---|
| 991 | ;===
|
---|
| 992 | ;
|
---|
| 993 | ;
|
---|
| 994 | RADDEG(X) Quit X*180/3.14159265358979
|
---|
| 995 | ;===
|
---|
| 996 | ;
|
---|
| 997 | ;
|
---|
| 998 | REPLACE(IN,SPEC) ;
|
---|
| 999 | NEW L,MASK,K,I,LT,F,VALUE
|
---|
| 1000 | SET L=$LENGTH(IN),MASK=$JUSTIFY("",L)
|
---|
| 1001 | SET I="" FOR SET I=$ORDER(SPEC(I)) QUIT:I="" DO
|
---|
| 1002 | . QUIT:'($DATA(SPEC(I,1))#2)
|
---|
| 1003 | . QUIT:SPEC(I,1)=""
|
---|
| 1004 | . QUIT:'($DATA(SPEC(I,2))#2)
|
---|
| 1005 | . SET LT=$LENGTH(SPEC(I,1))
|
---|
| 1006 | . SET F=0 FOR SET F=$FIND(IN,SPEC(I,1),F) QUIT:F<1 DO
|
---|
| 1007 | . . QUIT:$EXTRACT(MASK,F-LT,F-1)["X"
|
---|
| 1008 | . . SET VALUE(F-LT)=SPEC(I,2)
|
---|
| 1009 | . . SET $EXTRACT(MASK,F-LT,F-1)=$TRANSLATE($JUSTIFY("",LT)," ","X")
|
---|
| 1010 | . . QUIT
|
---|
| 1011 | . QUIT
|
---|
| 1012 | SET VALUE="" FOR K=1:1:L DO
|
---|
| 1013 | . IF $EXTRACT(MASK,K)=" " SET VALUE=VALUE_$EXTRACT(IN,K) QUIT
|
---|
| 1014 | . SET:$DATA(VALUE(K)) VALUE=VALUE_VALUE(K)
|
---|
| 1015 | . QUIT
|
---|
| 1016 | QUIT VALUE
|
---|
| 1017 | ;===
|
---|
| 1018 | ;
|
---|
| 1019 | ;
|
---|
| 1020 | SEC(X,PREC) ;
|
---|
| 1021 | New L,LIM,K,SIGN,VALUE
|
---|
| 1022 | ;
|
---|
| 1023 | ;;;" Set:X[":" X=$%DMSDEC^MATH(X,12) ;" Number ~~
|
---|
| 1024 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1025 | ;" Comment: The official description does not mention than
|
---|
| 1026 | ;" the function may also be called with the first
|
---|
| 1027 | ;" parameter in degrees, minutes and seconds.
|
---|
| 1028 | Set:X[":" X=$%DMSDEC^MATH(X)
|
---|
| 1029 | ;;;
|
---|
| 1030 | ;
|
---|
| 1031 | Set PREC=$Get(PREC,11)
|
---|
| 1032 | Set X=X#(2*$%PI^MATH())
|
---|
| 1033 | Set (VALUE,L)=1,SIGN=-1
|
---|
| 1034 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 1035 | For K=2:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 1036 | . Set L=L*X*X/(K-1*K),VALUE=VALUE+(SIGN*L)
|
---|
| 1037 | . Quit
|
---|
| 1038 | If 'VALUE Quit "INFINITE"
|
---|
| 1039 | Quit 1/VALUE
|
---|
| 1040 | ;===
|
---|
| 1041 | ;
|
---|
| 1042 | ;
|
---|
| 1043 | SECH(X,PREC) ;;;Quit 1/$%COSH^MATH(X,PREC) ;" Number ~~
|
---|
| 1044 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1045 | ;" Alan Frank (October 1995)
|
---|
| 1046 | Quit 1/$%COSH^MATH(X,$Get(PREC,11))
|
---|
| 1047 | ;;;
|
---|
| 1048 | ;===
|
---|
| 1049 | ;
|
---|
| 1050 | ;
|
---|
| 1051 | SIGN(X) Quit $SELECT(X<0:-1,X>0:1,1:0)
|
---|
| 1052 | ;===
|
---|
| 1053 | ;
|
---|
| 1054 | ;
|
---|
| 1055 | SIN(X,PREC) ;
|
---|
| 1056 | New L,LIM,K,SIGN,VALUE
|
---|
| 1057 | ;
|
---|
| 1058 | ;;;" Set:X[":" X=$%DMSDEC^MATH(X,12) ;" Number ~~
|
---|
| 1059 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1060 | ;" Comment: The official description does not mention than
|
---|
| 1061 | ;" the function may also be called with the first
|
---|
| 1062 | ;" parameter in degrees, minutes and seconds.
|
---|
| 1063 | Set:X[":" X=$%DMSDEC^MATH(X)
|
---|
| 1064 | ;;;
|
---|
| 1065 | ;
|
---|
| 1066 | Set PREC=$Get(PREC,11)
|
---|
| 1067 | Set X=X#(2*$%PI^MATH())
|
---|
| 1068 | Set (VALUE,L)=X,SIGN=-1
|
---|
| 1069 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 1070 | For K=3:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 1071 | . Set L=L/(K-1)*X/K*X,VALUE=VALUE+(SIGN*L)
|
---|
| 1072 | . Quit
|
---|
| 1073 | Quit VALUE
|
---|
| 1074 | ;===
|
---|
| 1075 | ;
|
---|
| 1076 | ;
|
---|
| 1077 | SIN(X) ;
|
---|
| 1078 | ;;;" ;" Number ~~
|
---|
| 1079 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1080 | ;" Comment: This version of the function is
|
---|
| 1081 | ;" optimized for speed, not for precision.
|
---|
| 1082 | ;" The 'precision' parameter is not supported,
|
---|
| 1083 | ;" and the precision is at best 1 in 10**-9.
|
---|
| 1084 | ;" Note that this function does not accept its
|
---|
| 1085 | ;" parameter in degrees, minutes and seconds.
|
---|
| 1086 | ;;;
|
---|
| 1087 | ;
|
---|
| 1088 | New A,N,PI,R,SIGN,XX
|
---|
| 1089 | ;
|
---|
| 1090 | ;" This approximation only works for 0 <= x <= pi/2
|
---|
| 1091 | ;" so reduce angle to correct quadrant.
|
---|
| 1092 | ;
|
---|
| 1093 | Set PI=$%PI^MATH(),X=X#(PI*2),SIGN=1
|
---|
| 1094 | Set:X>PI X=2*PI-X,SIGN=-1
|
---|
| 1095 | ;
|
---|
| 1096 | ;;;" Set:X*2<PI X=PI-X Set X=-PI/2+2 ;" Number ~~
|
---|
| 1097 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1098 | Set:X*2<PI X=PI-X
|
---|
| 1099 | ;;;
|
---|
| 1100 | ;
|
---|
| 1101 | ;
|
---|
| 1102 | Set XX=X*X,A(1)=-0.4999999963,A(2)=0.0416666418
|
---|
| 1103 | Set A(3)=-0.0013888397,A(4)=0.0000247609,A(5)=-0.0000002605
|
---|
| 1104 | Set (X,R)=1 For N=1:1:5 Set X=X*XX,R=A(N)*X+R
|
---|
| 1105 | Quit R*SIGN
|
---|
| 1106 | ;===
|
---|
| 1107 | ;
|
---|
| 1108 | ;
|
---|
| 1109 | SINH(X,PREC) ;
|
---|
| 1110 | ;
|
---|
| 1111 | ;;;" New F,I,P,R,T,XX ;" Number ~~
|
---|
| 1112 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1113 | ;" Eli Reidler (28 June 1996)
|
---|
| 1114 | New E,F,I,P,R,T,XX
|
---|
| 1115 | ;;;
|
---|
| 1116 | ;
|
---|
| 1117 | Set PREC=$Get(PREC,11)+1
|
---|
| 1118 | Set @("E=1E-"_PREC)
|
---|
| 1119 | Set XX=X*X,F=1,I=2,(P,R,T)=X
|
---|
| 1120 | For Set T=T*XX,F=I+1*I*F,R=T/F+R,P=P-R/R,I=I+2 If -E<P,P<E Quit
|
---|
| 1121 | Quit R
|
---|
| 1122 | ;===
|
---|
| 1123 | ;
|
---|
| 1124 | ;
|
---|
| 1125 | SQRT(X,PREC) ;
|
---|
| 1126 | If X<0 Set $Ecode=",M28,"
|
---|
| 1127 | If X=0 Quit 0
|
---|
| 1128 | ;
|
---|
| 1129 | ;;;" ;" Number ~~
|
---|
| 1130 | ;" Alan Frank (October 1995)
|
---|
| 1131 | Set PREC=$Get(PREC,11)
|
---|
| 1132 | ;;;
|
---|
| 1133 | ;
|
---|
| 1134 | ;
|
---|
| 1135 | ;;;" If X<1 Quit 1/$%SQRT^MATH(1/X) ;" Number ~~
|
---|
| 1136 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1137 | If X<1 Quit 1/$%SQRT^MATH(1/X,PREC)
|
---|
| 1138 | ;;;
|
---|
| 1139 | ;
|
---|
| 1140 | New P,R,E
|
---|
| 1141 | Set PREC=$Get(PREC,11)+1
|
---|
| 1142 | ;
|
---|
| 1143 | ;;;" Set @(E="1E-"_PREC) ;" Number ~~
|
---|
| 1144 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1145 | ;" Eli Reidler (28 June 1996)
|
---|
| 1146 | Set @("E=1E-"_PREC)
|
---|
| 1147 | ;;;
|
---|
| 1148 | ;
|
---|
| 1149 | Set R=X
|
---|
| 1150 | For Set P=R,R=X/R+R/2,P=P-R/R If -E<P,P<E Quit
|
---|
| 1151 | Quit R
|
---|
| 1152 | ;===
|
---|
| 1153 | ;
|
---|
| 1154 | ;
|
---|
| 1155 | TAN(X,PREC) ;
|
---|
| 1156 | New L,LIM,K,S,SIGN,VALUE
|
---|
| 1157 | ;
|
---|
| 1158 | ;;;" Set:X[":" X=$%DMSDEC^MATH(X,12) ;" Number ~~
|
---|
| 1159 | ;" Winfried Gerum (8 June 1995)
|
---|
| 1160 | ;" Comment: The official description does not mention than
|
---|
| 1161 | ;" the function may also be called with the first
|
---|
| 1162 | ;" parameter in degrees, minutes and seconds.
|
---|
| 1163 | Set:X[":" X=$%DMSDEC^MATH(X)
|
---|
| 1164 | ;;;
|
---|
| 1165 | ;
|
---|
| 1166 | Set PREC=$Get(PREC,11)
|
---|
| 1167 | Set X=X#(2*$%PI^MATH())
|
---|
| 1168 | Set (VALUE,L)=X,SIGN=-1
|
---|
| 1169 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 1170 | For K=3:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 1171 | . Set L=L/(K-1)*X/K*X,VALUE=VALUE+(SIGN*L)
|
---|
| 1172 | . Quit
|
---|
| 1173 | Set S=VALUE
|
---|
| 1174 | Set X=X#(2*$%PI^MATH())
|
---|
| 1175 | Set (VALUE,L)=1,SIGN=-1
|
---|
| 1176 | Set LIM=$Select((PREC+3)'>11:PREC+3,1:11),@("LIM=1E-"_LIM)
|
---|
| 1177 | For K=2:2 Do Quit:($Translate(L,"-")<LIM) Set SIGN=SIGN*-1
|
---|
| 1178 | . Set L=L*X*X/(K-1*K),VALUE=VALUE+(SIGN*L)
|
---|
| 1179 | . Quit
|
---|
| 1180 | If 'VALUE Quit "INFINITE"
|
---|
| 1181 | Quit S/VALUE
|
---|
| 1182 | ;===
|
---|
| 1183 | ;
|
---|
| 1184 | ;
|
---|
| 1185 | TANH(X,PREC) ;
|
---|
| 1186 | ;
|
---|
| 1187 | ;;;" ;" Number ~~
|
---|
| 1188 | ;" Alan Frank (October 1995)
|
---|
| 1189 | Set PREC=$Get(PREC,11)
|
---|
| 1190 | ;;;
|
---|
| 1191 | ;
|
---|
| 1192 | Quit $%SINH^MATH(X,PREC)/$%COSH^MATH(X,PREC)
|
---|
| 1193 | ;===
|
---|
| 1194 | ;
|
---|
| 1195 | ;
|
---|