PYTHIA8/pythia8170/src/PythiaStdlib.cxx

   1 // PythiaStdlib.cc is a part of the PYTHIA event generator.
   2 // Copyright (C) 2012 Torbjorn Sjostrand.
   3 // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
   4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
   5
   6 // Function definitions (not found in the header) for the Gamma function.
   7
   8 #include "PythiaStdlib.h"
   9
  10 namespace Pythia8 {
  11
  12 //==========================================================================
  13
  14 // The Gamma function for real arguments, using the Lanczos approximation.
  15 // Code based on http://en.wikipedia.org/wiki/Lanczos_approximation
  16
  17 double GammaCoef[9] = {
  18      0.99999999999980993,     676.5203681218851,   -1259.1392167224028,
  19       771.32342877765313,   -176.61502916214059,    12.507343278686905,
  20     -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7};
  21
  22 double GammaReal(double x) {
  23
  24   // Reflection formula (recursive!) for x < 0.5.
  25   if (x < 0.5) return M_PI / (sin(M_PI * x) * GammaReal(1 - x));
  26
  27   // Iterate through terms.
  28   double z = x - 1.;
  29   double gamma = GammaCoef[0];
  30   for (int i = 1; i < 9; ++i) gamma += GammaCoef[i] / (z + i);
  31
  32   // Answer.
  33   double t = z + 7.5;
  34   gamma *= sqrt(2. * M_PI) * pow(t, z + 0.5) * exp(-t);
  35   return gamma;
  36
  37 }
  38
  39 //==========================================================================
  40
  41 } // end namespace Pythia8