1*DECK ZUNI1 2 SUBROUTINE ZUNI1 (ZR, ZI, FNU, KODE, N, YR, YI, NZ, NLAST, FNUL, 3 + TOL, ELIM, ALIM) 4C***BEGIN PROLOGUE ZUNI1 5C***SUBSIDIARY 6C***PURPOSE Subsidiary to ZBESI and ZBESK 7C***LIBRARY SLATEC 8C***TYPE ALL (CUNI1-A, ZUNI1-A) 9C***AUTHOR Amos, D. E., (SNL) 10C***DESCRIPTION 11C 12C ZUNI1 COMPUTES I(FNU,Z) BY MEANS OF THE UNIFORM ASYMPTOTIC 13C EXPANSION FOR I(FNU,Z) IN -PI/3.LE.ARG Z.LE.PI/3. 14C 15C FNUL IS THE SMALLEST ORDER PERMITTED FOR THE ASYMPTOTIC 16C EXPANSION. NLAST=0 MEANS ALL OF THE Y VALUES WERE SET. 17C NLAST.NE.0 IS THE NUMBER LEFT TO BE COMPUTED BY ANOTHER 18C FORMULA FOR ORDERS FNU TO FNU+NLAST-1 BECAUSE FNU+NLAST-1.LT.FNUL. 19C Y(I)=CZERO FOR I=NLAST+1,N 20C 21C***SEE ALSO ZBESI, ZBESK 22C***ROUTINES CALLED D1MACH, ZABS, ZUCHK, ZUNIK, ZUOIK 23C***REVISION HISTORY (YYMMDD) 24C 830501 DATE WRITTEN 25C 910415 Prologue converted to Version 4.0 format. (BAB) 26C***END PROLOGUE ZUNI1 27C COMPLEX CFN,CONE,CRSC,CSCL,CSR,CSS,CWRK,CZERO,C1,C2,PHI,RZ,SUM,S1, 28C *S2,Y,Z,ZETA1,ZETA2 29 DOUBLE PRECISION ALIM, APHI, ASCLE, BRY, CONER, CRSC, 30 * CSCL, CSRR, CSSR, CWRKI, CWRKR, C1R, C2I, C2M, C2R, ELIM, FN, 31 * FNU, FNUL, PHII, PHIR, RAST, RS1, RZI, RZR, STI, STR, SUMI, 32 * SUMR, S1I, S1R, S2I, S2R, TOL, YI, YR, ZEROI, ZEROR, ZETA1I, 33 * ZETA1R, ZETA2I, ZETA2R, ZI, ZR, CYR, CYI, D1MACH, ZABS 34 INTEGER I, IFLAG, INIT, K, KODE, M, N, ND, NLAST, NN, NUF, NW, NZ 35 DIMENSION BRY(3), YR(N), YI(N), CWRKR(16), CWRKI(16), CSSR(3), 36 * CSRR(3), CYR(2), CYI(2) 37 EXTERNAL ZABS 38 DATA ZEROR,ZEROI,CONER / 0.0D0, 0.0D0, 1.0D0 / 39C***FIRST EXECUTABLE STATEMENT ZUNI1 40 NZ = 0 41 ND = N 42 NLAST = 0 43C----------------------------------------------------------------------- 44C COMPUTED VALUES WITH EXPONENTS BETWEEN ALIM AND ELIM IN MAG- 45C NITUDE ARE SCALED TO KEEP INTERMEDIATE ARITHMETIC ON SCALE, 46C EXP(ALIM)=EXP(ELIM)*TOL 47C----------------------------------------------------------------------- 48 CSCL = 1.0D0/TOL 49 CRSC = TOL 50 CSSR(1) = CSCL 51 CSSR(2) = CONER 52 CSSR(3) = CRSC 53 CSRR(1) = CRSC 54 CSRR(2) = CONER 55 CSRR(3) = CSCL 56 BRY(1) = 1.0D+3*D1MACH(1)/TOL 57C----------------------------------------------------------------------- 58C CHECK FOR UNDERFLOW AND OVERFLOW ON FIRST MEMBER 59C----------------------------------------------------------------------- 60 FN = MAX(FNU,1.0D0) 61 INIT = 0 62 CALL ZUNIK(ZR, ZI, FN, 1, 1, TOL, INIT, PHIR, PHII, ZETA1R, 63 * ZETA1I, ZETA2R, ZETA2I, SUMR, SUMI, CWRKR, CWRKI) 64 IF (KODE.EQ.1) GO TO 10 65 STR = ZR + ZETA2R 66 STI = ZI + ZETA2I 67 RAST = FN/ZABS(STR,STI) 68 STR = STR*RAST*RAST 69 STI = -STI*RAST*RAST 70 S1R = -ZETA1R + STR 71 S1I = -ZETA1I + STI 72 GO TO 20 73 10 CONTINUE 74 S1R = -ZETA1R + ZETA2R 75 S1I = -ZETA1I + ZETA2I 76 20 CONTINUE 77 RS1 = S1R 78 IF (ABS(RS1).GT.ELIM) GO TO 130 79 30 CONTINUE 80 NN = MIN(2,ND) 81 DO 80 I=1,NN 82 FN = FNU + (ND-I) 83 INIT = 0 84 CALL ZUNIK(ZR, ZI, FN, 1, 0, TOL, INIT, PHIR, PHII, ZETA1R, 85 * ZETA1I, ZETA2R, ZETA2I, SUMR, SUMI, CWRKR, CWRKI) 86 IF (KODE.EQ.1) GO TO 40 87 STR = ZR + ZETA2R 88 STI = ZI + ZETA2I 89 RAST = FN/ZABS(STR,STI) 90 STR = STR*RAST*RAST 91 STI = -STI*RAST*RAST 92 S1R = -ZETA1R + STR 93 S1I = -ZETA1I + STI + ZI 94 GO TO 50 95 40 CONTINUE 96 S1R = -ZETA1R + ZETA2R 97 S1I = -ZETA1I + ZETA2I 98 50 CONTINUE 99C----------------------------------------------------------------------- 100C TEST FOR UNDERFLOW AND OVERFLOW 101C----------------------------------------------------------------------- 102 RS1 = S1R 103 IF (ABS(RS1).GT.ELIM) GO TO 110 104 IF (I.EQ.1) IFLAG = 2 105 IF (ABS(RS1).LT.ALIM) GO TO 60 106C----------------------------------------------------------------------- 107C REFINE TEST AND SCALE 108C----------------------------------------------------------------------- 109 APHI = ZABS(PHIR,PHII) 110 RS1 = RS1 + LOG(APHI) 111 IF (ABS(RS1).GT.ELIM) GO TO 110 112 IF (I.EQ.1) IFLAG = 1 113 IF (RS1.LT.0.0D0) GO TO 60 114 IF (I.EQ.1) IFLAG = 3 115 60 CONTINUE 116C----------------------------------------------------------------------- 117C SCALE S1 IF ABS(S1).LT.ASCLE 118C----------------------------------------------------------------------- 119 S2R = PHIR*SUMR - PHII*SUMI 120 S2I = PHIR*SUMI + PHII*SUMR 121 STR = EXP(S1R)*CSSR(IFLAG) 122 S1R = STR*COS(S1I) 123 S1I = STR*SIN(S1I) 124 STR = S2R*S1R - S2I*S1I 125 S2I = S2R*S1I + S2I*S1R 126 S2R = STR 127 IF (IFLAG.NE.1) GO TO 70 128 CALL ZUCHK(S2R, S2I, NW, BRY(1), TOL) 129 IF (NW.NE.0) GO TO 110 130 70 CONTINUE 131 CYR(I) = S2R 132 CYI(I) = S2I 133 M = ND - I + 1 134 YR(M) = S2R*CSRR(IFLAG) 135 YI(M) = S2I*CSRR(IFLAG) 136 80 CONTINUE 137 IF (ND.LE.2) GO TO 100 138 RAST = 1.0D0/ZABS(ZR,ZI) 139 STR = ZR*RAST 140 STI = -ZI*RAST 141 RZR = (STR+STR)*RAST 142 RZI = (STI+STI)*RAST 143 BRY(2) = 1.0D0/BRY(1) 144 BRY(3) = D1MACH(2) 145 S1R = CYR(1) 146 S1I = CYI(1) 147 S2R = CYR(2) 148 S2I = CYI(2) 149 C1R = CSRR(IFLAG) 150 ASCLE = BRY(IFLAG) 151 K = ND - 2 152 FN = K 153 DO 90 I=3,ND 154 C2R = S2R 155 C2I = S2I 156 S2R = S1R + (FNU+FN)*(RZR*C2R-RZI*C2I) 157 S2I = S1I + (FNU+FN)*(RZR*C2I+RZI*C2R) 158 S1R = C2R 159 S1I = C2I 160 C2R = S2R*C1R 161 C2I = S2I*C1R 162 YR(K) = C2R 163 YI(K) = C2I 164 K = K - 1 165 FN = FN - 1.0D0 166 IF (IFLAG.GE.3) GO TO 90 167 STR = ABS(C2R) 168 STI = ABS(C2I) 169 C2M = MAX(STR,STI) 170 IF (C2M.LE.ASCLE) GO TO 90 171 IFLAG = IFLAG + 1 172 ASCLE = BRY(IFLAG) 173 S1R = S1R*C1R 174 S1I = S1I*C1R 175 S2R = C2R 176 S2I = C2I 177 S1R = S1R*CSSR(IFLAG) 178 S1I = S1I*CSSR(IFLAG) 179 S2R = S2R*CSSR(IFLAG) 180 S2I = S2I*CSSR(IFLAG) 181 C1R = CSRR(IFLAG) 182 90 CONTINUE 183 100 CONTINUE 184 RETURN 185C----------------------------------------------------------------------- 186C SET UNDERFLOW AND UPDATE PARAMETERS 187C----------------------------------------------------------------------- 188 110 CONTINUE 189 IF (RS1.GT.0.0D0) GO TO 120 190 YR(ND) = ZEROR 191 YI(ND) = ZEROI 192 NZ = NZ + 1 193 ND = ND - 1 194 IF (ND.EQ.0) GO TO 100 195 CALL ZUOIK(ZR, ZI, FNU, KODE, 1, ND, YR, YI, NUF, TOL, ELIM, ALIM) 196 IF (NUF.LT.0) GO TO 120 197 ND = ND - NUF 198 NZ = NZ + NUF 199 IF (ND.EQ.0) GO TO 100 200 FN = FNU + (ND-1) 201 IF (FN.GE.FNUL) GO TO 30 202 NLAST = ND 203 RETURN 204 120 CONTINUE 205 NZ = -1 206 RETURN 207 130 CONTINUE 208 IF (RS1.GT.0.0D0) GO TO 120 209 NZ = N 210 DO 140 I=1,N 211 YR(I) = ZEROR 212 YI(I) = ZEROI 213 140 CONTINUE 214 RETURN 215 END 216