GNU Octave: libcruft/amos/cacai.f Source File

Go to the documentation of this file.
00001       SUBROUTINE CACAI(Z, FNU, KODE, MR, N, Y, NZ, RL, TOL, ELIM, ALIM)
00002 C***BEGIN PROLOGUE  CACAI
00003 C***REFER TO  CAIRY
00004 C
00005 C     CACAI APPLIES THE ANALYTIC CONTINUATION FORMULA
00006 C
00007 C         K(FNU,ZN*EXP(MP))=K(FNU,ZN)*EXP(-MP*FNU) - MP*I(FNU,ZN)
00008 C                 MP=PI*MR*CMPLX(0.0,1.0)
00009 C
00010 C     TO CONTINUE THE K FUNCTION FROM THE RIGHT HALF TO THE LEFT
00011 C     HALF Z PLANE FOR USE WITH CAIRY WHERE FNU=1/3 OR 2/3 AND N=1.
00012 C     CACAI IS THE SAME AS CACON WITH THE PARTS FOR LARGER ORDERS AND
00013 C     RECURRENCE REMOVED. A RECURSIVE CALL TO CACON CAN RESULT IF CACON
00014 C     IS CALLED FROM CAIRY.
00015 C
00016 C***ROUTINES CALLED  CASYI,CBKNU,CMLRI,CSERI,CS1S2,R1MACH
00017 C***END PROLOGUE  CACAI
00018       COMPLEX CSGN, CSPN, C1, C2, Y, Z, ZN, CY
00019       REAL ALIM, ARG, ASCLE, AZ, CPN, DFNU, ELIM, FMR, FNU, PI, RL,
00020      * SGN, SPN, TOL, YY, R1MACH
00021       INTEGER INU, IUF, KODE, MR, N, NN, NW, NZ
00022       DIMENSION Y(N), CY(2)
00023       DATA PI / 3.14159265358979324E0 /
00024       NZ = 0
00025       ZN = -Z
00026       AZ = CABS(Z)
00027       NN = N
00028       DFNU = FNU + FLOAT(N-1)
00029       IF (AZ.LE.2.0E0) GO TO 10
00030       IF (AZ*AZ*0.25E0.GT.DFNU+1.0E0) GO TO 20
00031    10 CONTINUE
00032 C-----------------------------------------------------------------------
00033 C     POWER SERIES FOR THE I FUNCTION
00034 C-----------------------------------------------------------------------
00035       CALL CSERI(ZN, FNU, KODE, NN, Y, NW, TOL, ELIM, ALIM)
00036       GO TO 40
00037    20 CONTINUE
00038       IF (AZ.LT.RL) GO TO 30
00039 C-----------------------------------------------------------------------
00040 C     ASYMPTOTIC EXPANSION FOR LARGE Z FOR THE I FUNCTION
00041 C-----------------------------------------------------------------------
00042       CALL CASYI(ZN, FNU, KODE, NN, Y, NW, RL, TOL, ELIM, ALIM)
00043       IF (NW.LT.0) GO TO 70
00044       GO TO 40
00045    30 CONTINUE
00046 C-----------------------------------------------------------------------
00047 C     MILLER ALGORITHM NORMALIZED BY THE SERIES FOR THE I FUNCTION
00048 C-----------------------------------------------------------------------
00049       CALL CMLRI(ZN, FNU, KODE, NN, Y, NW, TOL)
00050       IF(NW.LT.0) GO TO 70
00051    40 CONTINUE
00052 C-----------------------------------------------------------------------
00053 C     ANALYTIC CONTINUATION TO THE LEFT HALF PLANE FOR THE K FUNCTION
00054 C-----------------------------------------------------------------------
00055       CALL CBKNU(ZN, FNU, KODE, 1, CY, NW, TOL, ELIM, ALIM)
00056       IF (NW.NE.0) GO TO 70
00057       FMR = FLOAT(MR)
00058       SGN = -SIGN(PI,FMR)
00059       CSGN = CMPLX(0.0E0,SGN)
00060       IF (KODE.EQ.1) GO TO 50
00061       YY = -AIMAG(ZN)
00062       CPN = COS(YY)
00063       SPN = SIN(YY)
00064       CSGN = CSGN*CMPLX(CPN,SPN)
00065    50 CONTINUE
00066 C-----------------------------------------------------------------------
00067 C     CALCULATE CSPN=EXP(FNU*PI*I) TO MINIMIZE LOSSES OF SIGNIFICANCE
00068 C     WHEN FNU IS LARGE
00069 C-----------------------------------------------------------------------
00070       INU = INT(FNU)
00071       ARG = (FNU-FLOAT(INU))*SGN
00072       CPN = COS(ARG)
00073       SPN = SIN(ARG)
00074       CSPN = CMPLX(CPN,SPN)
00075       IF (MOD(INU,2).EQ.1) CSPN = -CSPN
00076       C1 = CY(1)
00077       C2 = Y(1)
00078       IF (KODE.EQ.1) GO TO 60
00079       IUF = 0
00080       ASCLE = 1.0E+3*R1MACH(1)/TOL
00081       CALL CS1S2(ZN, C1, C2, NW, ASCLE, ALIM, IUF)
00082       NZ = NZ + NW
00083    60 CONTINUE
00084       Y(1) = CSPN*C1 + CSGN*C2
00085       RETURN
00086    70 CONTINUE
00087       NZ = -1
00088       IF(NW.EQ.(-2)) NZ=-2
00089       RETURN
00090       END