GNU Octave  3.8.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
dlbeta.f
Go to the documentation of this file.
1 *DECK DLBETA
2  DOUBLE PRECISION FUNCTION dlbeta (A, B)
3 C***BEGIN PROLOGUE DLBETA
4 C***PURPOSE Compute the natural logarithm of the complete Beta
5 C function.
6 C***LIBRARY SLATEC (FNLIB)
7 C***CATEGORY C7B
8 C***TYPE DOUBLE PRECISION (ALBETA-S, DLBETA-D, CLBETA-C)
9 C***KEYWORDS FNLIB, LOGARITHM OF THE COMPLETE BETA FUNCTION,
10 C SPECIAL FUNCTIONS
11 C***AUTHOR Fullerton, W., (LANL)
12 C***DESCRIPTION
13 C
14 C DLBETA(A,B) calculates the double precision natural logarithm of
15 C the complete beta function for double precision arguments
16 C A and B.
17 C
18 C***REFERENCES (NONE)
19 C***ROUTINES CALLED D9LGMC, DGAMMA, DLNGAM, DLNREL, XERMSG
20 C***REVISION HISTORY (YYMMDD)
21 C 770701 DATE WRITTEN
22 C 890531 Changed all specific intrinsics to generic. (WRB)
23 C 890531 REVISION DATE from Version 3.2
24 C 891214 Prologue converted to Version 4.0 format. (BAB)
25 C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ)
26 C 900727 Added EXTERNAL statement. (WRB)
27 C***END PROLOGUE DLBETA
28  DOUBLE PRECISION a, b, p, q, corr, sq2pil, d9lgmc, dgamma, dlngam,
29  1 dlnrel
30  EXTERNAL dgamma
31  SAVE sq2pil
32  DATA sq2pil / 0.9189385332 0467274178 0329736405 62 d0 /
33 C***FIRST EXECUTABLE STATEMENT DLBETA
34  p = min(a, b)
35  q = max(a, b)
36 C
37  IF (p .LE. 0.d0) CALL xermsg('SLATEC', 'DLBETA',
38  + 'BOTH ARGUMENTS MUST BE GT ZERO', 1, 2)
39 C
40  IF (p.GE.10.d0) go to 30
41  IF (q.GE.10.d0) go to 20
42 C
43 C P AND Q ARE SMALL.
44 C
45  dlbeta = log(dgamma(p) * (dgamma(q)/dgamma(p+q)) )
46  RETURN
47 C
48 C P IS SMALL, BUT Q IS BIG.
49 C
50  20 corr = d9lgmc(q) - d9lgmc(p+q)
51  dlbeta = dlngam(p) + corr + p - p*log(p+q)
52  1 + (q-0.5d0)*dlnrel(-p/(p+q))
53  RETURN
54 C
55 C P AND Q ARE BIG.
56 C
57  30 corr = d9lgmc(p) + d9lgmc(q) - d9lgmc(p+q)
58  dlbeta = -0.5d0*log(q) + sq2pil + corr + (p-0.5d0)*log(p/(p+q))
59  1 + q*dlnrel(-p/(p+q))
60  RETURN
61 C
62  END