1      SUBROUTINE CBUNI(Z, FNU, KODE, N, Y, NZ, NUI, NLAST, FNUL, TOL,
2     * ELIM, ALIM)
3C***BEGIN PROLOGUE  CBUNI
4C***REFER TO  CBESI,CBESK
5C
6C     CBUNI COMPUTES THE I BESSEL FUNCTION FOR LARGE CABS(Z).GT.
7C     FNUL AND FNU+N-1.LT.FNUL. THE ORDER IS INCREASED FROM
8C     FNU+N-1 GREATER THAN FNUL BY ADDING NUI AND COMPUTING
9C     ACCORDING TO THE UNIFORM ASYMPTOTIC EXPANSION FOR I(FNU,Z)
10C     ON IFORM=1 AND THE EXPANSION FOR J(FNU,Z) ON IFORM=2
11C
12C***ROUTINES CALLED  CUNI1,CUNI2,R1MACH
13C***END PROLOGUE  CBUNI
14      COMPLEX CSCL, CSCR, CY, RZ, ST, S1, S2, Y, Z
15      REAL ALIM, AX, AY, DFNU, ELIM, FNU, FNUI, FNUL, GNU, TOL, XX, YY,
16     * ASCLE, BRY, STR, STI, STM, R1MACH
17      INTEGER I, IFLAG, IFORM, K, KODE, N, NL, NLAST, NUI, NW, NZ
18      DIMENSION Y(N), CY(2), BRY(3)
19      NZ = 0
20      XX = REAL(Z)
21      YY = AIMAG(Z)
22      AX = ABS(XX)*1.7321E0
23      AY = ABS(YY)
24      IFORM = 1
25      IF (AY.GT.AX) IFORM = 2
26      IF (NUI.EQ.0) GO TO 60
27      FNUI = FLOAT(NUI)
28      DFNU = FNU + FLOAT(N-1)
29      GNU = DFNU + FNUI
30      IF (IFORM.EQ.2) GO TO 10
31C-----------------------------------------------------------------------
32C     ASYMPTOTIC EXPANSION FOR I(FNU,Z) FOR LARGE FNU APPLIED IN
33C     -PI/3.LE.ARG(Z).LE.PI/3
34C-----------------------------------------------------------------------
35      CALL CUNI1(Z, GNU, KODE, 2, CY, NW, NLAST, FNUL, TOL, ELIM, ALIM)
36      GO TO 20
37   10 CONTINUE
38C-----------------------------------------------------------------------
39C     ASYMPTOTIC EXPANSION FOR J(FNU,Z*EXP(M*HPI)) FOR LARGE FNU
40C     APPLIED IN PI/3.LT.ABS(ARG(Z)).LE.PI/2 WHERE M=+I OR -I
41C     AND HPI=PI/2
42C-----------------------------------------------------------------------
43      CALL CUNI2(Z, GNU, KODE, 2, CY, NW, NLAST, FNUL, TOL, ELIM, ALIM)
44   20 CONTINUE
45      IF (NW.LT.0) GO TO 50
46      IF (NW.NE.0) GO TO 90
47      AY = CABS(CY(1))
48C----------------------------------------------------------------------
49C     SCALE BACKWARD RECURRENCE, BRY(3) IS DEFINED BUT NEVER USED
50C----------------------------------------------------------------------
51      BRY(1) = 1.0E+3*R1MACH(1)/TOL
52      BRY(2) = 1.0E0/BRY(1)
53      BRY(3) = BRY(2)
54      IFLAG = 2
55      ASCLE = BRY(2)
56      AX = 1.0E0
57      CSCL = CMPLX(AX,0.0E0)
58      IF (AY.GT.BRY(1)) GO TO 21
59      IFLAG = 1
60      ASCLE = BRY(1)
61      AX = 1.0E0/TOL
62      CSCL = CMPLX(AX,0.0E0)
63      GO TO 25
64   21 CONTINUE
65      IF (AY.LT.BRY(2)) GO TO 25
66      IFLAG = 3
67      ASCLE = BRY(3)
68      AX = TOL
69      CSCL = CMPLX(AX,0.0E0)
70   25 CONTINUE
71      AY = 1.0E0/AX
72      CSCR = CMPLX(AY,0.0E0)
73      S1 = CY(2)*CSCL
74      S2 = CY(1)*CSCL
75      RZ = CMPLX(2.0E0,0.0E0)/Z
76      DO 30 I=1,NUI
77        ST = S2
78        S2 = CMPLX(DFNU+FNUI,0.0E0)*RZ*S2 + S1
79        S1 = ST
80        FNUI = FNUI - 1.0E0
81        IF (IFLAG.GE.3) GO TO 30
82        ST = S2*CSCR
83        STR = REAL(ST)
84        STI = AIMAG(ST)
85        STR = ABS(STR)
86        STI = ABS(STI)
87        STM = AMAX1(STR,STI)
88        IF (STM.LE.ASCLE) GO TO 30
89        IFLAG = IFLAG+1
90        ASCLE = BRY(IFLAG)
91        S1 = S1*CSCR
92        S2 = ST
93        AX = AX*TOL
94        AY = 1.0E0/AX
95        CSCL = CMPLX(AX,0.0E0)
96        CSCR = CMPLX(AY,0.0E0)
97        S1 = S1*CSCL
98        S2 = S2*CSCL
99   30 CONTINUE
100      Y(N) = S2*CSCR
101      IF (N.EQ.1) RETURN
102      NL = N - 1
103      FNUI = FLOAT(NL)
104      K = NL
105      DO 40 I=1,NL
106        ST = S2
107        S2 = CMPLX(FNU+FNUI,0.0E0)*RZ*S2 + S1
108        S1 = ST
109        ST = S2*CSCR
110        Y(K) = ST
111        FNUI = FNUI - 1.0E0
112        K = K - 1
113        IF (IFLAG.GE.3) GO TO 40
114        STR = REAL(ST)
115        STI = AIMAG(ST)
116        STR = ABS(STR)
117        STI = ABS(STI)
118        STM = AMAX1(STR,STI)
119        IF (STM.LE.ASCLE) GO TO 40
120        IFLAG = IFLAG+1
121        ASCLE = BRY(IFLAG)
122        S1 = S1*CSCR
123        S2 = ST
124        AX = AX*TOL
125        AY = 1.0E0/AX
126        CSCL = CMPLX(AX,0.0E0)
127        CSCR = CMPLX(AY,0.0E0)
128        S1 = S1*CSCL
129        S2 = S2*CSCL
130   40 CONTINUE
131      RETURN
132   50 CONTINUE
133      NZ = -1
134      IF(NW.EQ.(-2)) NZ=-2
135      RETURN
136   60 CONTINUE
137      IF (IFORM.EQ.2) GO TO 70
138C-----------------------------------------------------------------------
139C     ASYMPTOTIC EXPANSION FOR I(FNU,Z) FOR LARGE FNU APPLIED IN
140C     -PI/3.LE.ARG(Z).LE.PI/3
141C-----------------------------------------------------------------------
142      CALL CUNI1(Z, FNU, KODE, N, Y, NW, NLAST, FNUL, TOL, ELIM, ALIM)
143      GO TO 80
144   70 CONTINUE
145C-----------------------------------------------------------------------
146C     ASYMPTOTIC EXPANSION FOR J(FNU,Z*EXP(M*HPI)) FOR LARGE FNU
147C     APPLIED IN PI/3.LT.ABS(ARG(Z)).LE.PI/2 WHERE M=+I OR -I
148C     AND HPI=PI/2
149C-----------------------------------------------------------------------
150      CALL CUNI2(Z, FNU, KODE, N, Y, NW, NLAST, FNUL, TOL, ELIM, ALIM)
151   80 CONTINUE
152      IF (NW.LT.0) GO TO 50
153      NZ = NW
154      RETURN
155   90 CONTINUE
156      NLAST = N
157      RETURN
158      END
159