GNU Octave  6.2.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
d9lgmc.f
Go to the documentation of this file.
1 *DECK D9LGMC
2  DOUBLE PRECISION FUNCTION d9lgmc (X)
3 C***BEGIN PROLOGUE D9LGMC
4 C***SUBSIDIARY
5 C***PURPOSE Compute the log Gamma correction factor so that
6 C LOG(DGAMMA(X)) = LOG(SQRT(2*PI)) + (X-5.)*LOG(X) - X
7 C + D9LGMC(X).
8 C***LIBRARY SLATEC (FNLIB)
9 C***CATEGORY C7E
10 C***TYPE DOUBLE 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. so that
17 C LOG (DGAMMA(X)) = LOG(SQRT(2*PI)) + (X-.5)*LOG(X) - X + D9lGMC(X)
18 C
19 C Series for ALGM on the interval 0. to 1.00000E-02
20 C with weighted error 1.28E-31
21 C log weighted error 30.89
22 C significant figures required 29.81
23 C decimal places required 31.48
24 C
25 C***REFERENCES (NONE)
26 C***ROUTINES CALLED D1MACH, DCSEVL, INITDS, XERMSG
27 C***REVISION HISTORY (YYMMDD)
28 C 770601 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 D9LGMC
35  DOUBLE PRECISION x, algmcs(15), xbig, xmax, dcsevl, d1mach
36  LOGICAL first
37  SAVE algmcs, nalgm, xbig, xmax, first
38  DATA algmcs( 1) / +.1666389480 4518632472 0572965082 2 d+0 /
39  DATA algmcs( 2) / -.1384948176 0675638407 3298605913 5 d-4 /
40  DATA algmcs( 3) / +.9810825646 9247294261 5717154748 7 d-8 /
41  DATA algmcs( 4) / -.1809129475 5724941942 6330626671 9 d-10 /
42  DATA algmcs( 5) / +.6221098041 8926052271 2601554341 6 d-13 /
43  DATA algmcs( 6) / -.3399615005 4177219443 0333059966 6 d-15 /
44  DATA algmcs( 7) / +.2683181998 4826987489 5753884666 6 d-17 /
45  DATA algmcs( 8) / -.2868042435 3346432841 4462239999 9 d-19 /
46  DATA algmcs( 9) / +.3962837061 0464348036 7930666666 6 d-21 /
47  DATA algmcs( 10) / -.6831888753 9857668701 1199999999 9 d-23 /
48  DATA algmcs( 11) / +.1429227355 9424981475 7333333333 3 d-24 /
49  DATA algmcs( 12) / -.3547598158 1010705471 9999999999 9 d-26 /
50  DATA algmcs( 13) / +.1025680058 0104709120 0000000000 0 d-27 /
51  DATA algmcs( 14) / -.3401102254 3167487999 9999999999 9 d-29 /
52  DATA algmcs( 15) / +.1276642195 6300629333 3333333333 3 d-30 /
53  DATA first /.true./
54 C***FIRST EXECUTABLE STATEMENT D9LGMC
55  IF (first) THEN
56  nalgm = initds(algmcs, 15, real(d1mach(3)) )
57  xbig = 1.0d0/sqrt(d1mach(3))
58  xmax = exp(min(log(d1mach(2)/12.d0), -log(12.d0*d1mach(1))))
59  ENDIF
60  first = .false.
61 C
62  IF (x .LT. 10.d0) CALL xermsg ('SLATEC', 'D9LGMC',
63  + 'X MUST BE GE 10', 1, 2)
64  IF (x.GE.xmax) GO TO 20
65 C
66  d9lgmc = 1.d0/(12.d0*x)
67  IF (x.LT.xbig) d9lgmc = dcsevl(2.0d0*(10.d0/x)**2-1.d0, algmcs,
68  1 nalgm) / x
69  RETURN
70 C
71  20 d9lgmc = 0.d0
72  CALL xermsg ('SLATEC', 'D9LGMC', 'X SO BIG D9LGMC UNDERFLOWS', 2,
73  + 1)
74  RETURN
75 C
76  END
charNDArray min(char d, const charNDArray &m)
Definition: chNDArray.cc:207
double precision function d1mach(i)
Definition: d1mach.f:23
double precision function d9lgmc(X)
Definition: d9lgmc.f:3
ColumnVector real(const ComplexColumnVector &a)
Definition: dColVector.cc:137
double precision function dcsevl(X, CS, N)
Definition: dcsevl.f:3
function initds(OS, NOS, ETA)
Definition: initds.f:3
octave_int< T > xmax(const octave_int< T > &x, const octave_int< T > &y)
subroutine xermsg(LIBRAR, SUBROU, MESSG, NERR, LEVEL)
Definition: xermsg.f:3