1*DECK ZBKNU 2 SUBROUTINE ZBKNU (ZR, ZI, FNU, KODE, N, YR, YI, NZ, TOL, ELIM, 3 + ALIM) 4C***BEGIN PROLOGUE ZBKNU 5C***SUBSIDIARY 6C***PURPOSE Subsidiary to ZAIRY, ZBESH, ZBESI and ZBESK 7C***LIBRARY SLATEC 8C***TYPE ALL (CBKNU-A, ZBKNU-A) 9C***AUTHOR Amos, D. E., (SNL) 10C***DESCRIPTION 11C 12C ZBKNU COMPUTES THE K BESSEL FUNCTION IN THE RIGHT HALF Z PLANE. 13C 14C***SEE ALSO ZAIRY, ZBESH, ZBESI, ZBESK 15C***ROUTINES CALLED D1MACH, DGAMLN, I1MACH, ZABS, ZDIV, ZEXP, ZKSCL, 16C ZLOG, ZMLT, ZSHCH, ZSQRT, ZUCHK 17C***REVISION HISTORY (YYMMDD) 18C 830501 DATE WRITTEN 19C 910415 Prologue converted to Version 4.0 format. (BAB) 20C 930122 Added ZEXP, ZLOG and ZSQRT to EXTERNAL statement. (RWC) 21C***END PROLOGUE ZBKNU 22C 23 DOUBLE PRECISION AA, AK, ALIM, ASCLE, A1, A2, BB, BK, BRY, CAZ, 24 * CBI, CBR, CC, CCHI, CCHR, CKI, CKR, COEFI, COEFR, CONEI, CONER, 25 * CRSCR, CSCLR, CSHI, CSHR, CSI, CSR, CSRR, CSSR, CTWOR, 26 * CZEROI, CZEROR, CZI, CZR, DNU, DNU2, DPI, ELIM, ETEST, FC, FHS, 27 * FI, FK, FKS, FMUI, FMUR, FNU, FPI, FR, G1, G2, HPI, PI, PR, PTI, 28 * PTR, P1I, P1R, P2I, P2M, P2R, QI, QR, RAK, RCAZ, RTHPI, RZI, 29 * RZR, R1, S, SMUI, SMUR, SPI, STI, STR, S1I, S1R, S2I, S2R, TM, 30 * TOL, TTH, T1, T2, YI, YR, ZI, ZR, DGAMLN, D1MACH, ZABS, ELM, 31 * CELMR, ZDR, ZDI, AS, ALAS, HELIM, CYR, CYI 32 INTEGER I, IFLAG, INU, K, KFLAG, KK, KMAX, KODE, KODED, N, NZ, 33 * IDUM, I1MACH, J, IC, INUB, NW 34 DIMENSION YR(N), YI(N), CC(8), CSSR(3), CSRR(3), BRY(3), CYR(2), 35 * CYI(2) 36 EXTERNAL ZABS, ZEXP, ZLOG, ZSQRT 37C COMPLEX Z,Y,A,B,RZ,SMU,FU,FMU,F,FLRZ,CZ,S1,S2,CSH,CCH 38C COMPLEX CK,P,Q,COEF,P1,P2,CBK,PT,CZERO,CONE,CTWO,ST,EZ,CS,DK 39C 40 DATA KMAX / 30 / 41 DATA CZEROR,CZEROI,CONER,CONEI,CTWOR,R1/ 42 1 0.0D0 , 0.0D0 , 1.0D0 , 0.0D0 , 2.0D0 , 2.0D0 / 43 DATA DPI, RTHPI, SPI ,HPI, FPI, TTH / 44 1 3.14159265358979324D0, 1.25331413731550025D0, 45 2 1.90985931710274403D0, 1.57079632679489662D0, 46 3 1.89769999331517738D0, 6.66666666666666666D-01/ 47 DATA CC(1), CC(2), CC(3), CC(4), CC(5), CC(6), CC(7), CC(8)/ 48 1 5.77215664901532861D-01, -4.20026350340952355D-02, 49 2 -4.21977345555443367D-02, 7.21894324666309954D-03, 50 3 -2.15241674114950973D-04, -2.01348547807882387D-05, 51 4 1.13302723198169588D-06, 6.11609510448141582D-09/ 52C***FIRST EXECUTABLE STATEMENT ZBKNU 53 CAZ = ZABS(ZR,ZI) 54 CSCLR = 1.0D0/TOL 55 CRSCR = TOL 56 CSSR(1) = CSCLR 57 CSSR(2) = 1.0D0 58 CSSR(3) = CRSCR 59 CSRR(1) = CRSCR 60 CSRR(2) = 1.0D0 61 CSRR(3) = CSCLR 62 BRY(1) = 1.0D+3*D1MACH(1)/TOL 63 BRY(2) = 1.0D0/BRY(1) 64 BRY(3) = D1MACH(2) 65 NZ = 0 66 IFLAG = 0 67 KODED = KODE 68 RCAZ = 1.0D0/CAZ 69 STR = ZR*RCAZ 70 STI = -ZI*RCAZ 71 RZR = (STR+STR)*RCAZ 72 RZI = (STI+STI)*RCAZ 73 INU = FNU+0.5D0 74 DNU = FNU - INU 75 IF (ABS(DNU).EQ.0.5D0) GO TO 110 76 DNU2 = 0.0D0 77 IF (ABS(DNU).GT.TOL) DNU2 = DNU*DNU 78 IF (CAZ.GT.R1) GO TO 110 79C----------------------------------------------------------------------- 80C SERIES FOR ABS(Z).LE.R1 81C----------------------------------------------------------------------- 82 FC = 1.0D0 83 CALL ZLOG(RZR, RZI, SMUR, SMUI, IDUM) 84 FMUR = SMUR*DNU 85 FMUI = SMUI*DNU 86 CALL ZSHCH(FMUR, FMUI, CSHR, CSHI, CCHR, CCHI) 87 IF (DNU.EQ.0.0D0) GO TO 10 88 FC = DNU*DPI 89 FC = FC/SIN(FC) 90 SMUR = CSHR/DNU 91 SMUI = CSHI/DNU 92 10 CONTINUE 93 A2 = 1.0D0 + DNU 94C----------------------------------------------------------------------- 95C GAM(1-Z)*GAM(1+Z)=PI*Z/SIN(PI*Z), T1=1/GAM(1-DNU), T2=1/GAM(1+DNU) 96C----------------------------------------------------------------------- 97 T2 = EXP(-DGAMLN(A2,IDUM)) 98 T1 = 1.0D0/(T2*FC) 99 IF (ABS(DNU).GT.0.1D0) GO TO 40 100C----------------------------------------------------------------------- 101C SERIES FOR F0 TO RESOLVE INDETERMINACY FOR SMALL ABS(DNU) 102C----------------------------------------------------------------------- 103 AK = 1.0D0 104 S = CC(1) 105 DO 20 K=2,8 106 AK = AK*DNU2 107 TM = CC(K)*AK 108 S = S + TM 109 IF (ABS(TM).LT.TOL) GO TO 30 110 20 CONTINUE 111 30 G1 = -S 112 GO TO 50 113 40 CONTINUE 114 G1 = (T1-T2)/(DNU+DNU) 115 50 CONTINUE 116 G2 = (T1+T2)*0.5D0 117 FR = FC*(CCHR*G1+SMUR*G2) 118 FI = FC*(CCHI*G1+SMUI*G2) 119 CALL ZEXP(FMUR, FMUI, STR, STI) 120 PR = 0.5D0*STR/T2 121 PI = 0.5D0*STI/T2 122 CALL ZDIV(0.5D0, 0.0D0, STR, STI, PTR, PTI) 123 QR = PTR/T1 124 QI = PTI/T1 125 S1R = FR 126 S1I = FI 127 S2R = PR 128 S2I = PI 129 AK = 1.0D0 130 A1 = 1.0D0 131 CKR = CONER 132 CKI = CONEI 133 BK = 1.0D0 - DNU2 134 IF (INU.GT.0 .OR. N.GT.1) GO TO 80 135C----------------------------------------------------------------------- 136C GENERATE K(FNU,Z), 0.0D0 .LE. FNU .LT. 0.5D0 AND N=1 137C----------------------------------------------------------------------- 138 IF (CAZ.LT.TOL) GO TO 70 139 CALL ZMLT(ZR, ZI, ZR, ZI, CZR, CZI) 140 CZR = 0.25D0*CZR 141 CZI = 0.25D0*CZI 142 T1 = 0.25D0*CAZ*CAZ 143 60 CONTINUE 144 FR = (FR*AK+PR+QR)/BK 145 FI = (FI*AK+PI+QI)/BK 146 STR = 1.0D0/(AK-DNU) 147 PR = PR*STR 148 PI = PI*STR 149 STR = 1.0D0/(AK+DNU) 150 QR = QR*STR 151 QI = QI*STR 152 STR = CKR*CZR - CKI*CZI 153 RAK = 1.0D0/AK 154 CKI = (CKR*CZI+CKI*CZR)*RAK 155 CKR = STR*RAK 156 S1R = CKR*FR - CKI*FI + S1R 157 S1I = CKR*FI + CKI*FR + S1I 158 A1 = A1*T1*RAK 159 BK = BK + AK + AK + 1.0D0 160 AK = AK + 1.0D0 161 IF (A1.GT.TOL) GO TO 60 162 70 CONTINUE 163 YR(1) = S1R 164 YI(1) = S1I 165 IF (KODED.EQ.1) RETURN 166 CALL ZEXP(ZR, ZI, STR, STI) 167 CALL ZMLT(S1R, S1I, STR, STI, YR(1), YI(1)) 168 RETURN 169C----------------------------------------------------------------------- 170C GENERATE K(DNU,Z) AND K(DNU+1,Z) FOR FORWARD RECURRENCE 171C----------------------------------------------------------------------- 172 80 CONTINUE 173 IF (CAZ.LT.TOL) GO TO 100 174 CALL ZMLT(ZR, ZI, ZR, ZI, CZR, CZI) 175 CZR = 0.25D0*CZR 176 CZI = 0.25D0*CZI 177 T1 = 0.25D0*CAZ*CAZ 178 90 CONTINUE 179 FR = (FR*AK+PR+QR)/BK 180 FI = (FI*AK+PI+QI)/BK 181 STR = 1.0D0/(AK-DNU) 182 PR = PR*STR 183 PI = PI*STR 184 STR = 1.0D0/(AK+DNU) 185 QR = QR*STR 186 QI = QI*STR 187 STR = CKR*CZR - CKI*CZI 188 RAK = 1.0D0/AK 189 CKI = (CKR*CZI+CKI*CZR)*RAK 190 CKR = STR*RAK 191 S1R = CKR*FR - CKI*FI + S1R 192 S1I = CKR*FI + CKI*FR + S1I 193 STR = PR - FR*AK 194 STI = PI - FI*AK 195 S2R = CKR*STR - CKI*STI + S2R 196 S2I = CKR*STI + CKI*STR + S2I 197 A1 = A1*T1*RAK 198 BK = BK + AK + AK + 1.0D0 199 AK = AK + 1.0D0 200 IF (A1.GT.TOL) GO TO 90 201 100 CONTINUE 202 KFLAG = 2 203 A1 = FNU + 1.0D0 204 AK = A1*ABS(SMUR) 205 IF (AK.GT.ALIM) KFLAG = 3 206 STR = CSSR(KFLAG) 207 P2R = S2R*STR 208 P2I = S2I*STR 209 CALL ZMLT(P2R, P2I, RZR, RZI, S2R, S2I) 210 S1R = S1R*STR 211 S1I = S1I*STR 212 IF (KODED.EQ.1) GO TO 210 213 CALL ZEXP(ZR, ZI, FR, FI) 214 CALL ZMLT(S1R, S1I, FR, FI, S1R, S1I) 215 CALL ZMLT(S2R, S2I, FR, FI, S2R, S2I) 216 GO TO 210 217C----------------------------------------------------------------------- 218C IFLAG=0 MEANS NO UNDERFLOW OCCURRED 219C IFLAG=1 MEANS AN UNDERFLOW OCCURRED- COMPUTATION PROCEEDS WITH 220C KODED=2 AND A TEST FOR ON SCALE VALUES IS MADE DURING FORWARD 221C RECURSION 222C----------------------------------------------------------------------- 223 110 CONTINUE 224 CALL ZSQRT(ZR, ZI, STR, STI) 225 CALL ZDIV(RTHPI, CZEROI, STR, STI, COEFR, COEFI) 226 KFLAG = 2 227 IF (KODED.EQ.2) GO TO 120 228 IF (ZR.GT.ALIM) GO TO 290 229C BLANK LINE 230 STR = EXP(-ZR)*CSSR(KFLAG) 231 STI = -STR*SIN(ZI) 232 STR = STR*COS(ZI) 233 CALL ZMLT(COEFR, COEFI, STR, STI, COEFR, COEFI) 234 120 CONTINUE 235 IF (ABS(DNU).EQ.0.5D0) GO TO 300 236C----------------------------------------------------------------------- 237C MILLER ALGORITHM FOR ABS(Z).GT.R1 238C----------------------------------------------------------------------- 239 AK = COS(DPI*DNU) 240 AK = ABS(AK) 241 IF (AK.EQ.CZEROR) GO TO 300 242 FHS = ABS(0.25D0-DNU2) 243 IF (FHS.EQ.CZEROR) GO TO 300 244C----------------------------------------------------------------------- 245C COMPUTE R2=F(E). IF ABS(Z).GE.R2, USE FORWARD RECURRENCE TO 246C DETERMINE THE BACKWARD INDEX K. R2=F(E) IS A STRAIGHT LINE ON 247C 12.LE.E.LE.60. E IS COMPUTED FROM 2**(-E)=B**(1-I1MACH(14))= 248C TOL WHERE B IS THE BASE OF THE ARITHMETIC. 249C----------------------------------------------------------------------- 250 T1 = I1MACH(14)-1 251 T1 = T1*D1MACH(5)*3.321928094D0 252 T1 = MAX(T1,12.0D0) 253 T1 = MIN(T1,60.0D0) 254 T2 = TTH*T1 - 6.0D0 255 IF (ZR.NE.0.0D0) GO TO 130 256 T1 = HPI 257 GO TO 140 258 130 CONTINUE 259 T1 = DATAN(ZI/ZR) 260 T1 = ABS(T1) 261 140 CONTINUE 262 IF (T2.GT.CAZ) GO TO 170 263C----------------------------------------------------------------------- 264C FORWARD RECURRENCE LOOP WHEN ABS(Z).GE.R2 265C----------------------------------------------------------------------- 266 ETEST = AK/(DPI*CAZ*TOL) 267 FK = CONER 268 IF (ETEST.LT.CONER) GO TO 180 269 FKS = CTWOR 270 CKR = CAZ + CAZ + CTWOR 271 P1R = CZEROR 272 P2R = CONER 273 DO 150 I=1,KMAX 274 AK = FHS/FKS 275 CBR = CKR/(FK+CONER) 276 PTR = P2R 277 P2R = CBR*P2R - P1R*AK 278 P1R = PTR 279 CKR = CKR + CTWOR 280 FKS = FKS + FK + FK + CTWOR 281 FHS = FHS + FK + FK 282 FK = FK + CONER 283 STR = ABS(P2R)*FK 284 IF (ETEST.LT.STR) GO TO 160 285 150 CONTINUE 286 GO TO 310 287 160 CONTINUE 288 FK = FK + SPI*T1*SQRT(T2/CAZ) 289 FHS = ABS(0.25D0-DNU2) 290 GO TO 180 291 170 CONTINUE 292C----------------------------------------------------------------------- 293C COMPUTE BACKWARD INDEX K FOR ABS(Z).LT.R2 294C----------------------------------------------------------------------- 295 A2 = SQRT(CAZ) 296 AK = FPI*AK/(TOL*SQRT(A2)) 297 AA = 3.0D0*T1/(1.0D0+CAZ) 298 BB = 14.7D0*T1/(28.0D0+CAZ) 299 AK = (LOG(AK)+CAZ*COS(AA)/(1.0D0+0.008D0*CAZ))/COS(BB) 300 FK = 0.12125D0*AK*AK/CAZ + 1.5D0 301 180 CONTINUE 302C----------------------------------------------------------------------- 303C BACKWARD RECURRENCE LOOP FOR MILLER ALGORITHM 304C----------------------------------------------------------------------- 305 K = FK 306 FK = K 307 FKS = FK*FK 308 P1R = CZEROR 309 P1I = CZEROI 310 P2R = TOL 311 P2I = CZEROI 312 CSR = P2R 313 CSI = P2I 314 DO 190 I=1,K 315 A1 = FKS - FK 316 AK = (FKS+FK)/(A1+FHS) 317 RAK = 2.0D0/(FK+CONER) 318 CBR = (FK+ZR)*RAK 319 CBI = ZI*RAK 320 PTR = P2R 321 PTI = P2I 322 P2R = (PTR*CBR-PTI*CBI-P1R)*AK 323 P2I = (PTI*CBR+PTR*CBI-P1I)*AK 324 P1R = PTR 325 P1I = PTI 326 CSR = CSR + P2R 327 CSI = CSI + P2I 328 FKS = A1 - FK + CONER 329 FK = FK - CONER 330 190 CONTINUE 331C----------------------------------------------------------------------- 332C COMPUTE (P2/CS)=(P2/ABS(CS))*(CONJG(CS)/ABS(CS)) FOR BETTER 333C SCALING 334C----------------------------------------------------------------------- 335 TM = ZABS(CSR,CSI) 336 PTR = 1.0D0/TM 337 S1R = P2R*PTR 338 S1I = P2I*PTR 339 CSR = CSR*PTR 340 CSI = -CSI*PTR 341 CALL ZMLT(COEFR, COEFI, S1R, S1I, STR, STI) 342 CALL ZMLT(STR, STI, CSR, CSI, S1R, S1I) 343 IF (INU.GT.0 .OR. N.GT.1) GO TO 200 344 ZDR = ZR 345 ZDI = ZI 346 IF(IFLAG.EQ.1) GO TO 270 347 GO TO 240 348 200 CONTINUE 349C----------------------------------------------------------------------- 350C COMPUTE P1/P2=(P1/ABS(P2)*CONJG(P2)/ABS(P2) FOR SCALING 351C----------------------------------------------------------------------- 352 TM = ZABS(P2R,P2I) 353 PTR = 1.0D0/TM 354 P1R = P1R*PTR 355 P1I = P1I*PTR 356 P2R = P2R*PTR 357 P2I = -P2I*PTR 358 CALL ZMLT(P1R, P1I, P2R, P2I, PTR, PTI) 359 STR = DNU + 0.5D0 - PTR 360 STI = -PTI 361 CALL ZDIV(STR, STI, ZR, ZI, STR, STI) 362 STR = STR + 1.0D0 363 CALL ZMLT(STR, STI, S1R, S1I, S2R, S2I) 364C----------------------------------------------------------------------- 365C FORWARD RECURSION ON THE THREE TERM RECURSION WITH RELATION WITH 366C SCALING NEAR EXPONENT EXTREMES ON KFLAG=1 OR KFLAG=3 367C----------------------------------------------------------------------- 368 210 CONTINUE 369 STR = DNU + 1.0D0 370 CKR = STR*RZR 371 CKI = STR*RZI 372 IF (N.EQ.1) INU = INU - 1 373 IF (INU.GT.0) GO TO 220 374 IF (N.GT.1) GO TO 215 375 S1R = S2R 376 S1I = S2I 377 215 CONTINUE 378 ZDR = ZR 379 ZDI = ZI 380 IF(IFLAG.EQ.1) GO TO 270 381 GO TO 240 382 220 CONTINUE 383 INUB = 1 384 IF(IFLAG.EQ.1) GO TO 261 385 225 CONTINUE 386 P1R = CSRR(KFLAG) 387 ASCLE = BRY(KFLAG) 388 DO 230 I=INUB,INU 389 STR = S2R 390 STI = S2I 391 S2R = CKR*STR - CKI*STI + S1R 392 S2I = CKR*STI + CKI*STR + S1I 393 S1R = STR 394 S1I = STI 395 CKR = CKR + RZR 396 CKI = CKI + RZI 397 IF (KFLAG.GE.3) GO TO 230 398 P2R = S2R*P1R 399 P2I = S2I*P1R 400 STR = ABS(P2R) 401 STI = ABS(P2I) 402 P2M = MAX(STR,STI) 403 IF (P2M.LE.ASCLE) GO TO 230 404 KFLAG = KFLAG + 1 405 ASCLE = BRY(KFLAG) 406 S1R = S1R*P1R 407 S1I = S1I*P1R 408 S2R = P2R 409 S2I = P2I 410 STR = CSSR(KFLAG) 411 S1R = S1R*STR 412 S1I = S1I*STR 413 S2R = S2R*STR 414 S2I = S2I*STR 415 P1R = CSRR(KFLAG) 416 230 CONTINUE 417 IF (N.NE.1) GO TO 240 418 S1R = S2R 419 S1I = S2I 420 240 CONTINUE 421 STR = CSRR(KFLAG) 422 YR(1) = S1R*STR 423 YI(1) = S1I*STR 424 IF (N.EQ.1) RETURN 425 YR(2) = S2R*STR 426 YI(2) = S2I*STR 427 IF (N.EQ.2) RETURN 428 KK = 2 429 250 CONTINUE 430 KK = KK + 1 431 IF (KK.GT.N) RETURN 432 P1R = CSRR(KFLAG) 433 ASCLE = BRY(KFLAG) 434 DO 260 I=KK,N 435 P2R = S2R 436 P2I = S2I 437 S2R = CKR*P2R - CKI*P2I + S1R 438 S2I = CKI*P2R + CKR*P2I + S1I 439 S1R = P2R 440 S1I = P2I 441 CKR = CKR + RZR 442 CKI = CKI + RZI 443 P2R = S2R*P1R 444 P2I = S2I*P1R 445 YR(I) = P2R 446 YI(I) = P2I 447 IF (KFLAG.GE.3) GO TO 260 448 STR = ABS(P2R) 449 STI = ABS(P2I) 450 P2M = MAX(STR,STI) 451 IF (P2M.LE.ASCLE) GO TO 260 452 KFLAG = KFLAG + 1 453 ASCLE = BRY(KFLAG) 454 S1R = S1R*P1R 455 S1I = S1I*P1R 456 S2R = P2R 457 S2I = P2I 458 STR = CSSR(KFLAG) 459 S1R = S1R*STR 460 S1I = S1I*STR 461 S2R = S2R*STR 462 S2I = S2I*STR 463 P1R = CSRR(KFLAG) 464 260 CONTINUE 465 RETURN 466C----------------------------------------------------------------------- 467C IFLAG=1 CASES, FORWARD RECURRENCE ON SCALED VALUES ON UNDERFLOW 468C----------------------------------------------------------------------- 469 261 CONTINUE 470 HELIM = 0.5D0*ELIM 471 ELM = EXP(-ELIM) 472 CELMR = ELM 473 ASCLE = BRY(1) 474 ZDR = ZR 475 ZDI = ZI 476 IC = -1 477 J = 2 478 DO 262 I=1,INU 479 STR = S2R 480 STI = S2I 481 S2R = STR*CKR-STI*CKI+S1R 482 S2I = STI*CKR+STR*CKI+S1I 483 S1R = STR 484 S1I = STI 485 CKR = CKR+RZR 486 CKI = CKI+RZI 487 AS = ZABS(S2R,S2I) 488 ALAS = LOG(AS) 489 P2R = -ZDR+ALAS 490 IF(P2R.LT.(-ELIM)) GO TO 263 491 CALL ZLOG(S2R,S2I,STR,STI,IDUM) 492 P2R = -ZDR+STR 493 P2I = -ZDI+STI 494 P2M = EXP(P2R)/TOL 495 P1R = P2M*COS(P2I) 496 P1I = P2M*SIN(P2I) 497 CALL ZUCHK(P1R,P1I,NW,ASCLE,TOL) 498 IF(NW.NE.0) GO TO 263 499 J = 3 - J 500 CYR(J) = P1R 501 CYI(J) = P1I 502 IF(IC.EQ.(I-1)) GO TO 264 503 IC = I 504 GO TO 262 505 263 CONTINUE 506 IF(ALAS.LT.HELIM) GO TO 262 507 ZDR = ZDR-ELIM 508 S1R = S1R*CELMR 509 S1I = S1I*CELMR 510 S2R = S2R*CELMR 511 S2I = S2I*CELMR 512 262 CONTINUE 513 IF(N.NE.1) GO TO 270 514 S1R = S2R 515 S1I = S2I 516 GO TO 270 517 264 CONTINUE 518 KFLAG = 1 519 INUB = I+1 520 S2R = CYR(J) 521 S2I = CYI(J) 522 J = 3 - J 523 S1R = CYR(J) 524 S1I = CYI(J) 525 IF(INUB.LE.INU) GO TO 225 526 IF(N.NE.1) GO TO 240 527 S1R = S2R 528 S1I = S2I 529 GO TO 240 530 270 CONTINUE 531 YR(1) = S1R 532 YI(1) = S1I 533 IF(N.EQ.1) GO TO 280 534 YR(2) = S2R 535 YI(2) = S2I 536 280 CONTINUE 537 ASCLE = BRY(1) 538 CALL ZKSCL(ZDR,ZDI,FNU,N,YR,YI,NZ,RZR,RZI,ASCLE,TOL,ELIM) 539 INU = N - NZ 540 IF (INU.LE.0) RETURN 541 KK = NZ + 1 542 S1R = YR(KK) 543 S1I = YI(KK) 544 YR(KK) = S1R*CSRR(1) 545 YI(KK) = S1I*CSRR(1) 546 IF (INU.EQ.1) RETURN 547 KK = NZ + 2 548 S2R = YR(KK) 549 S2I = YI(KK) 550 YR(KK) = S2R*CSRR(1) 551 YI(KK) = S2I*CSRR(1) 552 IF (INU.EQ.2) RETURN 553 T2 = FNU + (KK-1) 554 CKR = T2*RZR 555 CKI = T2*RZI 556 KFLAG = 1 557 GO TO 250 558 290 CONTINUE 559C----------------------------------------------------------------------- 560C SCALE BY EXP(Z), IFLAG = 1 CASES 561C----------------------------------------------------------------------- 562 KODED = 2 563 IFLAG = 1 564 KFLAG = 2 565 GO TO 120 566C----------------------------------------------------------------------- 567C FNU=HALF ODD INTEGER CASE, DNU=-0.5 568C----------------------------------------------------------------------- 569 300 CONTINUE 570 S1R = COEFR 571 S1I = COEFI 572 S2R = COEFR 573 S2I = COEFI 574 GO TO 210 575C 576C 577 310 CONTINUE 578 NZ=-2 579 RETURN 580 END 581