1*DECK FAC
2      FUNCTION FAC (N)
3C***BEGIN PROLOGUE  FAC
4C***PURPOSE  Compute the factorial function.
5C***LIBRARY   SLATEC (FNLIB)
6C***CATEGORY  C1
7C***TYPE      SINGLE PRECISION (FAC-S, DFAC-D)
8C***KEYWORDS  FACTORIAL, FNLIB, SPECIAL FUNCTIONS
9C***AUTHOR  Fullerton, W., (LANL)
10C***DESCRIPTION
11C
12C FAC(N) evaluates the factorial function of N.  FAC is single
13C precision.  N must be an integer between 0 and 25 inclusive.
14C
15C***REFERENCES  (NONE)
16C***ROUTINES CALLED  GAMLIM, R9LGMC, XERMSG
17C***REVISION HISTORY  (YYMMDD)
18C   770601  DATE WRITTEN
19C   890531  Changed all specific intrinsics to generic.  (WRB)
20C   890531  REVISION DATE from Version 3.2
21C   891214  Prologue converted to Version 4.0 format.  (BAB)
22C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
23C***END PROLOGUE  FAC
24      DIMENSION FACN(26)
25      SAVE FACN, SQ2PIL, NMAX
26      DATA FACN( 1) / 1.0E0 /
27      DATA FACN( 2) / 1.0E0 /
28      DATA FACN( 3) / 2.0E0 /
29      DATA FACN( 4) / 6.0E0 /
30      DATA FACN( 5) / 24.0E0 /
31      DATA FACN( 6) / 120.0E0 /
32      DATA FACN( 7) / 720.0E0 /
33      DATA FACN( 8) / 5040.0E0 /
34      DATA FACN( 9) / 40320.0E0 /
35      DATA FACN(10) / 362880.0E0 /
36      DATA FACN(11) / 3628800.0E0 /
37      DATA FACN(12) / 39916800.0E0 /
38      DATA FACN(13) / 479001600.0E0 /
39      DATA FACN(14) / 6227020800.0E0 /
40      DATA FACN(15) / 87178291200.0E0 /
41      DATA FACN(16) / 1307674368000.0E0 /
42      DATA FACN(17) / 20922789888000.0E0 /
43      DATA FACN(18) / 355687428096000.0E0 /
44      DATA FACN(19) / 6402373705728000.0E0 /
45      DATA FACN(20) /  .12164510040883200E18 /
46      DATA FACN(21) /  .24329020081766400E19 /
47      DATA FACN(22) /  .51090942171709440E20 /
48      DATA FACN(23) /  .11240007277776077E22 /
49      DATA FACN(24) /  .25852016738884977E23 /
50      DATA FACN(25) /  .62044840173323944E24 /
51      DATA FACN(26) /  .15511210043330986E26 /
52      DATA SQ2PIL / 0.9189385332 0467274E0/
53      DATA NMAX / 0 /
54C***FIRST EXECUTABLE STATEMENT  FAC
55      IF (NMAX.NE.0) GO TO 10
56      CALL GAMLIM (XMIN, XMAX)
57      NMAX = XMAX - 1.
58C
59 10   IF (N .LT. 0) CALL XERMSG ('SLATEC', 'FAC',
60     +   'FACTORIAL OF NEGATIVE INTEGER UNDEFINED', 1, 2)
61C
62      IF (N.LE.25) FAC = FACN(N+1)
63      IF (N.LE.25) RETURN
64C
65      IF (N .GT. NMAX) CALL XERMSG ('SLATEC', 'FAC',
66     +   'N SO BIG FACTORIAL(N) OVERFLOWS', 2, 2)
67C
68      X = N + 1
69      FAC = EXP ( (X-0.5)*LOG(X) - X + SQ2PIL + R9LGMC(X) )
70C
71      RETURN
72      END
73