The accurate evaluation of the repeated integrals of the coerror function defined as
with is required in a number of important application areas in the physical sciences. The usual method for evaluating this function is based on three-term recurrence relationships, and controlling loss of accuracy is tricky.The method proposed here is based on numerical quadrature and gives a guaranteed accuracy over a predetermined domain. Assuming the use of IEEE standard, double precision arithmetic, a domain, D(n,x), is first determined in which the computation of the function will neither underflow nor overflow and 12 significant decimal digit accuracy may be obtained using 200 or less quadrature points. The domain is found using a combination of mathematical and computational techniques.
Further analysis using a MATLAB software implementation allows estimations of the number of quadrature points, N*, required for 12-digit accuracy to be made over subintervals of the positive x-range. These estimations may then be used to produce an expression for computing a suitable value for N* for a given x ≥ 0, thus providing a computationally efficient evaluation of the function. Test results are presented that confirm the predicted accuracy.
MATLAB functions are provided that implement these methods both in IEEE double precision and in extended precision. The software used to perform the accuracy and domain computations is also available.
This paper is a well-written account of the detail necessary to generate an accurate and efficient numerical method for evaluating the special function defined above. It will be of interest both to practitioners who require the accurate results that the associated software generates and to specialists in the numerical evaluation of special functions where it will act as a good example of a careful analysis of a new numerical approach.