|
3217
|
1 SUBROUTINE ZUNK1(ZR, ZI, FNU, KODE, MR, N, YR, YI, NZ, TOL, ELIM, |
|
|
2 * ALIM) |
|
|
3 C***BEGIN PROLOGUE ZUNK1 |
|
|
4 C***REFER TO ZBESK |
|
|
5 C |
|
|
6 C ZUNK1 COMPUTES K(FNU,Z) AND ITS ANALYTIC CONTINUATION FROM THE |
|
|
7 C RIGHT HALF PLANE TO THE LEFT HALF PLANE BY MEANS OF THE |
|
|
8 C UNIFORM ASYMPTOTIC EXPANSION. |
|
|
9 C MR INDICATES THE DIRECTION OF ROTATION FOR ANALYTIC CONTINUATION. |
|
|
10 C NZ=-1 MEANS AN OVERFLOW WILL OCCUR |
|
|
11 C |
|
|
12 C***ROUTINES CALLED ZKSCL,ZS1S2,ZUCHK,ZUNIK,D1MACH,XZABS |
|
|
13 C***END PROLOGUE ZUNK1 |
|
|
14 C COMPLEX CFN,CK,CONE,CRSC,CS,CSCL,CSGN,CSPN,CSR,CSS,CWRK,CY,CZERO, |
|
|
15 C *C1,C2,PHI,PHID,RZ,SUM,SUMD,S1,S2,Y,Z,ZETA1,ZETA1D,ZETA2,ZETA2D,ZR |
|
|
16 DOUBLE PRECISION ALIM, ANG, APHI, ASC, ASCLE, BRY, CKI, CKR, |
|
|
17 * CONER, CRSC, CSCL, CSGNI, CSPNI, CSPNR, CSR, CSRR, CSSR, |
|
|
18 * CWRKI, CWRKR, CYI, CYR, C1I, C1R, C2I, C2M, C2R, ELIM, FMR, FN, |
|
|
19 * FNF, FNU, PHIDI, PHIDR, PHII, PHIR, PI, RAST, RAZR, RS1, RZI, |
|
|
20 * RZR, SGN, STI, STR, SUMDI, SUMDR, SUMI, SUMR, S1I, S1R, S2I, |
|
|
21 * S2R, TOL, YI, YR, ZEROI, ZEROR, ZETA1I, ZETA1R, ZETA2I, ZETA2R, |
|
|
22 * ZET1DI, ZET1DR, ZET2DI, ZET2DR, ZI, ZR, ZRI, ZRR, D1MACH, XZABS |
|
|
23 INTEGER I, IB, IFLAG, IFN, IL, INIT, INU, IUF, K, KDFLG, KFLAG, |
|
|
24 * KK, KODE, MR, N, NW, NZ, INITD, IC, IPARD, J |
|
|
25 DIMENSION BRY(3), INIT(2), YR(N), YI(N), SUMR(2), SUMI(2), |
|
|
26 * ZETA1R(2), ZETA1I(2), ZETA2R(2), ZETA2I(2), CYR(2), CYI(2), |
|
|
27 * CWRKR(16,3), CWRKI(16,3), CSSR(3), CSRR(3), PHIR(2), PHII(2) |
|
|
28 DATA ZEROR,ZEROI,CONER / 0.0D0, 0.0D0, 1.0D0 / |
|
|
29 DATA PI / 3.14159265358979324D0 / |
|
|
30 C |
|
|
31 KDFLG = 1 |
|
|
32 NZ = 0 |
|
|
33 C----------------------------------------------------------------------- |
|
|
34 C EXP(-ALIM)=EXP(-ELIM)/TOL=APPROX. ONE PRECISION GREATER THAN |
|
|
35 C THE UNDERFLOW LIMIT |
|
|
36 C----------------------------------------------------------------------- |
|
|
37 CSCL = 1.0D0/TOL |
|
|
38 CRSC = TOL |
|
|
39 CSSR(1) = CSCL |
|
|
40 CSSR(2) = CONER |
|
|
41 CSSR(3) = CRSC |
|
|
42 CSRR(1) = CRSC |
|
|
43 CSRR(2) = CONER |
|
|
44 CSRR(3) = CSCL |
|
|
45 BRY(1) = 1.0D+3*D1MACH(1)/TOL |
|
|
46 BRY(2) = 1.0D0/BRY(1) |
|
|
47 BRY(3) = D1MACH(2) |
|
|
48 ZRR = ZR |
|
|
49 ZRI = ZI |
|
|
50 IF (ZR.GE.0.0D0) GO TO 10 |
|
|
51 ZRR = -ZR |
|
|
52 ZRI = -ZI |
|
|
53 10 CONTINUE |
|
|
54 J = 2 |
|
|
55 DO 70 I=1,N |
|
|
56 C----------------------------------------------------------------------- |
|
|
57 C J FLIP FLOPS BETWEEN 1 AND 2 IN J = 3 - J |
|
|
58 C----------------------------------------------------------------------- |
|
|
59 J = 3 - J |
|
|
60 FN = FNU + DBLE(FLOAT(I-1)) |
|
|
61 INIT(J) = 0 |
|
|
62 CALL ZUNIK(ZRR, ZRI, FN, 2, 0, TOL, INIT(J), PHIR(J), PHII(J), |
|
|
63 * ZETA1R(J), ZETA1I(J), ZETA2R(J), ZETA2I(J), SUMR(J), SUMI(J), |
|
|
64 * CWRKR(1,J), CWRKI(1,J)) |
|
|
65 IF (KODE.EQ.1) GO TO 20 |
|
|
66 STR = ZRR + ZETA2R(J) |
|
|
67 STI = ZRI + ZETA2I(J) |
|
|
68 RAST = FN/XZABS(STR,STI) |
|
|
69 STR = STR*RAST*RAST |
|
|
70 STI = -STI*RAST*RAST |
|
|
71 S1R = ZETA1R(J) - STR |
|
|
72 S1I = ZETA1I(J) - STI |
|
|
73 GO TO 30 |
|
|
74 20 CONTINUE |
|
|
75 S1R = ZETA1R(J) - ZETA2R(J) |
|
|
76 S1I = ZETA1I(J) - ZETA2I(J) |
|
|
77 30 CONTINUE |
|
|
78 RS1 = S1R |
|
|
79 C----------------------------------------------------------------------- |
|
|
80 C TEST FOR UNDERFLOW AND OVERFLOW |
|
|
81 C----------------------------------------------------------------------- |
|
|
82 IF (DABS(RS1).GT.ELIM) GO TO 60 |
|
|
83 IF (KDFLG.EQ.1) KFLAG = 2 |
|
|
84 IF (DABS(RS1).LT.ALIM) GO TO 40 |
|
|
85 C----------------------------------------------------------------------- |
|
|
86 C REFINE TEST AND SCALE |
|
|
87 C----------------------------------------------------------------------- |
|
|
88 APHI = XZABS(PHIR(J),PHII(J)) |
|
|
89 RS1 = RS1 + DLOG(APHI) |
|
|
90 IF (DABS(RS1).GT.ELIM) GO TO 60 |
|
|
91 IF (KDFLG.EQ.1) KFLAG = 1 |
|
|
92 IF (RS1.LT.0.0D0) GO TO 40 |
|
|
93 IF (KDFLG.EQ.1) KFLAG = 3 |
|
|
94 40 CONTINUE |
|
|
95 C----------------------------------------------------------------------- |
|
|
96 C SCALE S1 TO KEEP INTERMEDIATE ARITHMETIC ON SCALE NEAR |
|
|
97 C EXPONENT EXTREMES |
|
|
98 C----------------------------------------------------------------------- |
|
|
99 S2R = PHIR(J)*SUMR(J) - PHII(J)*SUMI(J) |
|
|
100 S2I = PHIR(J)*SUMI(J) + PHII(J)*SUMR(J) |
|
|
101 STR = DEXP(S1R)*CSSR(KFLAG) |
|
|
102 S1R = STR*DCOS(S1I) |
|
|
103 S1I = STR*DSIN(S1I) |
|
|
104 STR = S2R*S1R - S2I*S1I |
|
|
105 S2I = S1R*S2I + S2R*S1I |
|
|
106 S2R = STR |
|
|
107 IF (KFLAG.NE.1) GO TO 50 |
|
|
108 CALL ZUCHK(S2R, S2I, NW, BRY(1), TOL) |
|
|
109 IF (NW.NE.0) GO TO 60 |
|
|
110 50 CONTINUE |
|
|
111 CYR(KDFLG) = S2R |
|
|
112 CYI(KDFLG) = S2I |
|
|
113 YR(I) = S2R*CSRR(KFLAG) |
|
|
114 YI(I) = S2I*CSRR(KFLAG) |
|
|
115 IF (KDFLG.EQ.2) GO TO 75 |
|
|
116 KDFLG = 2 |
|
|
117 GO TO 70 |
|
|
118 60 CONTINUE |
|
|
119 IF (RS1.GT.0.0D0) GO TO 300 |
|
|
120 C----------------------------------------------------------------------- |
|
|
121 C FOR ZR.LT.0.0, THE I FUNCTION TO BE ADDED WILL OVERFLOW |
|
|
122 C----------------------------------------------------------------------- |
|
|
123 IF (ZR.LT.0.0D0) GO TO 300 |
|
|
124 KDFLG = 1 |
|
|
125 YR(I)=ZEROR |
|
|
126 YI(I)=ZEROI |
|
|
127 NZ=NZ+1 |
|
|
128 IF (I.EQ.1) GO TO 70 |
|
|
129 IF ((YR(I-1).EQ.ZEROR).AND.(YI(I-1).EQ.ZEROI)) GO TO 70 |
|
|
130 YR(I-1)=ZEROR |
|
|
131 YI(I-1)=ZEROI |
|
|
132 NZ=NZ+1 |
|
|
133 70 CONTINUE |
|
|
134 I = N |
|
|
135 75 CONTINUE |
|
|
136 RAZR = 1.0D0/XZABS(ZRR,ZRI) |
|
|
137 STR = ZRR*RAZR |
|
|
138 STI = -ZRI*RAZR |
|
|
139 RZR = (STR+STR)*RAZR |
|
|
140 RZI = (STI+STI)*RAZR |
|
|
141 CKR = FN*RZR |
|
|
142 CKI = FN*RZI |
|
|
143 IB = I + 1 |
|
|
144 IF (N.LT.IB) GO TO 160 |
|
|
145 C----------------------------------------------------------------------- |
|
|
146 C TEST LAST MEMBER FOR UNDERFLOW AND OVERFLOW. SET SEQUENCE TO ZERO |
|
|
147 C ON UNDERFLOW. |
|
|
148 C----------------------------------------------------------------------- |
|
|
149 FN = FNU + DBLE(FLOAT(N-1)) |
|
|
150 IPARD = 1 |
|
|
151 IF (MR.NE.0) IPARD = 0 |
|
|
152 INITD = 0 |
|
|
153 CALL ZUNIK(ZRR, ZRI, FN, 2, IPARD, TOL, INITD, PHIDR, PHIDI, |
|
|
154 * ZET1DR, ZET1DI, ZET2DR, ZET2DI, SUMDR, SUMDI, CWRKR(1,3), |
|
|
155 * CWRKI(1,3)) |
|
|
156 IF (KODE.EQ.1) GO TO 80 |
|
|
157 STR = ZRR + ZET2DR |
|
|
158 STI = ZRI + ZET2DI |
|
|
159 RAST = FN/XZABS(STR,STI) |
|
|
160 STR = STR*RAST*RAST |
|
|
161 STI = -STI*RAST*RAST |
|
|
162 S1R = ZET1DR - STR |
|
|
163 S1I = ZET1DI - STI |
|
|
164 GO TO 90 |
|
|
165 80 CONTINUE |
|
|
166 S1R = ZET1DR - ZET2DR |
|
|
167 S1I = ZET1DI - ZET2DI |
|
|
168 90 CONTINUE |
|
|
169 RS1 = S1R |
|
|
170 IF (DABS(RS1).GT.ELIM) GO TO 95 |
|
|
171 IF (DABS(RS1).LT.ALIM) GO TO 100 |
|
|
172 C---------------------------------------------------------------------------- |
|
|
173 C REFINE ESTIMATE AND TEST |
|
|
174 C------------------------------------------------------------------------- |
|
|
175 APHI = XZABS(PHIDR,PHIDI) |
|
|
176 RS1 = RS1+DLOG(APHI) |
|
|
177 IF (DABS(RS1).LT.ELIM) GO TO 100 |
|
|
178 95 CONTINUE |
|
|
179 IF (DABS(RS1).GT.0.0D0) GO TO 300 |
|
|
180 C----------------------------------------------------------------------- |
|
|
181 C FOR ZR.LT.0.0, THE I FUNCTION TO BE ADDED WILL OVERFLOW |
|
|
182 C----------------------------------------------------------------------- |
|
|
183 IF (ZR.LT.0.0D0) GO TO 300 |
|
|
184 NZ = N |
|
|
185 DO 96 I=1,N |
|
|
186 YR(I) = ZEROR |
|
|
187 YI(I) = ZEROI |
|
|
188 96 CONTINUE |
|
|
189 RETURN |
|
|
190 C--------------------------------------------------------------------------- |
|
|
191 C FORWARD RECUR FOR REMAINDER OF THE SEQUENCE |
|
|
192 C---------------------------------------------------------------------------- |
|
|
193 100 CONTINUE |
|
|
194 S1R = CYR(1) |
|
|
195 S1I = CYI(1) |
|
|
196 S2R = CYR(2) |
|
|
197 S2I = CYI(2) |
|
|
198 C1R = CSRR(KFLAG) |
|
|
199 ASCLE = BRY(KFLAG) |
|
|
200 DO 120 I=IB,N |
|
|
201 C2R = S2R |
|
|
202 C2I = S2I |
|
|
203 S2R = CKR*C2R - CKI*C2I + S1R |
|
|
204 S2I = CKR*C2I + CKI*C2R + S1I |
|
|
205 S1R = C2R |
|
|
206 S1I = C2I |
|
|
207 CKR = CKR + RZR |
|
|
208 CKI = CKI + RZI |
|
|
209 C2R = S2R*C1R |
|
|
210 C2I = S2I*C1R |
|
|
211 YR(I) = C2R |
|
|
212 YI(I) = C2I |
|
|
213 IF (KFLAG.GE.3) GO TO 120 |
|
|
214 STR = DABS(C2R) |
|
|
215 STI = DABS(C2I) |
|
|
216 C2M = DMAX1(STR,STI) |
|
|
217 IF (C2M.LE.ASCLE) GO TO 120 |
|
|
218 KFLAG = KFLAG + 1 |
|
|
219 ASCLE = BRY(KFLAG) |
|
|
220 S1R = S1R*C1R |
|
|
221 S1I = S1I*C1R |
|
|
222 S2R = C2R |
|
|
223 S2I = C2I |
|
|
224 S1R = S1R*CSSR(KFLAG) |
|
|
225 S1I = S1I*CSSR(KFLAG) |
|
|
226 S2R = S2R*CSSR(KFLAG) |
|
|
227 S2I = S2I*CSSR(KFLAG) |
|
|
228 C1R = CSRR(KFLAG) |
|
|
229 120 CONTINUE |
|
|
230 160 CONTINUE |
|
|
231 IF (MR.EQ.0) RETURN |
|
|
232 C----------------------------------------------------------------------- |
|
|
233 C ANALYTIC CONTINUATION FOR RE(Z).LT.0.0D0 |
|
|
234 C----------------------------------------------------------------------- |
|
|
235 NZ = 0 |
|
|
236 FMR = DBLE(FLOAT(MR)) |
|
|
237 SGN = -DSIGN(PI,FMR) |
|
|
238 C----------------------------------------------------------------------- |
|
|
239 C CSPN AND CSGN ARE COEFF OF K AND I FUNCTIONS RESP. |
|
|
240 C----------------------------------------------------------------------- |
|
|
241 CSGNI = SGN |
|
|
242 INU = INT(SNGL(FNU)) |
|
|
243 FNF = FNU - DBLE(FLOAT(INU)) |
|
|
244 IFN = INU + N - 1 |
|
|
245 ANG = FNF*SGN |
|
|
246 CSPNR = DCOS(ANG) |
|
|
247 CSPNI = DSIN(ANG) |
|
|
248 IF (MOD(IFN,2).EQ.0) GO TO 170 |
|
|
249 CSPNR = -CSPNR |
|
|
250 CSPNI = -CSPNI |
|
|
251 170 CONTINUE |
|
|
252 ASC = BRY(1) |
|
|
253 IUF = 0 |
|
|
254 KK = N |
|
|
255 KDFLG = 1 |
|
|
256 IB = IB - 1 |
|
|
257 IC = IB - 1 |
|
|
258 DO 270 K=1,N |
|
|
259 FN = FNU + DBLE(FLOAT(KK-1)) |
|
|
260 C----------------------------------------------------------------------- |
|
|
261 C LOGIC TO SORT OUT CASES WHOSE PARAMETERS WERE SET FOR THE K |
|
|
262 C FUNCTION ABOVE |
|
|
263 C----------------------------------------------------------------------- |
|
|
264 M=3 |
|
|
265 IF (N.GT.2) GO TO 175 |
|
|
266 172 CONTINUE |
|
|
267 INITD = INIT(J) |
|
|
268 PHIDR = PHIR(J) |
|
|
269 PHIDI = PHII(J) |
|
|
270 ZET1DR = ZETA1R(J) |
|
|
271 ZET1DI = ZETA1I(J) |
|
|
272 ZET2DR = ZETA2R(J) |
|
|
273 ZET2DI = ZETA2I(J) |
|
|
274 SUMDR = SUMR(J) |
|
|
275 SUMDI = SUMI(J) |
|
|
276 M = J |
|
|
277 J = 3 - J |
|
|
278 GO TO 180 |
|
|
279 175 CONTINUE |
|
|
280 IF ((KK.EQ.N).AND.(IB.LT.N)) GO TO 180 |
|
|
281 IF ((KK.EQ.IB).OR.(KK.EQ.IC)) GO TO 172 |
|
|
282 INITD = 0 |
|
|
283 180 CONTINUE |
|
|
284 CALL ZUNIK(ZRR, ZRI, FN, 1, 0, TOL, INITD, PHIDR, PHIDI, |
|
|
285 * ZET1DR, ZET1DI, ZET2DR, ZET2DI, SUMDR, SUMDI, |
|
|
286 * CWRKR(1,M), CWRKI(1,M)) |
|
|
287 IF (KODE.EQ.1) GO TO 200 |
|
|
288 STR = ZRR + ZET2DR |
|
|
289 STI = ZRI + ZET2DI |
|
|
290 RAST = FN/XZABS(STR,STI) |
|
|
291 STR = STR*RAST*RAST |
|
|
292 STI = -STI*RAST*RAST |
|
|
293 S1R = -ZET1DR + STR |
|
|
294 S1I = -ZET1DI + STI |
|
|
295 GO TO 210 |
|
|
296 200 CONTINUE |
|
|
297 S1R = -ZET1DR + ZET2DR |
|
|
298 S1I = -ZET1DI + ZET2DI |
|
|
299 210 CONTINUE |
|
|
300 C----------------------------------------------------------------------- |
|
|
301 C TEST FOR UNDERFLOW AND OVERFLOW |
|
|
302 C----------------------------------------------------------------------- |
|
|
303 RS1 = S1R |
|
|
304 IF (DABS(RS1).GT.ELIM) GO TO 260 |
|
|
305 IF (KDFLG.EQ.1) IFLAG = 2 |
|
|
306 IF (DABS(RS1).LT.ALIM) GO TO 220 |
|
|
307 C----------------------------------------------------------------------- |
|
|
308 C REFINE TEST AND SCALE |
|
|
309 C----------------------------------------------------------------------- |
|
|
310 APHI = XZABS(PHIDR,PHIDI) |
|
|
311 RS1 = RS1 + DLOG(APHI) |
|
|
312 IF (DABS(RS1).GT.ELIM) GO TO 260 |
|
|
313 IF (KDFLG.EQ.1) IFLAG = 1 |
|
|
314 IF (RS1.LT.0.0D0) GO TO 220 |
|
|
315 IF (KDFLG.EQ.1) IFLAG = 3 |
|
|
316 220 CONTINUE |
|
|
317 STR = PHIDR*SUMDR - PHIDI*SUMDI |
|
|
318 STI = PHIDR*SUMDI + PHIDI*SUMDR |
|
|
319 S2R = -CSGNI*STI |
|
|
320 S2I = CSGNI*STR |
|
|
321 STR = DEXP(S1R)*CSSR(IFLAG) |
|
|
322 S1R = STR*DCOS(S1I) |
|
|
323 S1I = STR*DSIN(S1I) |
|
|
324 STR = S2R*S1R - S2I*S1I |
|
|
325 S2I = S2R*S1I + S2I*S1R |
|
|
326 S2R = STR |
|
|
327 IF (IFLAG.NE.1) GO TO 230 |
|
|
328 CALL ZUCHK(S2R, S2I, NW, BRY(1), TOL) |
|
|
329 IF (NW.EQ.0) GO TO 230 |
|
|
330 S2R = ZEROR |
|
|
331 S2I = ZEROI |
|
|
332 230 CONTINUE |
|
|
333 CYR(KDFLG) = S2R |
|
|
334 CYI(KDFLG) = S2I |
|
|
335 C2R = S2R |
|
|
336 C2I = S2I |
|
|
337 S2R = S2R*CSRR(IFLAG) |
|
|
338 S2I = S2I*CSRR(IFLAG) |
|
|
339 C----------------------------------------------------------------------- |
|
|
340 C ADD I AND K FUNCTIONS, K SEQUENCE IN Y(I), I=1,N |
|
|
341 C----------------------------------------------------------------------- |
|
|
342 S1R = YR(KK) |
|
|
343 S1I = YI(KK) |
|
|
344 IF (KODE.EQ.1) GO TO 250 |
|
|
345 CALL ZS1S2(ZRR, ZRI, S1R, S1I, S2R, S2I, NW, ASC, ALIM, IUF) |
|
|
346 NZ = NZ + NW |
|
|
347 250 CONTINUE |
|
|
348 YR(KK) = S1R*CSPNR - S1I*CSPNI + S2R |
|
|
349 YI(KK) = CSPNR*S1I + CSPNI*S1R + S2I |
|
|
350 KK = KK - 1 |
|
|
351 CSPNR = -CSPNR |
|
|
352 CSPNI = -CSPNI |
|
|
353 IF (C2R.NE.0.0D0 .OR. C2I.NE.0.0D0) GO TO 255 |
|
|
354 KDFLG = 1 |
|
|
355 GO TO 270 |
|
|
356 255 CONTINUE |
|
|
357 IF (KDFLG.EQ.2) GO TO 275 |
|
|
358 KDFLG = 2 |
|
|
359 GO TO 270 |
|
|
360 260 CONTINUE |
|
|
361 IF (RS1.GT.0.0D0) GO TO 300 |
|
|
362 S2R = ZEROR |
|
|
363 S2I = ZEROI |
|
|
364 GO TO 230 |
|
|
365 270 CONTINUE |
|
|
366 K = N |
|
|
367 275 CONTINUE |
|
|
368 IL = N - K |
|
|
369 IF (IL.EQ.0) RETURN |
|
|
370 C----------------------------------------------------------------------- |
|
|
371 C RECUR BACKWARD FOR REMAINDER OF I SEQUENCE AND ADD IN THE |
|
|
372 C K FUNCTIONS, SCALING THE I SEQUENCE DURING RECURRENCE TO KEEP |
|
|
373 C INTERMEDIATE ARITHMETIC ON SCALE NEAR EXPONENT EXTREMES. |
|
|
374 C----------------------------------------------------------------------- |
|
|
375 S1R = CYR(1) |
|
|
376 S1I = CYI(1) |
|
|
377 S2R = CYR(2) |
|
|
378 S2I = CYI(2) |
|
|
379 CSR = CSRR(IFLAG) |
|
|
380 ASCLE = BRY(IFLAG) |
|
|
381 FN = DBLE(FLOAT(INU+IL)) |
|
|
382 DO 290 I=1,IL |
|
|
383 C2R = S2R |
|
|
384 C2I = S2I |
|
|
385 S2R = S1R + (FN+FNF)*(RZR*C2R-RZI*C2I) |
|
|
386 S2I = S1I + (FN+FNF)*(RZR*C2I+RZI*C2R) |
|
|
387 S1R = C2R |
|
|
388 S1I = C2I |
|
|
389 FN = FN - 1.0D0 |
|
|
390 C2R = S2R*CSR |
|
|
391 C2I = S2I*CSR |
|
|
392 CKR = C2R |
|
|
393 CKI = C2I |
|
|
394 C1R = YR(KK) |
|
|
395 C1I = YI(KK) |
|
|
396 IF (KODE.EQ.1) GO TO 280 |
|
|
397 CALL ZS1S2(ZRR, ZRI, C1R, C1I, C2R, C2I, NW, ASC, ALIM, IUF) |
|
|
398 NZ = NZ + NW |
|
|
399 280 CONTINUE |
|
|
400 YR(KK) = C1R*CSPNR - C1I*CSPNI + C2R |
|
|
401 YI(KK) = C1R*CSPNI + C1I*CSPNR + C2I |
|
|
402 KK = KK - 1 |
|
|
403 CSPNR = -CSPNR |
|
|
404 CSPNI = -CSPNI |
|
|
405 IF (IFLAG.GE.3) GO TO 290 |
|
|
406 C2R = DABS(CKR) |
|
|
407 C2I = DABS(CKI) |
|
|
408 C2M = DMAX1(C2R,C2I) |
|
|
409 IF (C2M.LE.ASCLE) GO TO 290 |
|
|
410 IFLAG = IFLAG + 1 |
|
|
411 ASCLE = BRY(IFLAG) |
|
|
412 S1R = S1R*CSR |
|
|
413 S1I = S1I*CSR |
|
|
414 S2R = CKR |
|
|
415 S2I = CKI |
|
|
416 S1R = S1R*CSSR(IFLAG) |
|
|
417 S1I = S1I*CSSR(IFLAG) |
|
|
418 S2R = S2R*CSSR(IFLAG) |
|
|
419 S2I = S2I*CSSR(IFLAG) |
|
|
420 CSR = CSRR(IFLAG) |
|
|
421 290 CONTINUE |
|
|
422 RETURN |
|
|
423 300 CONTINUE |
|
|
424 NZ = -1 |
|
|
425 RETURN |
|
|
426 END |
