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