1 SUBROUTINE zbiry(ZR, ZI, ID, KODE, BIR, BII, IERR)
124 DOUBLE PRECISION AA, AD, AK, ALIM, ATRM, AZ, AZ3, BB, BII, BIR,
125 * bk, cc, ck, coef, conei, coner, csqi, csqr, cyi, cyr, c1, c2,
126 * dig, dk, d1, d2, eaa, elim, fid, fmr, fnu, fnul, pi, rl, r1m5,
127 * sfac, sti,
str, s1i, s1r, s2i, s2r, tol, trm1i, trm1r, trm2i,
128 * trm2r, tth, zi, zr, ztai, ztar, z3i, z3r,
d1mach, xzabs
129 INTEGER ID, IERR, K, KODE, K1, K2, NZ, I1MACH
131 DATA tth, c1, c2, coef, pi /6.66666666666666667
d-01,
132 * 6.14926627446000736
d-01,4.48288357353826359
d-01,
133 * 5.77350269189625765
d-01,3.14159265358979324
d+00/
134 DATA coner, conei /1.0d0,0.0d0/
138 IF (id.LT.0 .OR. id.GT.1) ierr=1
139 IF (kode.LT.1 .OR. kode.GT.2) ierr=1
140 IF (ierr.NE.0)
RETURN
143 fid = dble(float(id))
144 IF (az.GT.1.0e0) go to 70
152 IF (az.LT.tol) go to 130
154 IF (aa.LT.tol/az) go to 40
162 z3r =
str*zr - sti*zi
163 z3i =
str*zi + sti*zr
166 bk = 3.0d0 - fid - fid
168 dk = 3.0d0 + fid + fid
172 ak = 24.0d0 + 9.0d0*fid
173 bk = 30.0d0 - 9.0d0*fid
175 str = (trm1r*z3r-trm1i*z3i)/d1
176 trm1i = (trm1r*z3i+trm1i*z3r)/d1
180 str = (trm2r*z3r-trm2i*z3i)/d2
181 trm2i = (trm2r*z3i+trm2i*z3r)/d2
189 IF (atrm.LT.tol*ad) go to 40
194 IF (id.EQ.1) go to 50
195 bir = c1*s1r + c2*(zr*s2r-zi*s2i)
196 bii = c1*s1i + c2*(zr*s2i+zi*s2r)
197 IF (kode.EQ.1)
RETURN
199 ztar = tth*(zr*
str-zi*sti)
200 ztai = tth*(zr*sti+zi*
str)
210 IF (az.LE.tol) go to 60
212 str = s1r*zr - s1i*zi
213 sti = s1r*zi + s1i*zr
214 bir = bir + cc*(
str*zr-sti*zi)
215 bii = bii + cc*(
str*zi+sti*zr)
217 IF (kode.EQ.1)
RETURN
219 ztar = tth*(zr*
str-zi*sti)
220 ztai = tth*(zr*sti+zi*
str)
231 fnu = (1.0d0+fid)/3.0d0
246 k = min0(iabs(k1),iabs(k2))
247 elim = 2.303d0*(dble(float(k))*r1m5-3.0d0)
249 aa = r1m5*dble(float(k1))
250 dig = dmin1(aa,18.0d0)
252 alim = elim + dmax1(-aa,-41.45d0)
253 rl = 1.2d0*dig + 3.0d0
254 fnul = 10.0d0 + 6.0d0*(dig-3.0d0)
259 bb=dble(float(i1mach(9)))*0.5d0
262 IF (az.GT.aa) go to 260
265 CALL
xzsqrt(zr, zi, csqr, csqi)
266 ztar = tth*(zr*csqr-zi*csqi)
267 ztai = tth*(zr*csqi+zi*csqr)
273 IF (zr.GE.0.0d0) go to 80
279 IF (zi.NE.0.0d0 .OR. zr.GT.0.0d0) go to 90
284 IF (kode.EQ.2) go to 100
289 IF (bb.LT.alim) go to 100
290 bb = bb + 0.25d0*dlog(az)
292 IF (bb.GT.elim) go to 190
295 IF (aa.GE.0.0d0 .AND. zr.GT.0.0d0) go to 110
297 IF (zi.LT.0.0d0) fmr = -pi
305 CALL
zbinu(ztar, ztai, fnu, kode, 1, cyr, cyi, nz, rl, fnul, tol,
307 IF (nz.LT.0) go to 200
312 s1r = (
str*cyr(1)-sti*cyi(1))*z3r
313 s1i = (
str*cyi(1)+sti*cyr(1))*z3r
314 fnu = (2.0d0-fid)/3.0d0
315 CALL
zbinu(ztar, ztai, fnu, kode, 2, cyr, cyi, nz, rl, fnul, tol,
324 CALL
zdiv(cyr(1), cyi(1), ztar, ztai,
str, sti)
325 s2r = (fnu+fnu)*
str + cyr(2)
326 s2i = (fnu+fnu)*sti + cyi(2)
330 s1r = coef*(s1r+s2r*
str-s2i*sti)
331 s1i = coef*(s1i+s2r*sti+s2i*
str)
332 IF (id.EQ.1) go to 120
333 str = csqr*s1r - csqi*s1i
334 s1i = csqr*s1i + csqi*s1r
340 str = zr*s1r - zi*s1i
341 s1i = zr*s1i + zi*s1r
347 aa = c1*(1.0d0-fid) + fid*c2
356 IF(nz.EQ.(-1)) go to 190
subroutine zbinu(ZR, ZI, FNU, KODE, N, CYR, CYI, NZ, RL, FNUL, TOL, ELIM, ALIM)
std::string str(char sep= 'x') const
std::string dimension(void) const
double precision function d1mach(i)
F77_RET_T const double const double double * d
subroutine zdiv(AR, AI, BR, BI, CR, CI)
subroutine zbiry(ZR, ZI, ID, KODE, BIR, BII, IERR)
subroutine xzsqrt(AR, AI, BR, BI)