1 SUBROUTINE zunk1(ZR, ZI, FNU, KODE, MR, N, YR, YI, NZ, TOL, ELIM,
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 /
50 IF (zr.GE.0.0d0) go to 10
60 fn = fnu + dble(float(i-1))
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
68 rast = fn/xzabs(
str,sti)
75 s1r = zeta1r(j) - zeta2r(j)
76 s1i = zeta1i(j) - zeta2i(j)
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
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
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)
104 str = s2r*s1r - s2i*s1i
105 s2i = s1r*s2i + s2r*s1i
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
113 yr(i) = s2r*csrr(kflag)
114 yi(i) = s2i*csrr(kflag)
115 IF (kdflg.EQ.2) go to 75
119 IF (rs1.GT.0.0d0) go to 300
123 IF (zr.LT.0.0d0) go to 300
129 IF ((yr(i-1).EQ.zeror).AND.(yi(i-1).EQ.zeroi)) go to 70
136 razr = 1.0d0/xzabs(zrr,zri)
144 IF (n.LT.ib) go to 160
149 fn = fnu + dble(float(n-1))
151 IF (mr.NE.0) ipard = 0
153 CALL
zunik(zrr, zri, fn, 2, ipard, tol, initd, phidr, phidi,
154 * zet1dr, zet1di, zet2dr, zet2di, sumdr, sumdi, cwrkr(1,3),
156 IF (kode.EQ.1) go to 80
159 rast = fn/xzabs(
str,sti)
166 s1r = zet1dr - zet2dr
167 s1i = zet1di - zet2di
170 IF (dabs(rs1).GT.elim) go to 95
171 IF (dabs(rs1).LT.alim) go to 100
175 aphi = xzabs(phidr,phidi)
177 IF (dabs(rs1).LT.elim) go to 100
179 IF (dabs(rs1).GT.0.0d0) go to 300
183 IF (zr.LT.0.0d0) go to 300
203 s2r = ckr*c2r - cki*c2i + s1r
204 s2i = ckr*c2i + cki*c2r + s1i
213 IF (kflag.GE.3) go to 120
217 IF (c2m.LE.ascle) go to 120
224 s1r = s1r*cssr(kflag)
225 s1i = s1i*cssr(kflag)
226 s2r = s2r*cssr(kflag)
227 s2i = s2i*cssr(kflag)
236 fmr = dble(float(mr))
243 fnf = fnu - dble(float(inu))
248 IF (
mod(ifn,2).EQ.0) go to 170
259 fn = fnu + dble(float(kk-1))
265 IF (n.GT.2) go to 175
280 IF ((kk.EQ.n).AND.(ib.LT.n)) go to 180
281 IF ((kk.EQ.ib).OR.(kk.EQ.ic)) go to 172
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
290 rast = fn/xzabs(
str,sti)
297 s1r = -zet1dr + zet2dr
298 s1i = -zet1di + zet2di
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
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
317 str = phidr*sumdr - phidi*sumdi
318 sti = phidr*sumdi + phidi*sumdr
321 str = dexp(s1r)*cssr(iflag)
324 str = s2r*s1r - s2i*s1i
325 s2i = s2r*s1i + s2i*s1r
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
337 s2r = s2r*csrr(iflag)
338 s2i = s2i*csrr(iflag)
344 IF (kode.EQ.1) go to 250
345 CALL
zs1s2(zrr, zri, s1r, s1i, s2r, s2i, nw, asc, alim, iuf)
348 yr(kk) = s1r*cspnr - s1i*cspni + s2r
349 yi(kk) = cspnr*s1i + cspni*s1r + s2i
353 IF (c2r.NE.0.0d0 .OR. c2i.NE.0.0d0) go to 255
357 IF (kdflg.EQ.2) go to 275
361 IF (rs1.GT.0.0d0) go to 300
381 fn = dble(float(inu+il))
385 s2r = s1r + (fn+fnf)*(rzr*c2r-rzi*c2i)
386 s2i = s1i + (fn+fnf)*(rzr*c2i+rzi*c2r)
396 IF (kode.EQ.1) go to 280
397 CALL
zs1s2(zrr, zri, c1r, c1i, c2r, c2i, nw, asc, alim, iuf)
400 yr(kk) = c1r*cspnr - c1i*cspni + c2r
401 yi(kk) = c1r*cspni + c1i*cspnr + c2i
405 IF (iflag.GE.3) go to 290
409 IF (c2m.LE.ascle) go to 290
416 s1r = s1r*cssr(iflag)
417 s1i = s1i*cssr(iflag)
418 s2r = s2r*cssr(iflag)
419 s2i = s2i*cssr(iflag)
std::string str(char sep= 'x') const
std::string dimension(void) const
octave_int< T > mod(const octave_int< T > &x, const octave_int< T > &y)
double precision function d1mach(i)
F77_RET_T const double const double double * d
subroutine zs1s2(ZRR, ZRI, S1R, S1I, S2R, S2I, NZ, ASCLE, ALIM, IUF)
subroutine zunik(ZRR, ZRI, FNU, IKFLG, IPMTR, TOL, INIT, PHIR, PHII, ZETA1R, ZETA1I, ZETA2R, ZETA2I, SUMR, SUMI, CWRKR, CWRKI)
subroutine zuchk(YR, YI, NZ, ASCLE, TOL)
subroutine zunk1(ZR, ZI, FNU, KODE, MR, N, YR, YI, NZ, TOL, ELIM, ALIM)