Just to recall what this is all about, the following formula is known as Simpson's Rule. It approximates the integral of a function #f:[a,b] to mathbb{R}#, where #[a,b]# is a closed bound interval.
#int_a^b f(x) dx approx {b-a}/6 [f(a)+4 f({a+b}/2)+f(b)]#.
The error committed in this approximation process is given by #1/90 ({b−a}/2)^5 f^{(4)}(zeta)# for some #zeta in [a,b]#. Proving that this is really the expression that describes the error committed by Simpson's approximation requires a lot of space and theorems, so I will skip it. Interested readers could find it on page 6 of the paper Simpson’s Rule Is Exact for Quintics by Louis A. Talman.
So it has to exist some #zeta in [a,b]# such that #E=1/90 ({b−a}/2)^5 f^{(4)}(zeta)#, and this is the exact value of the error. Notice that this formula works only in the hypothesis that #f# admits the fourth derivative for every point #zeta in (a,b)#.
When dealing with approximation of integrals, one is not interested in finding exact values. So, we don't either require to know the error value exactly. In fact, if we compute the exact error, then we can retrieve the exact value of the integral by simply adding (subtracting) the error to (from) the approximated value. In a few words, searching for the exact error is the same thing as searching for the exact value.
The interesting fact about error (actually the reason why approximation is so widely used and useful, especially in applications) is that it can be estimated, often considering the worst case. In this manner we don't have exact information, but we are sure that the error on the approximated value is bounded. Once this bound is known, we could decide if the precision accomplished by the approximation is enough or not, even in the worst case (and maybe search for a way to improve the approximation...).
In our case, we can think of #E# as a function of #zeta# (I adjust it with the absolute value because error is expressed as a positive number):
#E(zeta)=1/90 ({b−a}/2)^5 |f^{(4)}(zeta)|#
The structure of this function depends on #f#. If #f# is differentiable five times in #(a,b)#, then the maximum error can be found using the First Derivative Test on #E(zeta)# and comparing the values. Otherwise, other techniques has to be employed.
Example 1
We want to estimate #int_0^2 e^{x} dx#. Let's try to approximate the integral value via Simpson's Rule:
#int_0^2 e^{x} dx approx (2-0)/6 [e^0+4e^1+e^2]=1/3(1+4e+e^2) approx 6.421#
Now we want to express a measure of "how much wrong" can this result be: an upper bound for the error. We know that the fourth derivative of #e^x# is #e^x#. By the formula given above, we get #E(zeta)=1/90 ({2−0}/2)^5 |e^{zeta}|=e^{zeta}/90#
So #E'(zeta)=e^{zeta}/90>0# for all #zeta in [0,2]#. The function is strictly increasing on #[0,2]#, so #E(2)=e^2/90 < 0.083# is the maximum value of #E(zeta)# on #[0,2]#. So #E<0.083#.
We conclude that
#6.421-0.083 < int_0^2 e^{x} dx < 6.421+0.083#
Relative error is #< 1,3 %#.
Thanks to the Fundamental Theorem of Calculus we know that the exact value for this integral is #e^2-1 approx 6.389#, but to get this value we need to be able to compute antiderivatives (with #e^x# the computation is very simple, but in other cases it can be very hard or even impossible).
Now, we can compute the exact error because we know the exact value of the integral:
#|[1/3(1+4e+e^2)] - [e^2-1]| approx 0.032#
The upper bound on #E# that we found previously is consistent with this result, which is a smaller number. In this case, Simpson's formula is more than twice preciser than we could estimate. And the error committed is quite small!
Example 2
Now we want to approximate the value of #int_{-3}^3 e^(-x^2) dx#.
graph{e^(-x^2) [-3.5, 3.5, -0.5, 1.5]}
We are not able to express the antiderivative of the function #e^(-x^2)# in terms of finite compositions of functions. So there's no way to evaluate this integral precisely. We have to approximate!
By Simpson's Rule:
#int_{-3}^3 e^(-x^2) dx approx {3+3}/6 [e^{-9}+4e^0+e^{-9}]=2e^{-9}+4 approx 4.000#
For some #zeta in [-3,3]# the error is given by
#E(zeta)=1/90((3+3)/2)^5 |4 e^(-zeta^2) (4 zeta^4-12 zeta^2+3)|=27/10 |4 e^(-zeta^2) (4 zeta^4-12 zeta^2+3)|#
graph{|54/5e^(-x^2)(4x^4-12x^2+3)| [-4, 4, -2, 35]}
The global maximum on #[-3,3]# of #E(zeta)# is #E(0)=27/10 |4e^0 *3|=162/5=32.4#, so the value of t
This error is very high and the result is the following
#-28.400 <= int_{-3}^3 e^{-x^2} dx <= 36.400#
which is not very accurate! The relative error can be bounded only by #902 %#... not really a satisfying result.
To improve the precision we can split the integral as follows:
#int_{-3}^3 e^{-x^2} dx=int_{-3}^{-2} e^{-x^2} dx+int_{-2}^{-1} e^{-x^2} dx+int_{-1}^0 e^{-x^2} dx+int_{0}^1 e^{-x^2} dx+int_{1}^2 e^{-x^2} dx+int_{2}^3 e^{-x^2} dx#
and use Simpson's Rule on each of them. By symmetry, this reduces to
#int_{-3}^3 e^{-x^2} dx=2(int_{0}^1 e^{-x^2} dx+int_{1}^2 e^{-x^2} dx+int_{2}^3 e^{-x^2} dx) approx 2(0.747+0.135+0.004)=1.772#
Error on each interval is given by #tilde{E}(zeta)=1/2880 |4 e^(-zeta^2) (4 zeta^4-12 zeta^2+3)|# and we have to maximize it. The global maximum on #[0,1]# is at #x=0#, on #[1,2]# at #x=1# and on #[2,3]# at #x=sqrt{5/2+sqrt{5/2}}#. So the error bound is given by
#2*(tilde{E}(0)+tilde{E}(1)+tilde{E}(sqrt{5/2+sqrt{5/2}})) <2*(0.005+0.003+0.001)=0.018#
The relative error is #< 1,1%#.
The result is that we get a much preciser value of the integral:
#1.754 < int_{-3}^3 e^{-x^2} dx < 1.790#
