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
r9lgmc.f
Go to the documentation of this file.
1 *DECK R9LGMC
2  FUNCTION r9lgmc (X)
3 C***BEGIN PROLOGUE R9LGMC
4 C***SUBSIDIARY
5 C***PURPOSE Compute the log Gamma correction factor so that
6 C LOG(GAMMA(X)) = LOG(SQRT(2*PI)) + (X-.5)*LOG(X) - X
7 C + R9LGMC(X).
8 C***LIBRARY SLATEC (FNLIB)
9 C***CATEGORY C7E
10 C***TYPE SINGLE PRECISION (R9LGMC-S, D9LGMC-D, C9LGMC-C)
11 C***KEYWORDS COMPLETE GAMMA FUNCTION, CORRECTION TERM, FNLIB,
12 C LOG GAMMA, LOGARITHM, SPECIAL FUNCTIONS
13 C***AUTHOR Fullerton, W., (LANL)
14 C***DESCRIPTION
15 C
16 C Compute the log gamma correction factor for X .GE. 10.0 so that
17 C LOG (GAMMA(X)) = LOG(SQRT(2*PI)) + (X-.5)*LOG(X) - X + R9LGMC(X)
18 C
19 C Series for ALGM on the interval 0. to 1.00000D-02
20 C with weighted error 3.40E-16
21 C log weighted error 15.47
22 C significant figures required 14.39
23 C decimal places required 15.86
24 C
25 C***REFERENCES (NONE)
26 C***ROUTINES CALLED CSEVL, INITS, R1MACH, XERMSG
27 C***REVISION HISTORY (YYMMDD)
28 C 770801 DATE WRITTEN
29 C 890531 Changed all specific intrinsics to generic. (WRB)
30 C 890531 REVISION DATE from Version 3.2
31 C 891214 Prologue converted to Version 4.0 format. (BAB)
32 C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ)
33 C 900720 Routine changed from user-callable to subsidiary. (WRB)
34 C***END PROLOGUE R9LGMC
35  dimension algmcs(6)
36  LOGICAL first
37  SAVE algmcs, nalgm, xbig, xmax, first
38  DATA algmcs( 1) / .1666389480 45186e0 /
39  DATA algmcs( 2) / -.0000138494 817606e0 /
40  DATA algmcs( 3) / .0000000098 108256e0 /
41  DATA algmcs( 4) / -.0000000000 180912e0 /
42  DATA algmcs( 5) / .0000000000 000622e0 /
43  DATA algmcs( 6) / -.0000000000 000003e0 /
44  DATA first /.true./
45 C***FIRST EXECUTABLE STATEMENT R9LGMC
46  IF (first) THEN
47  nalgm = inits(algmcs, 6, r1mach(3))
48  xbig = 1.0/sqrt(r1mach(3))
49  xmax = exp(min(log(r1mach(2)/12.0), -log(12.0*r1mach(1))) )
50  ENDIF
51  first = .false.
52 C
53  IF (x .LT. 10.0) CALL xermsg('SLATEC', 'R9LGMC',
54  + 'X MUST BE GE 10', 1, 2)
55  IF (x.GE.xmax) go to 20
56 C
57  r9lgmc = 1.0/(12.0*x)
58  IF (x.LT.xbig) r9lgmc = csevl(2.0*(10./x)**2-1., algmcs, nalgm)/x
59  RETURN
60 C
61  20 r9lgmc = 0.0
62  CALL xermsg('SLATEC', 'R9LGMC', 'X SO BIG R9LGMC UNDERFLOWS', 2,
63  + 1)
64  RETURN
65 C
66  END