albeta.f

Go to the documentation of this file.

*DECK ALBETA
      FUNCTION albeta (A, B)
C***BEGIN PROLOGUE  ALBETA
C***PURPOSE  Compute the natural logarithm of the complete Beta
C            function.
C***LIBRARY   SLATEC (FNLIB)
C***CATEGORY  C7B
C***TYPE      SINGLE PRECISION (ALBETA-S, DLBETA-D, CLBETA-C)
C***KEYWORDS  FNLIB, LOGARITHM OF THE COMPLETE BETA FUNCTION,
C             SPECIAL FUNCTIONS
C***AUTHOR  Fullerton, W., (LANL)
C***DESCRIPTION
C
C ALBETA computes the natural log of the complete beta function.
C
C Input Parameters:
C       A   real and positive
C       B   real and positive
C
C***REFERENCES  (NONE)
C***ROUTINES CALLED  ALNGAM, ALNREL, GAMMA, R9LGMC, XERMSG
C***REVISION HISTORY  (YYMMDD)
C   770701  DATE WRITTEN
C   890531  Changed all specific intrinsics to generic.  (WRB)
C   890531  REVISION DATE from Version 3.2
C   891214  Prologue converted to Version 4.0 format.  (BAB)
C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
C   900326  Removed duplicate information from DESCRIPTION section.
C           (WRB)
C   900727  Added EXTERNAL statement.  (WRB)
C***END PROLOGUE  ALBETA
      EXTERNAL gamma
      SAVE sq2pil
      DATA sq2pil / 0.9189385332 0467274 e0 /
C***FIRST EXECUTABLE STATEMENT  ALBETA
      p = min(a, b)
      q = max(a, b)
C
      IF (p .LE. 0.0) CALL xermsg('SLATEC', 'ALBETA',
     +   'BOTH ARGUMENTS MUST BE GT ZERO', 1, 2)
      IF (p.GE.10.0) go to 30
      IF (q.GE.10.0) go to 20
C
C P AND Q ARE SMALL.
C
      albeta = log(gamma(p) * (gamma(q)/gamma(p+q)) )
      RETURN
C
C P IS SMALL, BUT Q IS BIG.
C
 20   corr = r9lgmc(q) - r9lgmc(p+q)
      albeta = alngam(p) + corr + p - p*log(p+q) +
     1  (q-0.5)*alnrel(-p/(p+q))
      RETURN
C
C P AND Q ARE BIG.
C
 30   corr = r9lgmc(p) + r9lgmc(q) - r9lgmc(p+q)
      albeta = -0.5*log(q) + sq2pil + corr + (p-0.5)*log(p/(p+q))
     1  + q*alnrel(-p/(p+q))
      RETURN
C
      END