Back
 APM  Vol.11 No.1 , January 2021
Evaluation of Exponential Integral by Means of Fast-Converging Power Series
Abstract: Exponential integral for real arguments is evaluated by employing a fast-converging power series originally developed for the resolution of Grandi’s paradox. Laguerre’s historic solution is first recapitulated and then the new solution method is described in detail. Numerical results obtained from the present series solution are compared with the tabulated values correct to nine decimal places. Finally, comments are made for the further use of the present approach for integrals involving definite functions in denominator.

1. Introduction

Exponential integral is encountered in various physics problems such as heat transfer, scattering, and neutron transport processes. Standard solutions to the integral are typically given in series forms with restrictions on the magnitudes of arguments or asymptotic expansions [1]. Historically, Laguerre (1834-1886) gave a continued fraction expansion for its solution in [2], which was reprinted in [3]. Numerical treatment of a generalized exponential integral based on rational minimax approximations and a comprehensive account of the relevant literature can be found in Milgram [4]. Based on fast-converging series expansions the present work constructs a solution uniformly valid for the entire domain of positive real arguments.

2. Laguerre’s Solution

Laguerre [2] [3] investigated the integral x t 1 e t d t and derived the following continued fraction representation as a solution.

e x x + 1 1 x + 3 1 x + 5 4 1 4 x + 7 9 1 9 x + 9 16 1 16 x + (1)

which gives the following expressions according to the orders kept

e x x + 1 , ( x + 3 ) e x x 2 + 4 x + 2 , ( x 2 + 8 x + 11 ) e x x 3 + 9 x 2 + 18 x + 6 , ( x 3 + 15 x 2 + 58 x + 50 ) e x x 4 + 16 x 3 + 72 x 2 + 96 x + 24 . (2)

Setting the lower limit of integration x = 1 in the above expressions gives respectively,

e 1 2 0.18394, 4 e 1 7 0.21022, 20 e 1 34 0.21640, 124 e 1 209 0.21826, (3)

which obviously converges to the nine decimal place correct value 0.219383934 given in the tables in [ [1], p. 240]. Laguerre also gave the second-order approximation result for x = 4 as 7 / 34 e 4 0.003770866 , which is accurate to five decimal places when compared to the nine decimal place correct value 0.003779352. For x = 4 the present approach described in §3 attains the five decimal place accuracy with n = 9 terms while the nine decimal place accuracy requires n = 27 terms.

3. Integration by Use of a Novel Series Expansion

The exponential integral for real arguments is defined as [ [1], p. 262]

Γ ( 0, x ) = E 1 ( x ) = x t 1 e t d t (4)

where Γ ( 0, x ) is the incomplete gamma function and x > 0 is the argument or arbitrary lower limit of integration. By introducing a change of variable in the form t = x + u the integral in (4) transforms to

E 1 ( x ) = 0 e ( x + u ) ( x + u ) d u = e x x [ 0 x e u ( 1 + u / x ) d u + x e u ( 1 + u / x ) d u ] (5)

Evaluations of the above two-part integrals are carried out separately by employing relatively fast-converging series expansions developed in resolving Grandi’s paradox [5]. Detailed derivation of the rectified and convergence-improved expansions of 1 / ( 1 + u / x ) for u / x 1 and u / x 1 can be found in Beji [5]. Accordingly, when u / x 1 or u x

1 1 + u / x = 1 + ( 1 ) 1 2 n [ ( n 1 ) + ( n 2 ) + + ( n n 1 ) + ( n n ) ] ( u x ) + ( 1 ) 2 2 n [ ( n 2 ) + ( n 3 ) + + ( n n 1 ) + ( n n ) ] ( u x ) 2 + + ( 1 ) n 1 2 n [ ( n n 1 ) + ( n n ) ] ( u x ) n 1 + ( 1 ) n 2 n ( n n ) ( u x ) n = 1 + 1 2 n p = 1 n [ m = p n ( n m ) ] ( 1 ) p ( u x ) p = C 0 + p = 1 n C p ( u x ) p (6)

For u / x 1 or u x Equation (6) is manipulated into the following convergent form [5]:

1 1 + u / x = x u + ( 1 ) 1 2 n [ ( n 1 ) + ( n 2 ) + + ( n n 1 ) + ( n n ) ] ( x u ) 2 + ( 1 ) 2 2 n [ ( n 2 ) + ( n 3 ) + + ( n n 1 ) + ( n n ) ] ( x u ) 3 + + ( 1 ) n 1 2 n [ ( n n 1 ) + ( n n ) ] ( x u ) n + ( 1 ) n 2 n ( n n ) ( x u ) n + 1 = x u + 1 2 n p = 1 n [ m = p n ( n m ) ] ( 1 ) p ( x u ) p + 1 = C 0 ( x u ) + p = 1 n C p ( x u ) p + 1 (7)

where the coefficients for both (6) and (7) are defined as

C 0 = 1 and C p = ( 1 ) p 2 n m = p n ( n m ) for p = 1 , 2 , , n (8)

We remark that unlike the standard expansion 1 / ( 1 + u / x ) = 1 ( u / x ) + ( u / x ) 2 ( u / x ) 3 + , which is unsatisfactory in many aspects, Expansions (6) and (7) are quite advantageous in terms of accurate representation and convergence. The most remarkable feature of these series is the change of coefficients in such a way that the series provides a highly accurate expansion of the generating function for the selected truncation order. Integrals involving similar but complex functions can be evaluated as well by employing the extended forms of the above type series as presented in [6].

The series given in Equation (6) is used in the first part of the integration in (5) where u x .

0 x e u ( 1 + u / x ) d u = 0 x [ C 0 + p = 1 n C p ( u x ) p ] e u d u (9)

Integrating by parts it can be shown that for an arbitrary p 1 value

C p 0 x ( u x ) p e u d u = C p p ! x p [ ( u p p ! + u p 1 ( p 1 ) ! + u p 2 ( p 2 ) ! + + u 2 2 ! + u 1 ! + 1 0 ! ) e u ] 0 x

(10)

Establishing p = 1,2, , n 1, n separately and then adding all the terms give

0 x [ C 0 + p = 1 n C p ( u x ) p ] e u d u = [ ( C 0 + 1 ! 0 ! C 1 x + 2 ! 0 ! C 2 x 2 + 3 ! 0 ! C 3 x 3 + + n ! 0 ! C n x n ) e u ] 0 x [ ( 1 ! 1 ! C 1 x + 2 ! 1 ! C 2 x 2 + 3 ! 1 ! C 3 x 3 + + n ! 1 ! C n x n ) u e u ] 0 x [ ( 2 ! 2 ! C 2 x 2 + 3 ! 2 ! C 3 x 3 + 4 ! 2 ! C 4 x 4 + + n ! 2 ! C n x n ) u 2 e u ] 0 x [ ( n ! n ! C n x n ) u n e u ] 0 x (11)

Equation (11) can be expressed in a compact form as double summations.

0 x [ C 0 + p = 1 n C p ( u x ) p ] e u d u = [ p = 0 n ( m = p n m ! p ! C m x m ) u p e u ] 0 x (12)

The first part of integration in (5) has thus been evaluated in terms of the series introduced in (6). The second part of the integration is carried out by employing (7), which is convergent for u x .

x e u ( 1 + u / x ) d u = x [ C 0 ( x u ) + p = 1 n C p ( x u ) p + 1 ] e u d u (13)

Again integrating by parts successively we obtain for an arbitrary p 1 value

C p x p + 1 x e u u p + 1 d u = C p x p + 1 p ! [ ( ( p 1 ) ! u p + ( p 2 ) ! u p 1 + ( 1 ) p ( p p ) ! u ) e u ] x + ( 1 ) p C p x p + 1 p ! x e u u d u (14)

Starting from p = 0 , which simply produces C 0 x x ( e u / u ) d u , we add the terms evaluated individually according to (14) for p = 1 , 2 , , n and obtain

x [ C 0 ( x u ) + p = 1 n C p ( x u ) p + 1 ] e u d u = [ 1 0 ! C 0 x 1 1 ! C 1 x 2 + 1 2 ! C 2 x 3 1 3 ! C 3 x 4 + + ( 1 ) n 1 n ! C n x n + 1 ] x e u u d u + [ ( 0 ! 1 ! C 1 x 2 + 0 ! 2 ! C 2 x 3 0 ! 3 ! C 3 x 4 + + ( 1 ) n 0 ! n ! C n x n + 1 ) e u u ] x + [ ( 1 ! 2 ! C 2 x 3 + 1 ! 3 ! C 3 x 4 1 ! 4 ! C 4 x 5 + + ( 1 ) n 1 1 ! n ! C n x n + 1 ) e u u 2 ] x + [ ( 2 ! 3 ! C 3 x 4 + 2 ! 4 ! C 4 x 5 2 ! 5 ! C 5 x 6 + + ( 1 ) n 2 2 ! n ! C n x n + 1 ) e u u 3 ] x + + [ ( ( 1 ) n n + 1 ( n 1 ) ! n ! C n x n + 1 ) e u u n ] x (15)

Noting that the integral on the right is E 1 ( x ) with integration variable u and expressing the summations in compact forms give

x [ C 0 ( x u ) + p = 1 n C p ( x u ) p + 1 ] e u d u = [ p = 0 n ( 1 ) p p ! C p x p + 1 ] E 1 ( x ) + [ p = 1 n ( m = p n ( 1 ) m p + 1 ( p 1 ) ! m ! C m x m + 1 ) e u u p ] x (16)

We now consider Equation (5) and make use of (12) and (16) so that

E 1 ( x ) = e x x [ 0 x e u ( 1 + u / x ) d u + x e u ( 1 + u / x ) d u ] = e x x [ p = 0 n ( m = p n m ! p ! C m x m ) u p e u ] 0 x + e x x [ p = 0 n ( 1 ) p p ! C p x p + 1 ] E 1 ( x ) + e x x [ p = 1 n ( m = p n ( 1 ) m p + 1 ( p 1 ) ! m ! C m x m + 1 ) e u u p ] x (17)

Evaluating the integral limits, rearranging, introducing a new running index k = m p + 1 , and gathering the terms multiplying E 1 ( x ) together give

[ x e x p = 0 n ( 1 ) p p ! C p x p + 1 ] E 1 ( x ) = ( 1 e x ) p = 0 n p ! C p x p e x p = 1 n [ m = p n [ ( 1 ) k ( m k ) ! m ! x k + ( x p ) m ! ( m k ) ! x k ] C m ] (18)

which in turn yields the series solution uniformly valid in the entire range of 0 < x < for the exponential integral E 1 ( x ) = x t 1 e t d t .

4. Computed Results

Equation (18) is now used for sample computations made for different arguments or values of the lower limit of the integration x. Figure 1 shows the plot of E 1 ( x ) versus x as computed from Equation (18) for the range 0.25 x 5.00 . Table 1 compares the numerically computed and tabulated values taken from Abramowitz and Stegun [ [1], pp. 238-242] correct for nine decimal places. The relative error percentage is computed as 100 | 1 E 1 c o m p . / E 1 e x a c t | . Obviously, the error decreases with increasing x and truncation order n. The latter is expected and normally increasing the truncation order for smaller x values decreases the error percentage but this was not possible due to computational capability limitations. For smaller x values the numbers computed in Equation (18) become large hence n cannot be increased further. For instance, for x = 0.25 the maximum possible n is 10 and the corresponding relative error percentage is 1.56%1. This point is the only drawback of the present approach. Machines that can handle big numbers would have no such problem but the present code using double precision arithmetic (16 decimal points) is insufficient for desired accuracy.

Figure 1. E 1 ( x ) as computed from Equation (18) for the range 0.25 x 5.00 .

Table 1. E 1 ( x ) values as given in Abramowitz and Stegun [1] and as computed by employing Equation (18) for a range of x values between 0.25 and 5.00.

Note however that for x 3.50 the nine decimal place accuracy is achieved and virtually maintained without need to increase n = 27 further. The FORTRAN code used for computations is included in the appendix.

5. Concluding Remarks

Exponential integral is evaluated by employing a novel series expansion, which has much better convergence characteristics compared to the corresponding standard expansion. While different means of evaluating the exponential integral are available, new contributions of the present approach may be attributed to its originality, uniform validity in the entire range of positive real numbers, and its potential of extension to various integrals of similar type with real or complex polynomials in the denominator as treated in [5] and [6].

Appendix

NOTES

1For the same case the fourth-order approximation of Laguerre in (2) has 7.47% error.

Cite this paper: Beji, S. (2021) Evaluation of Exponential Integral by Means of Fast-Converging Power Series. Advances in Pure Mathematics, 11, 101-108. doi: 10.4236/apm.2021.111006.
References

[1]   Abramowitz, M. and Stegun, I.A. (1972) Handbook of Mathematical Functions. Dover Publications, New York.

[2]   Laguerre, E.N. (1879) Sur l’intégrale x(e-x/x)dx. Bulletin de la Societé Mathématique de France, 7, 72-81.
https://doi.org/10.24033/bsmf.157

[3]   Laguerre, E.N. (1898) Œ uvres De Laguerre. Du Bureau des Longitudes, De L’École Polytechnique, Paris.

[4]   Milgram, M.S. (1985) The Generalized Integro-Exponential Function. Mathematics of Computation, 44, 443-458.
https://doi.org/10.1090/S0025-5718-1985-0777276-4

[5]   Beji, S. (2020) Resolution of Grandi’s Paradox and Investigations on Related Series. Applied Mathematics E-Notes, 20, 265-277.

[6]   Beji, S. (2020) Resolution of Grandi’s Paradox as Extended to Complex Valued Functions. Advances in Pure Mathematics, 10, 447-463.
https://doi.org/10.4236/apm.2020.108027

 
 
Top