The Simpson method uses quadratic approximations of the
function on each subinterval
(three function values, at ,
and are used to determine uniquely the quadratic polynomial).
We add up all the integrals of the quadratic polynomials on all the
subintervals, to obtain an approximation of the function integral
over . This approximation is
(the Simpson formula)
An efficient implementation of the algorithm that computes
for different values of
might use 3 variables:
To implement the adaptive Simpson method, we need to compute efficiently, re-using the function values that were calculated for . The key is to observe that, when doubling the number of subintervals, the old node points are a subset of the new node points become are a subset of the node points .