Directory of image this file is from
This file as a plain text file

C
C     ..................................................................
C
C        SUBROUTINE CEL1
C
C        PURPOSE
C           CALCULATE COMPLETE ELLIPTIC INTEGRAL OF FIRST KIND
C
C        USAGE
C           CALL CEL1(RES,AK,IER)
C
C        DESCRIPTION OF PARAMETERS
C           RES   - RESULT VALUE
C           AK    - MODULUS (INPUT)
C           IER   - RESULTANT ERROR CODE WHERE
C                   IER=0  NO ERROR
C                   IER=1  AK NOT IN RANGE -1 TO +1
C
C        REMARKS
C           THE RESULT IS SET TO 1.7E38 IF ABS(AK) GE 1                         0
C           FOR MODULUS AK AND COMPLEMENTARY MODULUS CK,
C           EQUATION AK*AK+CK*CK=1.0 IS USED.
C           AK MUST BE IN THE RANGE -1 TO +1
C
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
C           NONE
C
C        METHOD
C           DEFINITION
C           CEL1(AK)=INTEGRAL(1/SQRT((1+T*T)*(1+(CK*T)**2)), SUMMED
C           OVER T FROM 0 TO INFINITY).
C           EQUIVALENT ARE THE DEFINITIONS
C           CEL1(AK)=INTEGRAL(1/(COS(T)SQRT(1+(CK*TAN(T))**2)),SUMMED
C           OVER T FROM 0 TO PI/2),
C           CEL1(AK)=INTEGRAL(1/SQRT(1-(AK*SIN(T))**2),SUMMED OVER T
C           FROM 0 TO PI/2), WHERE K=SQRT(1.-CK*CK).
C           EVALUATION
C           LANDENS TRANSFORMATION IS USED FOR CALCULATION.
C           REFERENCE
C           R.BULIRSCH, 'NUMERICAL CALCULATION OF ELLIPTIC INTEGRALS
C           AND ELLIPTIC FUNCTIONS', HANDBOOK SERIES SPECIAL FUNCTIONS,
C           NUMERISCHE MATHEMATIK VOL. 7, 1965, PP. 78-90.
C
C     ..................................................................
C
      SUBROUTINE CEL1(RES,AK,IER)
      IER=0
      ARI=2.
      GEO=(0.5-AK)+0.5
      GEO=GEO+GEO*AK
      RES=0.5
      IF(GEO)1,2,4
    1 IER=1
    2 RES=1.7E38                                                                0
      RETURN
    3 GEO=GEO*AARI
    4 GEO=SQRT(GEO)
      GEO=GEO+GEO
      AARI=ARI
      ARI=ARI+GEO
      RES=RES+RES
      IF(GEO/AARI-0.9999)3,5,5
    5 RES=RES/ARI*6.283185E0
      RETURN
      END