r9lgit.f

Go to the documentation of this file.
00001 *DECK R9LGIT
00002       FUNCTION R9LGIT (A, X, ALGAP1)
00003 C***BEGIN PROLOGUE  R9LGIT
00004 C***SUBSIDIARY
00005 C***PURPOSE  Compute the logarithm of Tricomi's incomplete Gamma
00006 C            function with Perron's continued fraction for large X and
00007 C            A .GE. X.
00008 C***LIBRARY   SLATEC (FNLIB)
00009 C***CATEGORY  C7E
00010 C***TYPE      SINGLE PRECISION (R9LGIT-S, D9LGIT-D)
00011 C***KEYWORDS  FNLIB, INCOMPLETE GAMMA FUNCTION, LOGARITHM,
00012 C             PERRON'S CONTINUED FRACTION, SPECIAL FUNCTIONS, TRICOMI
00013 C***AUTHOR  Fullerton, W., (LANL)
00014 C***DESCRIPTION
00015 C
00016 C Compute the log of Tricomi's incomplete gamma function with Perron's
00017 C continued fraction for large X and for A .GE. X.
00018 C
00019 C***REFERENCES  (NONE)
00020 C***ROUTINES CALLED  R1MACH, XERMSG
00021 C***REVISION HISTORY  (YYMMDD)
00022 C   770701  DATE WRITTEN
00023 C   890531  Changed all specific intrinsics to generic.  (WRB)
00024 C   890531  REVISION DATE from Version 3.2
00025 C   891214  Prologue converted to Version 4.0 format.  (BAB)
00026 C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
00027 C   900720  Routine changed from user-callable to subsidiary.  (WRB)
00028 C***END PROLOGUE  R9LGIT
00029       SAVE EPS, SQEPS
00030       DATA EPS, SQEPS / 2*0.0 /
00031 C***FIRST EXECUTABLE STATEMENT  R9LGIT
00032       IF (EPS.EQ.0.0) EPS = 0.5*R1MACH(3)
00033       IF (SQEPS.EQ.0.0) SQEPS = SQRT(R1MACH(4))
00034 C
00035       IF (X .LE. 0.0 .OR. A .LT. X) CALL XERMSG ('SLATEC', 'R9LGIT',
00036      +   'X SHOULD BE GT 0.0 AND LE A', 2, 2)
00037 C
00038       AX = A + X
00039       A1X = AX + 1.0
00040       R = 0.0
00041       P = 1.0
00042       S = P
00043       DO 20 K=1,200
00044         FK = K
00045         T = (A+FK)*X*(1.0+R)
00046         R = T/((AX+FK)*(A1X+FK)-T)
00047         P = R*P
00048         S = S + P
00049         IF (ABS(P).LT.EPS*S) GO TO 30
00050  20   CONTINUE
00051       CALL XERMSG ('SLATEC', 'R9LGIT',
00052      +   'NO CONVERGENCE IN 200 TERMS OF CONTINUED FRACTION', 3, 2)
00053 C
00054  30   HSTAR = 1.0 - X*S/A1X
00055       IF (HSTAR .LT. SQEPS) CALL XERMSG ('SLATEC', 'R9LGIT',
00056      +   'RESULT LESS THAN HALF PRECISION', 1, 1)
00057 C
00058       R9LGIT = -X - ALGAP1 - LOG(HSTAR)
00059 C
00060       RETURN
00061       END
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Defines