Solution 11577: Algorithm Used in Computing Numeric Integration on the TI-82, TI-83 Family, TI-84 Family, TI-85, and TI-86 - INTERNAL ONLY
How do the TI-82, TI-83 family, TI-84 family, TI-85, and TI-86 perform numeric integration?
The method used in the TI-82, TI-83 family, TI-84 family, TI-85, and TI-86 is known as Gauss-Kronrod integration. It is a good deal more sophisticated than commonly known methods like Simpson's rule or even Romberg integration. In fact, the widely used mainframe library program QUADPACK uses this technique. A detailed discussion is much better left to the literature on this subject, which is voluminous. A few comments, however, may help the user in utilizing this function.
The concept begins with Gaussian quadrature rules, which have the property that by sampling a function at "n" points, an integral can be estimated that would be exact if the function were a polynomial of degree "2n-1". For comparison, the trapezoid rule is of polynomial degree one and Simpson's rule is of polynomial degree 3. But while Simpson's rule gives an exact integral for cubic polynomials based on three points, a Gauss quadrature rule with three points will give the exact integral of polynomials up to order 5. The next concept is to use a pair of Gauss quadrature rules of different order to provide both an integral result and a corresponding error estimate. The development of optimal extension pairs (pairs for which the lower order sample points are a subset of the higher order ones so that extra function evaluations are avoided) was done by A.S. Kronrod in 1965. The use of these pairs allow the method to be "adaptive". The integral is developed by beginning with the whole interval, computing the integral and the error estimate. If the error is too large, the interval is bisected and the process repeated on each half. Anytime a segment passes with regard to its "share" of the total error budget, that integral is added to the total integral. Other segments that don't pass are further divided. Note that this is quite different than some "adaptive" methods that iterate with smaller intervals until the result changes by less than some tolerance. The method computes error from two integral estimates and only takes new function samples in rapidly changing regions that require it.
This means that the "tolerance" needs to be specified to achieve for the integral. In the implementation this is an "absolute" error bound. This means the algorithm will quit when the absolute value of the difference in its two integrals is less than the tolerance you specify. Hence, it may not be wise to ask for a tolerance of .00001 if computing an integral that has a value of order 10,000 because 10 accurate digits must be given. Another thing to watch out for is the fact the function is sampled a finite number of times on the first integration attempt. These points are not equally spaced and are in fact clustered toward the end points. However, it is not uncommon to have a function that is virtually constant or perhaps zero throughout much of the range for integration. In this case, the integrator may quit early and fail to detect the existence of function behavior in other regions. An example of this is: fnInt(e^(-t/1e-6),t,0,1). This may seem appropriate to an electrical engineer who wants to compute all the "energy" in this waveform with a microsecond order fall time. However, just over a value of 0.002 in the interval 0 to 1, the value of this expression becomes less than 1e-999 and underflows to zero. So the algorithm returns a value of zero for this integral. A range of (0,.0001) would be more appropriate for the problem and also because the integral is small, a tol of 1e-10 is not excessive.
References:
"Numerical Methods and Software", David Kahaner, Cleve Moler and Stephen Nash, Prentice Hall, 1989, Chapter 5.
"QUADPACK, A Subroutine for Automatic Integration", R. Piessens, E. de Doncker-Kapenga, C.W. Uberhuber, D.K. Kahaner, Springer-Verlag, 1983.
Numerous articles covering this technique may also be found in journals such as:
Mathematics of Computation
ACM Transactions on Mathematical Software
SIAM Journal on Numerical Analysis