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 . Historically, Laguerre (1834-1886) gave a continued fraction expansion for its solution in , which was reprinted in . 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 . 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   investigated the integral and derived the following continued fraction representation as a solution.
which gives the following expressions according to the orders kept
Setting the lower limit of integration in the above expressions gives respectively,
which obviously converges to the nine decimal place correct value 0.219383934 given in the tables in [ , p. 240]. Laguerre also gave the second-order approximation result for as , which is accurate to five decimal places when compared to the nine decimal place correct value 0.003779352. For the present approach described in §3 attains the five decimal place accuracy with terms while the nine decimal place accuracy requires terms.
3. Integration by Use of a Novel Series Expansion
The exponential integral for real arguments is defined as [ , p. 262]
where is the incomplete gamma function and is the argument or arbitrary lower limit of integration. By introducing a change of variable in the form the integral in (4) transforms to
Evaluations of the above two-part integrals are carried out separately by employing relatively fast-converging series expansions developed in resolving Grandi’s paradox . Detailed derivation of the rectified and convergence-improved expansions of for and can be found in Beji . Accordingly, when or
For or Equation (6) is manipulated into the following convergent form :
where the coefficients for both (6) and (7) are defined as
We remark that unlike the standard expansion , 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 .
The series given in Equation (6) is used in the first part of the integration in (5) where .
Integrating by parts it can be shown that for an arbitrary value
Establishing separately and then adding all the terms give
Equation (11) can be expressed in a compact form as double summations.
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 .
Again integrating by parts successively we obtain for an arbitrary value
Starting from , which simply produces , we add the terms evaluated individually according to (14) for and obtain
Noting that the integral on the right is with integration variable u and expressing the summations in compact forms give
We now consider Equation (5) and make use of (12) and (16) so that
Evaluating the integral limits, rearranging, introducing a new running index , and gathering the terms multiplying together give
which in turn yields the series solution uniformly valid in the entire range of for the exponential integral .
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 versus x as computed from Equation (18) for the range . Table 1 compares the numerically computed and tabulated values taken from Abramowitz and Stegun [ , pp. 238-242] correct for nine decimal places. The relative error percentage is computed as . 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 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. as computed from Equation (18) for the range .
Table 1. values as given in Abramowitz and Stegun  and as computed by employing Equation (18) for a range of x values between 0.25 and 5.00.
Note however that for the nine decimal place accuracy is achieved and virtually maintained without need to increase 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  and .
1For the same case the fourth-order approximation of Laguerre in (2) has 7.47% error.