Application of Continued Fractions to the Computation of the Gamma Function

V. E. Kornilov

(Presented by the Department of Higher Mathematics)

Along with the existing fractional-rational approximations for the gamma function [1], pp. 171–174, the present article studies a new type of approximation of this function.

In the article the following abbreviated notation is introduced:

\[ \left(\frac{x}{e}\right)^x \equiv z;\qquad (a)_\kappa = a(a+1)\ldots(a+\kappa-1),\qquad (a)_0=1. \]

1. The gamma function is usually represented in the form of the following integral:

\[ \Gamma(x)=\int_0^\infty e^{-t}t^{x-1}\,dt,\qquad x>0. \tag{1} \]

Or

\[ \Gamma(x)=\int_y^\infty e^{-t}t^{x-1}\,dt+\int_0^y e^{-t}t^{x-1}\,dt. \tag{2} \]

The integrals situated on the right-hand side of equality 2 are expanded into the following continued fractions [1], pp. 140, 145:

\[ \int_y^\infty e^{-t}t^{x-1}\,dt = \frac{y^x e^{-y}}{ y+\dfrac{1-x}{ 1+\ldots+\dfrac{n}{ y+\dfrac{n+1-x}{1+\ldots} } } }; \tag{3} \]

\[ \int_0^y e^{-t}t^{x-1}\,dt = \frac{y^x e^{-y}}{ x-\dfrac{yx}{ 1+x+\ldots+\dfrac{ny}{ 2n+x-\dfrac{(n+x)y}{ 2n+1+x+\ldots } } } }. \tag{4} \]

The gamma function can be approximately replaced by the sum of suitable fractions, which are computed according to equalities (3) and (4). By choos-

the continued fractions of the continued fractions (3) and (4) take their simplest form in the case when the variable limit of the integrals (3) and (4) is \(y=x\), since:

a) the convergents of the continued fraction (3) are then such that the greatest powers of \(x\) are half as large as in the case \(y\ne x\);

b) the convergents of the continued fraction (4) contain simpler denominators, although the numerators are of higher degree than in the case \(y\ne x\).

In equalities (3) and (4) we set \(y=x\) and successively compute the convergents, which are computed below from the 1st to the 20th convergent for integrals of the form (3) and from the 1st to the 12th convergent for integrals of the form (4):

\[ \left\{ \begin{aligned} \int_x^\infty e^{-t}t^{z-1}\,dt \approx{}& 1)\ \frac{z}{x};\quad 2)\ z;\quad 3)\ z\frac{x+1}{x};\quad 4)\ \frac{3z}{x+2}; \\[6pt] &5)\ z\frac{5x+2}{x(x+6)};\quad 6)\ z\frac{2x+11}{7x+6};\quad 7)\ z\frac{x^2+13x+3}{x(5x+12)}; \\[6pt] &8)\ z\frac{23x+50}{3x^2+46x+24};\quad 9)\ z\frac{31x^2+154x+24}{x(3x^2+86x+120)}; \\[6pt] &10)\ z\frac{8x^2+219x+274}{55x^2+326x+120};\quad 11)\ z\frac{4x^3+187x^2+522x+60}{x(35x^2+378x+360)}; \\[6pt] &12)\ \frac{z(203x^2+2084x+1764)}{15x^3+760x^2+2556x+720}; \\[6pt] &13)\ \frac{z(251x^3+4328x^2+8028x+720)}{x(15x^3+1180x^2+7092x+5040)}; \\[6pt] &14)\ z\frac{48x^3+3665x^2+20852x+13068}{525x^3+9856x^2+22212x+5040}; \\[6pt] &15)\ z\frac{24x^4+2711x^3+25574x^2+34632x+2520}{x(315x^3+9058x^2+35928x+20160)}; \\[6pt] &16)\ \frac{z(2141x^3+59616x^2+223012x+109584)}{105x^4+12460x^3+128492x^2+212976x+40320}; \\[6pt] &17)\ z\frac{505x^4+20598{,}4x^3+126439{,}2x^2+132739{,}2x+8064}{x(21x^4+3500x^3+54684x^2+157564{,}8x+72576)}; \\[6pt] &18)\ \frac{z(384x^4+62645x^3+945728x^2+2561220x+1026576)}{5985x^4+257068x^3+1731276x^2+2239344x+362880}; \\[6pt] &19)\ \frac{192x^5+42685x^4+936328x^3+4125492x^2+3499920x+181440}{z^{-1}(3465x^5+207284x^4+2096028x^3+4664880x^2+1814400x)}; \\[6pt] &20)\ \frac{z(5313x^4+310675{,}6x^3+3029408{,}8x^2+6317092{,}8x+2125728)}{189x^5+43470x^4+1006292x^3+4880635{,}2x^2+5131872x+725760}. \end{aligned} \right. \tag{5} \]

(convergents 7), 11), 15), 19) are reduced by 2, and fractions 17) and 20) are reduced by 5).

\[ \left\{ \begin{aligned} \int_{0}^{x} e^{-t} t^{x-1}\,dt \simeq {}& 1)\ \frac{z}{x};\qquad 2)\ z\frac{x+1}{x};\qquad 3)\ z\frac{x^{2}+4x+2}{2x(x+1)};\\[6pt] &4)\ z\frac{(x+2)(5x+3)}{(x)_2(x+6)};\qquad 5)\ z\frac{7x^{3}+41x^{2}+62x+24}{(x)_3(x+12)};\\[6pt] &6)\ z\frac{2x^{3}+47x^{2}+92x+40}{(x+3)^{-1}(x)_3(11x+60)};\\[6pt] &7)\ \frac{x^{5}+43x^{4}+338x^{3}+959x^{2}1062x+360}{z^{-1}(x)_4(7x+60)};\\[6pt] &8)\ \frac{41x^{4}+720x^{3}+2743x^{2}+3492x+1260}{z^{-1}(x+4)^{-1}(x)_4(3x^{2}+158x+840)};\\[6pt] &9)\ \frac{49x^{6}+1556x^{5}+15399x^{4}+67120x^{3}+139436x^{2}+129744x+40320}{z^{-1}(x)_5(3x^{2}+226x+1680)};\\[6pt] &10)\ \frac{8x^{6}+909x^{5}+15699x^{4}+90260x^{3}+219132x^{2}+223056x+72576}{z^{-1}(x+5)^{-1}(x)_5(95x^{2}+2874x+15120)};\\[6pt] &11)\ \frac{ \displaystyle \frac{x^{8}}{3780}+\frac{91x^{7}}{2160}+\frac{18757x^{6}}{15120}+\frac{44819x^{5}}{3024}+\frac{134679,1}{1512x^{-4}}+ \frac{109002,7}{378x^{-3}}+\frac{29809}{60x^{-2}}+\frac{8641}{21x^{-1}}+120 }{ \displaystyle z^{-1}15120^{-1}(x)_6\bigl(55x^{2}+2202x+15120\bigr) };\\[6pt] &12)\ \frac{ \displaystyle \frac{37,3x^{7}}{33264}+\frac{24301x^{6}}{332640}+\frac{445063x^{5}}{332640}+\frac{343908,7}{33264x^{-4}}+ \frac{93173}{2376x^{-3}}+\frac{622483,9}{8316x^{-2}}+\frac{91691}{1386x^{-1}}+20 }{ \displaystyle z^{-1}(x+6)^{-1}332640^{-1}(x)_6\bigl(15x^{3}+2740x^{2}+63564x+332640\bigr) }; \end{aligned} \right. \tag{6} \]

(the convergents 7) and 11) have been reduced by 2).

The gamma function, according to equality (2), is approximately replaced by the sum of the fractional-rational approximations (5) and (6). The most advantageous choice of approximations consists in making the errors of these approximations approximately equal.

1. Taking into account the circumstance that for real \(x>0\) the remainder terms \(zR_m(x)\) and \(zr_k(x)\) for the convergents (5) and (6), respectively, can be represented in the form of infinite functional series of Leibniz type, namely ([1], pp. 11, 12):

\[ R_m(x)=\sum_{n=m}^{\infty}(-1)^n \frac{(1-x)^{\left[\frac{n+1}{2}\right]}(1)^{\left[\frac{n}{2}\right]}} {Q_n(x)Q_{n+1}(x)}, \tag{7} \]

\[ r_\kappa(x)=\sum_{n=0}^{\infty}(-1)^{\left[\frac{\kappa}{2}\right]+n} \frac{x^{\kappa+2n}(x)\left[\frac{\kappa+1}{2}\right]+n(1)\left[\frac{\kappa}{2}\right]+n} {(\kappa+1+2n+x)^{-1}q_{\kappa+2n}(x)q_{\kappa+2+2n}(x)}, \tag{8} \]

where \([v]\) denotes the greatest integer not exceeding \(v\); it is not difficult to conclude that the remainder terms, in absolute value, are less than the first terms of the right-hand sides of equalities (7) and (8). That is,

\[ |R_m(x)|< \frac{|(1-x)|^{\left[\frac{m+1}{2}\right]} |(1)|^{\left[\frac{m}{2}\right]}} {Q_m(x)Q_{m+1}(x)},\qquad \left[\frac{m+1}{2}\right]>x>0. \tag{9} \]

\[ |r_\kappa(x)|< \frac{(\kappa+1+x)x^\kappa(x)\left[\frac{\kappa+1}{2}\right](1)\left[\frac{\kappa}{2}\right]} {q_\kappa(x)q_{\kappa+2}(x)},\qquad x>0. \tag{10} \]

On the basis of the foregoing, the gamma function is represented by the following equality:

\[ \Gamma(x)=z\left[\frac{P_m(x)}{Q_m(x)}+\frac{p_\kappa(x)}{q_\kappa(x)}+R_m(x)+r_\kappa(x)\right], \tag{11} \]

where the terms in square brackets are determined, respectively, according to equalities (5), (6) and inequalities (9) and (10).

Let us determine an estimate of the remainder term of the gamma function if, for its computation, we use the 19th convergent from equalities (5) and the 11th convergent from equalities (6) on the interval \(2\le x\le 3\); at the ends of the interval \(R_{19}(2)=R_{19}(3)=0\). We denote the remainder term of the gamma function by \(\rho_{15}(x)\):

\[ |\rho_{15}(x)|\le z\left(|R_{19}(x)|+|r_{11}(x)|\right). \tag{12} \]

After computing \(|\rho_{15}(x)|\) with step \(0.1\), we have

\[ |\rho_{15}(2)|<0.00000000133=0.0^8 133, \]

\[ \begin{aligned} |\rho_{15}(2.1)|&<0.0^7 364, & |\rho_{15}(2.6)|&<0.0^7 425,\\ |\rho_{15}(2.2)|&<0.0^7 516, & |\rho_{15}(2.7)|&<0.0^7 393,\\ |\rho_{15}(2.3)|&<0.0^7 553, & |\rho_{15}(2.8)|&<0.0^7 392,\\ |\rho_{15}(2.4)|&<0.0^7 529, & |\rho_{15}(2.9)|&<0.0^7 412,\\ |\rho_{15}(2.5)|&<0.0^7 475, & |\rho_{15}(3)|&<0.0^7 530. \end{aligned} \tag{13} \]

Applying formula (11), we compute the gamma function for \(x=2.5\):

\[ \Gamma(2.5)\simeq \sqrt{\frac{2.5^3}{e^5}} \left( \frac{163301833}{95843005}+ \frac{43399497}{18135117} \right) = \]

\[ =0.32446944673\cdot 4.0969663890=1.3293404175. \]

Bearing in mind that

\[ \Gamma(2.5)=0.75\sqrt{\pi}=1.3293403882, \]

we obtain \(|\rho_{15}(2.5)|<0.0^7 293\).

Conclusion

Using the 19th fractional-rational approximation of integral (3) and the 11th approximation of integral (4) on the interval \([2,3]\), one can compute the gamma function with relative error \(\delta_x<0.0^7 6\) ([2], p. 115).

References

1. A. N. Khovanskii. Application of continued fractions and their generalizations to problems of approximate analysis. Moscow, Gostekhizdat, 1956.

2. I. N. Bronshtein, K. A. Semendyaev. Handbook of Mathematics. Moscow, GITTL, 1956.