The Gamma Function

The techniques used to create logarithm and hyperbolic trig methods can be applied to model more complicated mathematical functions as well. As an example, let's write a method that will compute the gamma function. The gamma function is one of the higher transcendental functions and can be used to solve certain types of second-order differential equations. It is defined by the integral expression in Eq. (17.7).

Equation 17.7

An interesting quality of the gamma function is it satisfies the following recurrence expression (Eq. (17.8)).

Equation 17.8

To numerically solve for the gamma function, an approximate form must be derived for Eq. (17.7). A number of methods to approximate the gamma function have been developed over the years . The approximate form we will use was derived by Lanczos ¹ and approximates the gamma function according to Eq. (17.9).

Equation 17.9

Lanczos found that the error of this approximation was minimized when g = 5 and N = 6. The implementation of this method is actually quite simple. The c _i coefficients and the constant are defined as static members of the Math2 class. The gamma() method only requires one call each to the Math.pow() and Math.exp() methods. The code listing for the gamma() method is ”

 public static double gamma(double z) {   double grp1, grp2, grp3;   int j;   grp1 = z + 0.5;   grp2 = z + 5.5;   grp3 = c[0];   for(j=1; j<7; ++j) {     grp3 += c[j]/(z+j);   }   return Math.pow(grp2,grp1)*          Math.exp(-grp2)*TwoPi*grp3/z; } 

Despite the complicated look of Eq. (17.9), the gamma() function requires relatively few lines of code. Since the gamma() function increases as a function of ( z “ 1)! most computers will experience an arithmetic overflow at z > 140 or so. For this reason, gamma functions are sometimes written to return the logarithm of the gamma function rather than the gamma function itself.