SUBROUTINE ZUOIK(ZR, ZI, FNU, KODE, IKFLG, N, YR, YI, NUF, TOL, * ELIM, ALIM) C***BEGIN PROLOGUE ZUOIK C***REFER TO ZBESI,ZBESK,ZBESH C C ZUOIK COMPUTES THE LEADING TERMS OF THE UNIFORM ASYMPTOTIC C EXPANSIONS FOR THE I AND K FUNCTIONS AND COMPARES THEM C (IN LOGARITHMIC FORM) TO ALIM AND ELIM FOR OVER AND UNDERFLOW C WHERE ALIM.LT.ELIM. IF THE MAGNITUDE, BASED ON THE LEADING C EXPONENTIAL, IS LESS THAN ALIM OR GREATER THAN -ALIM, THEN C THE RESULT IS ON SCALE. IF NOT, THEN A REFINED TEST USING OTHER C MULTIPLIERS (IN LOGARITHMIC FORM) IS MADE BASED ON ELIM. HERE C EXP(-ELIM)=SMALLEST MACHINE NUMBER*1.0E+3 AND EXP(-ALIM)= C EXP(-ELIM)/TOL C C IKFLG=1 MEANS THE I SEQUENCE IS TESTED C =2 MEANS THE K SEQUENCE IS TESTED C NUF = 0 MEANS THE LAST MEMBER OF THE SEQUENCE IS ON SCALE C =-1 MEANS AN OVERFLOW WOULD OCCUR C IKFLG=1 AND NUF.GT.0 MEANS THE LAST NUF Y VALUES WERE SET TO ZERO C THE FIRST N-NUF VALUES MUST BE SET BY ANOTHER ROUTINE C IKFLG=2 AND NUF.EQ.N MEANS ALL Y VALUES WERE SET TO ZERO C IKFLG=2 AND 0.LT.NUF.LT.N NOT CONSIDERED. Y MUST BE SET BY C ANOTHER ROUTINE C C***ROUTINES CALLED ZUCHK,ZUNHJ,ZUNIK,D1MACH,AZABS,AZLOG C***END PROLOGUE ZUOIK C COMPLEX ARG,ASUM,BSUM,CWRK,CZ,CZERO,PHI,SUM,Y,Z,ZB,ZETA1,ZETA2,ZN, C *ZR DOUBLE PRECISION AARG, AIC, ALIM, APHI, ARGI, ARGR, ASUMI, ASUMR, * ASCLE, AX, AY, BSUMI, BSUMR, CWRKI, CWRKR, CZI, CZR, ELIM, FNN, * FNU, GNN, GNU, PHII, PHIR, RCZ, STR, STI, SUMI, SUMR, TOL, YI, * YR, ZBI, ZBR, ZEROI, ZEROR, ZETA1I, ZETA1R, ZETA2I, ZETA2R, ZI, * ZNI, ZNR, ZR, ZRI, ZRR, D1MACH, AZABS INTEGER I, IDUM, IFORM, IKFLG, INIT, KODE, N, NN, NUF, NW DIMENSION YR(N), YI(N), CWRKR(16), CWRKI(16) DATA ZEROR,ZEROI / 0.0D0, 0.0D0 / DATA AIC / 1.265512123484645396D+00 / NUF = 0 NN = N ZRR = ZR ZRI = ZI IF (ZR.GE.0.0D0) GO TO 10 ZRR = -ZR ZRI = -ZI 10 CONTINUE ZBR = ZRR ZBI = ZRI AX = DABS(ZR)*1.7321D0 AY = DABS(ZI) IFORM = 1 IF (AY.GT.AX) IFORM = 2 GNU = DMAX1(FNU,1.0D0) IF (IKFLG.EQ.1) GO TO 20 FNN = DBLE(FLOAT(NN)) GNN = FNU + FNN - 1.0D0 GNU = DMAX1(GNN,FNN) 20 CONTINUE C----------------------------------------------------------------------- C ONLY THE MAGNITUDE OF ARG AND PHI ARE NEEDED ALONG WITH THE C REAL PARTS OF ZETA1, ZETA2 AND ZB. NO ATTEMPT IS MADE TO GET C THE SIGN OF THE IMAGINARY PART CORRECT. C----------------------------------------------------------------------- IF (IFORM.EQ.2) GO TO 30 INIT = 0 CALL ZUNIK(ZRR, ZRI, GNU, IKFLG, 1, TOL, INIT, PHIR, PHII, * ZETA1R, ZETA1I, ZETA2R, ZETA2I, SUMR, SUMI, CWRKR, CWRKI) CZR = -ZETA1R + ZETA2R CZI = -ZETA1I + ZETA2I GO TO 50 30 CONTINUE ZNR = ZRI ZNI = -ZRR IF (ZI.GT.0.0D0) GO TO 40 ZNR = -ZNR 40 CONTINUE CALL ZUNHJ(ZNR, ZNI, GNU, 1, TOL, PHIR, PHII, ARGR, ARGI, ZETA1R, * ZETA1I, ZETA2R, ZETA2I, ASUMR, ASUMI, BSUMR, BSUMI) CZR = -ZETA1R + ZETA2R CZI = -ZETA1I + ZETA2I AARG = AZABS(ARGR,ARGI) 50 CONTINUE IF (KODE.EQ.1) GO TO 60 CZR = CZR - ZBR CZI = CZI - ZBI 60 CONTINUE IF (IKFLG.EQ.1) GO TO 70 CZR = -CZR CZI = -CZI 70 CONTINUE APHI = AZABS(PHIR,PHII) RCZ = CZR C----------------------------------------------------------------------- C OVERFLOW TEST C----------------------------------------------------------------------- IF (RCZ.GT.ELIM) GO TO 210 IF (RCZ.LT.ALIM) GO TO 80 RCZ = RCZ + DLOG(APHI) IF (IFORM.EQ.2) RCZ = RCZ - 0.25D0*DLOG(AARG) - AIC IF (RCZ.GT.ELIM) GO TO 210 GO TO 130 80 CONTINUE C----------------------------------------------------------------------- C UNDERFLOW TEST C----------------------------------------------------------------------- IF (RCZ.LT.(-ELIM)) GO TO 90 IF (RCZ.GT.(-ALIM)) GO TO 130 RCZ = RCZ + DLOG(APHI) IF (IFORM.EQ.2) RCZ = RCZ - 0.25D0*DLOG(AARG) - AIC IF (RCZ.GT.(-ELIM)) GO TO 110 90 CONTINUE DO 100 I=1,NN YR(I) = ZEROR YI(I) = ZEROI 100 CONTINUE NUF = NN RETURN 110 CONTINUE ASCLE = 1.0D+3*D1MACH(1)/TOL CALL AZLOG(PHIR, PHII, STR, STI, IDUM) CZR = CZR + STR CZI = CZI + STI IF (IFORM.EQ.1) GO TO 120 CALL AZLOG(ARGR, ARGI, STR, STI, IDUM) CZR = CZR - 0.25D0*STR - AIC CZI = CZI - 0.25D0*STI 120 CONTINUE AX = DEXP(RCZ)/TOL AY = CZI CZR = AX*DCOS(AY) CZI = AX*DSIN(AY) CALL ZUCHK(CZR, CZI, NW, ASCLE, TOL) IF (NW.NE.0) GO TO 90 130 CONTINUE IF (IKFLG.EQ.2) RETURN IF (N.EQ.1) RETURN C----------------------------------------------------------------------- C SET UNDERFLOWS ON I SEQUENCE C----------------------------------------------------------------------- 140 CONTINUE GNU = FNU + DBLE(FLOAT(NN-1)) IF (IFORM.EQ.2) GO TO 150 INIT = 0 CALL ZUNIK(ZRR, ZRI, GNU, IKFLG, 1, TOL, INIT, PHIR, PHII, * ZETA1R, ZETA1I, ZETA2R, ZETA2I, SUMR, SUMI, CWRKR, CWRKI) CZR = -ZETA1R + ZETA2R CZI = -ZETA1I + ZETA2I GO TO 160 150 CONTINUE CALL ZUNHJ(ZNR, ZNI, GNU, 1, TOL, PHIR, PHII, ARGR, ARGI, ZETA1R, * ZETA1I, ZETA2R, ZETA2I, ASUMR, ASUMI, BSUMR, BSUMI) CZR = -ZETA1R + ZETA2R CZI = -ZETA1I + ZETA2I AARG = AZABS(ARGR,ARGI) 160 CONTINUE IF (KODE.EQ.1) GO TO 170 CZR = CZR - ZBR CZI = CZI - ZBI 170 CONTINUE APHI = AZABS(PHIR,PHII) RCZ = CZR IF (RCZ.LT.(-ELIM)) GO TO 180 IF (RCZ.GT.(-ALIM)) RETURN RCZ = RCZ + DLOG(APHI) IF (IFORM.EQ.2) RCZ = RCZ - 0.25D0*DLOG(AARG) - AIC IF (RCZ.GT.(-ELIM)) GO TO 190 180 CONTINUE YR(NN) = ZEROR YI(NN) = ZEROI NN = NN - 1 NUF = NUF + 1 IF (NN.EQ.0) RETURN GO TO 140 190 CONTINUE ASCLE = 1.0D+3*D1MACH(1)/TOL CALL AZLOG(PHIR, PHII, STR, STI, IDUM) CZR = CZR + STR CZI = CZI + STI IF (IFORM.EQ.1) GO TO 200 CALL AZLOG(ARGR, ARGI, STR, STI, IDUM) CZR = CZR - 0.25D0*STR - AIC CZI = CZI - 0.25D0*STI 200 CONTINUE AX = DEXP(RCZ)/TOL AY = CZI CZR = AX*DCOS(AY) CZI = AX*DSIN(AY) CALL ZUCHK(CZR, CZI, NW, ASCLE, TOL) IF (NW.NE.0) GO TO 180 RETURN 210 CONTINUE NUF = -1 RETURN END