GNU Octave  9.1.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
zunk2.f
Go to the documentation of this file.
1  SUBROUTINE zunk2(ZR, ZI, FNU, KODE, MR, N, YR, YI, NZ, TOL, ELIM,
2  * ALIM)
3 C***BEGIN PROLOGUE ZUNK2
4 C***REFER TO ZBESK
5 C
6 C ZUNK2 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 EXPANSIONS FOR H(KIND,FNU,ZN) AND J(FNU,ZN)
9 C WHERE ZN IS IN THE RIGHT HALF PLANE, KIND=(3-MR)/2, MR=+1 OR
10 C -1. HERE ZN=ZR*I OR -ZR*I WHERE ZR=Z IF Z IS IN THE RIGHT
11 C HALF PLANE OR ZR=-Z IF Z IS IN THE LEFT HALF PLANE. MR INDIC-
12 C ATES THE DIRECTION OF ROTATION FOR ANALYTIC CONTINUATION.
13 C NZ=-1 MEANS AN OVERFLOW WILL OCCUR
14 C
15 C***ROUTINES CALLED ZAIRY,ZKSCL,ZS1S2,ZUCHK,ZUNHJ,D1MACH,XZABS
16 C***END PROLOGUE ZUNK2
17 C COMPLEX AI,ARG,ARGD,ASUM,ASUMD,BSUM,BSUMD,CFN,CI,CIP,CK,CONE,CRSC,
18 C *CR1,CR2,CS,CSCL,CSGN,CSPN,CSR,CSS,CY,CZERO,C1,C2,DAI,PHI,PHID,RZ,
19 C *S1,S2,Y,Z,ZB,ZETA1,ZETA1D,ZETA2,ZETA2D,ZN,ZR
20  DOUBLE PRECISION AARG, AIC, AII, AIR, ALIM, ANG, APHI, ARGDI,
21  * argdr, argi, argr, asc, ascle, asumdi, asumdr, asumi, asumr,
22  * bry, bsumdi, bsumdr, bsumi, bsumr, car, cipi, cipr, cki, ckr,
23  * coner, crsc, cr1i, cr1r, cr2i, cr2r, cscl, csgni, csi,
24  * cspni, cspnr, csr, csrr, cssr, cyi, cyr, c1i, c1r, c2i, c2m,
25  * c2r, daii, dair, elim, fmr, fn, fnf, fnu, hpi, phidi, phidr,
26  * phii, phir, pi, pti, ptr, rast, razr, rs1, rzi, rzr, sar, sgn,
27  * sti, str, s1i, s1r, s2i, s2r, tol, yi, yr, yy, zbi, zbr, zeroi,
28  * zeror, zeta1i, zeta1r, zeta2i, zeta2r, zet1di, zet1dr, zet2di,
29  * zet2dr, zi, zni, znr, zr, zri, zrr, d1mach, xzabs
30  INTEGER I, IB, IFLAG, IFN, IL, IN, INU, IUF, K, KDFLG, KFLAG, KK,
31  * kode, mr, n, nai, ndai, nw, nz, idum, j, ipard, ic
32  dimension bry(3), yr(n), yi(n), asumr(2), asumi(2), bsumr(2),
33  * bsumi(2), phir(2), phii(2), argr(2), argi(2), zeta1r(2),
34  * zeta1i(2), zeta2r(2), zeta2i(2), cyr(2), cyi(2), cipr(4),
35  * cipi(4), cssr(3), csrr(3)
36  DATA zeror,zeroi,coner,cr1r,cr1i,cr2r,cr2i /
37  1 0.0d0, 0.0d0, 1.0d0,
38  1 1.0d0,1.73205080756887729d0 , -0.5d0,-8.66025403784438647d-01 /
39  DATA hpi, pi, aic /
40  1 1.57079632679489662d+00, 3.14159265358979324d+00,
41  1 1.26551212348464539d+00/
42  DATA cipr(1),cipi(1),cipr(2),cipi(2),cipr(3),cipi(3),cipr(4),
43  * cipi(4) /
44  1 1.0d0,0.0d0 , 0.0d0,-1.0d0 , -1.0d0,0.0d0 , 0.0d0,1.0d0 /
45 C
46  kdflg = 1
47  nz = 0
48 C-----------------------------------------------------------------------
49 C EXP(-ALIM)=EXP(-ELIM)/TOL=APPROX. ONE PRECISION GREATER THAN
50 C THE UNDERFLOW LIMIT
51 C-----------------------------------------------------------------------
52  cscl = 1.0d0/tol
53  crsc = tol
54  cssr(1) = cscl
55  cssr(2) = coner
56  cssr(3) = crsc
57  csrr(1) = crsc
58  csrr(2) = coner
59  csrr(3) = cscl
60  bry(1) = 1.0d+3*d1mach(1)/tol
61  bry(2) = 1.0d0/bry(1)
62  bry(3) = d1mach(2)
63  zrr = zr
64  zri = zi
65  IF (zr.GE.0.0d0) GO TO 10
66  zrr = -zr
67  zri = -zi
68  10 CONTINUE
69  yy = zri
70  znr = zri
71  zni = -zrr
72  zbr = zrr
73  zbi = zri
74  inu = int(sngl(fnu))
75  fnf = fnu - dble(float(inu))
76  ang = -hpi*fnf
77  car = dcos(ang)
78  sar = dsin(ang)
79  c2r = hpi*sar
80  c2i = -hpi*car
81  kk = mod(inu,4) + 1
82  str = c2r*cipr(kk) - c2i*cipi(kk)
83  sti = c2r*cipi(kk) + c2i*cipr(kk)
84  csr = cr1r*str - cr1i*sti
85  csi = cr1r*sti + cr1i*str
86  IF (yy.GT.0.0d0) GO TO 20
87  znr = -znr
88  zbi = -zbi
89  20 CONTINUE
90 C-----------------------------------------------------------------------
91 C K(FNU,Z) IS COMPUTED FROM H(2,FNU,-I*Z) WHERE Z IS IN THE FIRST
92 C QUADRANT. FOURTH QUADRANT VALUES (YY.LE.0.0E0) ARE COMPUTED BY
93 C CONJUGATION SINCE THE K FUNCTION IS REAL ON THE POSITIVE REAL AXIS
94 C-----------------------------------------------------------------------
95  j = 2
96  DO 80 i=1,n
97 C-----------------------------------------------------------------------
98 C J FLIP FLOPS BETWEEN 1 AND 2 IN J = 3 - J
99 C-----------------------------------------------------------------------
100  j = 3 - j
101  fn = fnu + dble(float(i-1))
102  CALL zunhj(znr, zni, fn, 0, tol, phir(j), phii(j), argr(j),
103  * argi(j), zeta1r(j), zeta1i(j), zeta2r(j), zeta2i(j), asumr(j),
104  * asumi(j), bsumr(j), bsumi(j))
105  IF (kode.EQ.1) GO TO 30
106  str = zbr + zeta2r(j)
107  sti = zbi + zeta2i(j)
108  rast = fn/xzabs(str,sti)
109  str = str*rast*rast
110  sti = -sti*rast*rast
111  s1r = zeta1r(j) - str
112  s1i = zeta1i(j) - sti
113  GO TO 40
114  30 CONTINUE
115  s1r = zeta1r(j) - zeta2r(j)
116  s1i = zeta1i(j) - zeta2i(j)
117  40 CONTINUE
118 C-----------------------------------------------------------------------
119 C TEST FOR UNDERFLOW AND OVERFLOW
120 C-----------------------------------------------------------------------
121  rs1 = s1r
122  IF (dabs(rs1).GT.elim) GO TO 70
123  IF (kdflg.EQ.1) kflag = 2
124  IF (dabs(rs1).LT.alim) GO TO 50
125 C-----------------------------------------------------------------------
126 C REFINE TEST AND SCALE
127 C-----------------------------------------------------------------------
128  aphi = xzabs(phir(j),phii(j))
129  aarg = xzabs(argr(j),argi(j))
130  rs1 = rs1 + dlog(aphi) - 0.25d0*dlog(aarg) - aic
131  IF (dabs(rs1).GT.elim) GO TO 70
132  IF (kdflg.EQ.1) kflag = 1
133  IF (rs1.LT.0.0d0) GO TO 50
134  IF (kdflg.EQ.1) kflag = 3
135  50 CONTINUE
136 C-----------------------------------------------------------------------
137 C SCALE S1 TO KEEP INTERMEDIATE ARITHMETIC ON SCALE NEAR
138 C EXPONENT EXTREMES
139 C-----------------------------------------------------------------------
140  c2r = argr(j)*cr2r - argi(j)*cr2i
141  c2i = argr(j)*cr2i + argi(j)*cr2r
142  CALL zairy(c2r, c2i, 0, 2, air, aii, nai, idum)
143  CALL zairy(c2r, c2i, 1, 2, dair, daii, ndai, idum)
144  str = dair*bsumr(j) - daii*bsumi(j)
145  sti = dair*bsumi(j) + daii*bsumr(j)
146  ptr = str*cr2r - sti*cr2i
147  pti = str*cr2i + sti*cr2r
148  str = ptr + (air*asumr(j)-aii*asumi(j))
149  sti = pti + (air*asumi(j)+aii*asumr(j))
150  ptr = str*phir(j) - sti*phii(j)
151  pti = str*phii(j) + sti*phir(j)
152  s2r = ptr*csr - pti*csi
153  s2i = ptr*csi + pti*csr
154  str = dexp(s1r)*cssr(kflag)
155  s1r = str*dcos(s1i)
156  s1i = str*dsin(s1i)
157  str = s2r*s1r - s2i*s1i
158  s2i = s1r*s2i + s2r*s1i
159  s2r = str
160  IF (kflag.NE.1) GO TO 60
161  CALL zuchk(s2r, s2i, nw, bry(1), tol)
162  IF (nw.NE.0) GO TO 70
163  60 CONTINUE
164  IF (yy.LE.0.0d0) s2i = -s2i
165  cyr(kdflg) = s2r
166  cyi(kdflg) = s2i
167  yr(i) = s2r*csrr(kflag)
168  yi(i) = s2i*csrr(kflag)
169  str = csi
170  csi = -csr
171  csr = str
172  IF (kdflg.EQ.2) GO TO 85
173  kdflg = 2
174  GO TO 80
175  70 CONTINUE
176  IF (rs1.GT.0.0d0) GO TO 320
177 C-----------------------------------------------------------------------
178 C FOR ZR.LT.0.0, THE I FUNCTION TO BE ADDED WILL OVERFLOW
179 C-----------------------------------------------------------------------
180  IF (zr.LT.0.0d0) GO TO 320
181  kdflg = 1
182  yr(i)=zeror
183  yi(i)=zeroi
184  nz=nz+1
185  str = csi
186  csi =-csr
187  csr = str
188  IF (i.EQ.1) GO TO 80
189  IF ((yr(i-1).EQ.zeror).AND.(yi(i-1).EQ.zeroi)) GO TO 80
190  yr(i-1)=zeror
191  yi(i-1)=zeroi
192  nz=nz+1
193  80 CONTINUE
194  i = n
195  85 CONTINUE
196  razr = 1.0d0/xzabs(zrr,zri)
197  str = zrr*razr
198  sti = -zri*razr
199  rzr = (str+str)*razr
200  rzi = (sti+sti)*razr
201  ckr = fn*rzr
202  cki = fn*rzi
203  ib = i + 1
204  IF (n.LT.ib) GO TO 180
205 C-----------------------------------------------------------------------
206 C TEST LAST MEMBER FOR UNDERFLOW AND OVERFLOW. SET SEQUENCE TO ZERO
207 C ON UNDERFLOW.
208 C-----------------------------------------------------------------------
209  fn = fnu + dble(float(n-1))
210  ipard = 1
211  IF (mr.NE.0) ipard = 0
212  CALL zunhj(znr, zni, fn, ipard, tol, phidr, phidi, argdr, argdi,
213  * zet1dr, zet1di, zet2dr, zet2di, asumdr, asumdi, bsumdr, bsumdi)
214  IF (kode.EQ.1) GO TO 90
215  str = zbr + zet2dr
216  sti = zbi + zet2di
217  rast = fn/xzabs(str,sti)
218  str = str*rast*rast
219  sti = -sti*rast*rast
220  s1r = zet1dr - str
221  s1i = zet1di - sti
222  GO TO 100
223  90 CONTINUE
224  s1r = zet1dr - zet2dr
225  s1i = zet1di - zet2di
226  100 CONTINUE
227  rs1 = s1r
228  IF (dabs(rs1).GT.elim) GO TO 105
229  IF (dabs(rs1).LT.alim) GO TO 120
230 C----------------------------------------------------------------------------
231 C REFINE ESTIMATE AND TEST
232 C-------------------------------------------------------------------------
233  aphi = xzabs(phidr,phidi)
234  rs1 = rs1+dlog(aphi)
235  IF (dabs(rs1).LT.elim) GO TO 120
236  105 CONTINUE
237  IF (rs1.GT.0.0d0) GO TO 320
238 C-----------------------------------------------------------------------
239 C FOR ZR.LT.0.0, THE I FUNCTION TO BE ADDED WILL OVERFLOW
240 C-----------------------------------------------------------------------
241  IF (zr.LT.0.0d0) GO TO 320
242  nz = n
243  DO 106 i=1,n
244  yr(i) = zeror
245  yi(i) = zeroi
246  106 CONTINUE
247  RETURN
248  120 CONTINUE
249  s1r = cyr(1)
250  s1i = cyi(1)
251  s2r = cyr(2)
252  s2i = cyi(2)
253  c1r = csrr(kflag)
254  ascle = bry(kflag)
255  DO 130 i=ib,n
256  c2r = s2r
257  c2i = s2i
258  s2r = ckr*c2r - cki*c2i + s1r
259  s2i = ckr*c2i + cki*c2r + s1i
260  s1r = c2r
261  s1i = c2i
262  ckr = ckr + rzr
263  cki = cki + rzi
264  c2r = s2r*c1r
265  c2i = s2i*c1r
266  yr(i) = c2r
267  yi(i) = c2i
268  IF (kflag.GE.3) GO TO 130
269  str = dabs(c2r)
270  sti = dabs(c2i)
271  c2m = dmax1(str,sti)
272  IF (c2m.LE.ascle) GO TO 130
273  kflag = kflag + 1
274  ascle = bry(kflag)
275  s1r = s1r*c1r
276  s1i = s1i*c1r
277  s2r = c2r
278  s2i = c2i
279  s1r = s1r*cssr(kflag)
280  s1i = s1i*cssr(kflag)
281  s2r = s2r*cssr(kflag)
282  s2i = s2i*cssr(kflag)
283  c1r = csrr(kflag)
284  130 CONTINUE
285  180 CONTINUE
286  IF (mr.EQ.0) RETURN
287 C-----------------------------------------------------------------------
288 C ANALYTIC CONTINUATION FOR RE(Z).LT.0.0D0
289 C-----------------------------------------------------------------------
290  nz = 0
291  fmr = dble(float(mr))
292  sgn = -dsign(pi,fmr)
293 C-----------------------------------------------------------------------
294 C CSPN AND CSGN ARE COEFF OF K AND I FUNCIONS RESP.
295 C-----------------------------------------------------------------------
296  csgni = sgn
297  IF (yy.LE.0.0d0) csgni = -csgni
298  ifn = inu + n - 1
299  ang = fnf*sgn
300  cspnr = dcos(ang)
301  cspni = dsin(ang)
302  IF (mod(ifn,2).EQ.0) GO TO 190
303  cspnr = -cspnr
304  cspni = -cspni
305  190 CONTINUE
306 C-----------------------------------------------------------------------
307 C CS=COEFF OF THE J FUNCTION TO GET THE I FUNCTION. I(FNU,Z) IS
308 C COMPUTED FROM EXP(I*FNU*HPI)*J(FNU,-I*Z) WHERE Z IS IN THE FIRST
309 C QUADRANT. FOURTH QUADRANT VALUES (YY.LE.0.0E0) ARE COMPUTED BY
310 C CONJUGATION SINCE THE I FUNCTION IS REAL ON THE POSITIVE REAL AXIS
311 C-----------------------------------------------------------------------
312  csr = sar*csgni
313  csi = car*csgni
314  in = mod(ifn,4) + 1
315  c2r = cipr(in)
316  c2i = cipi(in)
317  str = csr*c2r + csi*c2i
318  csi = -csr*c2i + csi*c2r
319  csr = str
320  asc = bry(1)
321  iuf = 0
322  kk = n
323  kdflg = 1
324  ib = ib - 1
325  ic = ib - 1
326  DO 290 k=1,n
327  fn = fnu + dble(float(kk-1))
328 C-----------------------------------------------------------------------
329 C LOGIC TO SORT OUT CASES WHOSE PARAMETERS WERE SET FOR THE K
330 C FUNCTION ABOVE
331 C-----------------------------------------------------------------------
332  IF (n.GT.2) GO TO 175
333  172 CONTINUE
334  phidr = phir(j)
335  phidi = phii(j)
336  argdr = argr(j)
337  argdi = argi(j)
338  zet1dr = zeta1r(j)
339  zet1di = zeta1i(j)
340  zet2dr = zeta2r(j)
341  zet2di = zeta2i(j)
342  asumdr = asumr(j)
343  asumdi = asumi(j)
344  bsumdr = bsumr(j)
345  bsumdi = bsumi(j)
346  j = 3 - j
347  GO TO 210
348  175 CONTINUE
349  IF ((kk.EQ.n).AND.(ib.LT.n)) GO TO 210
350  IF ((kk.EQ.ib).OR.(kk.EQ.ic)) GO TO 172
351  CALL zunhj(znr, zni, fn, 0, tol, phidr, phidi, argdr,
352  * argdi, zet1dr, zet1di, zet2dr, zet2di, asumdr,
353  * asumdi, bsumdr, bsumdi)
354  210 CONTINUE
355  IF (kode.EQ.1) GO TO 220
356  str = zbr + zet2dr
357  sti = zbi + zet2di
358  rast = fn/xzabs(str,sti)
359  str = str*rast*rast
360  sti = -sti*rast*rast
361  s1r = -zet1dr + str
362  s1i = -zet1di + sti
363  GO TO 230
364  220 CONTINUE
365  s1r = -zet1dr + zet2dr
366  s1i = -zet1di + zet2di
367  230 CONTINUE
368 C-----------------------------------------------------------------------
369 C TEST FOR UNDERFLOW AND OVERFLOW
370 C-----------------------------------------------------------------------
371  rs1 = s1r
372  IF (dabs(rs1).GT.elim) GO TO 280
373  IF (kdflg.EQ.1) iflag = 2
374  IF (dabs(rs1).LT.alim) GO TO 240
375 C-----------------------------------------------------------------------
376 C REFINE TEST AND SCALE
377 C-----------------------------------------------------------------------
378  aphi = xzabs(phidr,phidi)
379  aarg = xzabs(argdr,argdi)
380  rs1 = rs1 + dlog(aphi) - 0.25d0*dlog(aarg) - aic
381  IF (dabs(rs1).GT.elim) GO TO 280
382  IF (kdflg.EQ.1) iflag = 1
383  IF (rs1.LT.0.0d0) GO TO 240
384  IF (kdflg.EQ.1) iflag = 3
385  240 CONTINUE
386  CALL zairy(argdr, argdi, 0, 2, air, aii, nai, idum)
387  CALL zairy(argdr, argdi, 1, 2, dair, daii, ndai, idum)
388  str = dair*bsumdr - daii*bsumdi
389  sti = dair*bsumdi + daii*bsumdr
390  str = str + (air*asumdr-aii*asumdi)
391  sti = sti + (air*asumdi+aii*asumdr)
392  ptr = str*phidr - sti*phidi
393  pti = str*phidi + sti*phidr
394  s2r = ptr*csr - pti*csi
395  s2i = ptr*csi + pti*csr
396  str = dexp(s1r)*cssr(iflag)
397  s1r = str*dcos(s1i)
398  s1i = str*dsin(s1i)
399  str = s2r*s1r - s2i*s1i
400  s2i = s2r*s1i + s2i*s1r
401  s2r = str
402  IF (iflag.NE.1) GO TO 250
403  CALL zuchk(s2r, s2i, nw, bry(1), tol)
404  IF (nw.EQ.0) GO TO 250
405  s2r = zeror
406  s2i = zeroi
407  250 CONTINUE
408  IF (yy.LE.0.0d0) s2i = -s2i
409  cyr(kdflg) = s2r
410  cyi(kdflg) = s2i
411  c2r = s2r
412  c2i = s2i
413  s2r = s2r*csrr(iflag)
414  s2i = s2i*csrr(iflag)
415 C-----------------------------------------------------------------------
416 C ADD I AND K FUNCTIONS, K SEQUENCE IN Y(I), I=1,N
417 C-----------------------------------------------------------------------
418  s1r = yr(kk)
419  s1i = yi(kk)
420  IF (kode.EQ.1) GO TO 270
421  CALL zs1s2(zrr, zri, s1r, s1i, s2r, s2i, nw, asc, alim, iuf)
422  nz = nz + nw
423  270 CONTINUE
424  yr(kk) = s1r*cspnr - s1i*cspni + s2r
425  yi(kk) = s1r*cspni + s1i*cspnr + s2i
426  kk = kk - 1
427  cspnr = -cspnr
428  cspni = -cspni
429  str = csi
430  csi = -csr
431  csr = str
432  IF (c2r.NE.0.0d0 .OR. c2i.NE.0.0d0) GO TO 255
433  kdflg = 1
434  GO TO 290
435  255 CONTINUE
436  IF (kdflg.EQ.2) GO TO 295
437  kdflg = 2
438  GO TO 290
439  280 CONTINUE
440  IF (rs1.GT.0.0d0) GO TO 320
441  s2r = zeror
442  s2i = zeroi
443  GO TO 250
444  290 CONTINUE
445  k = n
446  295 CONTINUE
447  il = n - k
448  IF (il.EQ.0) RETURN
449 C-----------------------------------------------------------------------
450 C RECUR BACKWARD FOR REMAINDER OF I SEQUENCE AND ADD IN THE
451 C K FUNCTIONS, SCALING THE I SEQUENCE DURING RECURRENCE TO KEEP
452 C INTERMEDIATE ARITHMETIC ON SCALE NEAR EXPONENT EXTREMES.
453 C-----------------------------------------------------------------------
454  s1r = cyr(1)
455  s1i = cyi(1)
456  s2r = cyr(2)
457  s2i = cyi(2)
458  csr = csrr(iflag)
459  ascle = bry(iflag)
460  fn = dble(float(inu+il))
461  DO 310 i=1,il
462  c2r = s2r
463  c2i = s2i
464  s2r = s1r + (fn+fnf)*(rzr*c2r-rzi*c2i)
465  s2i = s1i + (fn+fnf)*(rzr*c2i+rzi*c2r)
466  s1r = c2r
467  s1i = c2i
468  fn = fn - 1.0d0
469  c2r = s2r*csr
470  c2i = s2i*csr
471  ckr = c2r
472  cki = c2i
473  c1r = yr(kk)
474  c1i = yi(kk)
475  IF (kode.EQ.1) GO TO 300
476  CALL zs1s2(zrr, zri, c1r, c1i, c2r, c2i, nw, asc, alim, iuf)
477  nz = nz + nw
478  300 CONTINUE
479  yr(kk) = c1r*cspnr - c1i*cspni + c2r
480  yi(kk) = c1r*cspni + c1i*cspnr + c2i
481  kk = kk - 1
482  cspnr = -cspnr
483  cspni = -cspni
484  IF (iflag.GE.3) GO TO 310
485  c2r = dabs(ckr)
486  c2i = dabs(cki)
487  c2m = dmax1(c2r,c2i)
488  IF (c2m.LE.ascle) GO TO 310
489  iflag = iflag + 1
490  ascle = bry(iflag)
491  s1r = s1r*csr
492  s1i = s1i*csr
493  s2r = ckr
494  s2i = cki
495  s1r = s1r*cssr(iflag)
496  s1i = s1i*cssr(iflag)
497  s2r = s2r*cssr(iflag)
498  s2i = s2i*cssr(iflag)
499  csr = csrr(iflag)
500  310 CONTINUE
501  RETURN
502  320 CONTINUE
503  nz = -1
504  RETURN
505  END
d1mach
double precision function d1mach(i)
Definition: d1mach.f:23
mod
T mod(T x, T y)
Definition: lo-mappers.h:294
xzabs
double precision function xzabs(ZR, ZI)
Definition: xzabs.f:2
zairy
subroutine zairy(ZR, ZI, ID, KODE, AIR, AII, NZ, IERR)
Definition: zairy.f:2
zs1s2
subroutine zs1s2(ZRR, ZRI, S1R, S1I, S2R, S2I, NZ, ASCLE, ALIM, IUF)
Definition: zs1s2.f:3
zuchk
subroutine zuchk(YR, YI, NZ, ASCLE, TOL)
Definition: zuchk.f:2
zunhj
subroutine zunhj(ZR, ZI, FNU, IPMTR, TOL, PHIR, PHII, ARGR, ARGI, ZETA1R, ZETA1I, ZETA2R, ZETA2I, ASUMR, ASUMI, BSUMR, BSUMI)
Definition: zunhj.f:3
zunk2
subroutine zunk2(ZR, ZI, FNU, KODE, MR, N, YR, YI, NZ, TOL, ELIM, ALIM)
Definition: zunk2.f:3