/*
 * stirling.c
 * This program checks Stirling's formula.
 */

#include <stdio.h>
#include <math.h>

#define PI 3.14159265
#define N  20

/* Prototypes: */
double log_fact (int n);

int main ()
{
  int    n;
  double exact_res;
  double stirling_res;

  for (n = 1; n <= N; n++)
  {
    /* Compute the exact value of ln(n!): */
    exact_res = log_fact(n);

    /*                                            n^(n+0.5)
     * Use Stirling's formula: n! ~ sqrt(2*PI) * -----------
     *                                              e^n
     *
     * Therefore: ln(n!) ~ 0.5*ln(2*PI) + (n+0.5)*ln(n) - n
     */
    stirling_res = 0.5*log(2*PI) + (n + 0.5)*log(n) - n;

    /* Print the results. */
    printf ("ln(%d!) = %.6lf\t", n, exact_res);
    printf ("Stirling's apporximation is %.6lf\n", stirling_res);
  }

  return (0);
}

/* ------------------------------------------------------------------------
 * Function: log_fact
 * Purpose : Compute the natural logarithm of a factorial.
 * Input   : n - The factorial base.
 * Returns : The value of ln(n!).
 */
double log_fact (int n)
{
  if (n <= 1)
    /* Since 0! = 1! = 1, we can simple return ln(1) = 0: */
    return (0);

  /* Use the equality: ln(n!) = ln((n-1)!) + ln(n) */
  return (log_fact(n - 1) + log (n));
}

