Gamma special function

[Jim-215-Fisher]

During implementing the gamma probability distribution(#278), I realized there is a need for gamma function and log_gamma with complex argument as well as incomplete gamma functions. So here is the PR for gamma special function. The gamma and log_gamma functions with complex argument was coded for single and double precision. All implemented functions were tested against available multi-precision package version with satisfactory results. This PR includes the following:

1. an extension of gamma function to integer/complex argument,
2. an extension of log_gamma function to integer/complex argument,
3. a log_factorial function for integer argument,
4. lower/upper incomplete gamma functions,
5. normalized lower/upper incomplete gamma functions (P and Q).

Please test it and provide any comment and suggestion.

[jvdp1]

Thank you. It looks already good to me, while I am not a user of such functions. A main question is why not all functions support all kinds of reals.

[jvdp1]

Explain first the aim of this function.

[jvdp1]

I would suggest to move these sentences to the section Returned value.

[Jim-215-Fisher]

Thank you. It looks already good to me, while I am not a user of such functions. A main question is why not all functions support all kinds of reals.

Very good question. Short answer is limits of precision. Detailed answer is as follows. During the calculation of most special functions such as Gamma, Bessel functions, there are no finite analytical formulae to exactly calculate the function values. Those functions have to be approximated by using various techniques, which introducing truncation error and rounding error. These errors are more pronouncing during complex value multiplication and division arithmetic of the special functions. Because of the propagation of these errors, the calculated results will have much less significant digits or even incorrect value if the same precision as calculated result is used during calculation. By using higher precision during calculation, the desired precision digits are guaranteed to be correct. So single precision uses double precision for calculation, double precision uses quadruple precision for calculation. Extended precision calculation may need higher precision a little over quadruple, while quadruple precision calculation need much higher precision not available right now. If there is a package for arbitrary precision, all kinds of reals can be supported.

[jvdp1]

LGTM (only a minor comment). Thank you.

[milancurcic]

Thank you!