Romberg Integration
1. Information about
Romberg Integration
Romberg integration uses the Composite Trapezoidal rule to give preliminary approximations and than applies Richardson extrapolation process [1] to obtain improvements of the approximations.
To begin the explanation of Romberg
Integration recall that the Composite trapezoidal rule for approximating the
integral of a function f on an interval [a,b] using m subintervals is:
,
where
, 
and 
for each j=0,1,………,m.
We first obtain Composite Trapezoidal
rule approximation with 
,
……,and 
, where n is positive integer. The values of the step size 
corresponding to 
are 
, and with this notation the Trapezoidal rule becomes
where
is a number in (a,b).
At this point we can introduce the
notation 
as follows:
by
following this porcedure
than
the general expression is as follows:
.
The
following figure shows 
and 
.
Figure.1:and
of f(x)=(x2+x+1)cos(x) on the interval [0,pi/2].
If Richardson extraploation is applied to the above
formulas, (for detailed derivation see [1]) the following expression can be
obtained:
The following table shows the results of the above
expression.
The Romberg technique has the additional desirable
property that it allows an entire new row in the table to be calculated by
simply doing the one application of the Composite Trapezoidal rule and than
using the previously calculated values to obatin the succeeding entries in the
row. The method used to construct the a table of this type calculates the
entries row by row. The implementation of this technique can be seen in the
following algorithm:
INPUT endpoints a,b; integer n
OUTPUT an array R (
is the approximation to the value of the integral.)
Step 1: Set h = b-a;
.
Step 2: OUTPUT (
).
Step 3: For i =2,3,…….,n do Steps 4-8.
Step
4: Set 
Step
5: For j=2,3,……,i
Set
(Extrapolation)
Step
6: OUTPUT (
for j=1,2,……,i)
Step 7: Set h=h/2
Step 8: For j = 1,2,……,i set 
(Update row 1 of R)
Step 9: STOP.
The above algorithm requires the use of integer n
to determine the number of rows to be generated. It is often more useful to
define an error tolerance instead of this integer as a stopping criteria. This
tolerance is the difference of 
and 
. To guard against possibility that two consecutive row
elements agree with each other but not with the value of the integral being
approximated. Then the difference between 
and 
should also be checked. If both satisfy the tolerance, the
integral can be approximated by 
. This is not a universal safeguard but it provides
sufficiently accurate results.
2. Numerical Results and
Conclusions
f(x)=sin(x) on the interval [0,p]
The following figure is plot of this
function on the provided interval.
Figure.2: The plot of f(x)=sin(x)
on the provided interval.
The following table gives the results of the
integration of this function on the provided interval with different methods.
|
Trapezoidal |
Simpson’s |
Gaussian Quad. |
Romberg |
Epsilon |
|
1.9741510 |
2.0045765 |
1.9999906 |
1.9999999 |
0.01 |
|
1.9969667 |
2.0026127 |
1.9999999 |
1.9999999 |
0.001 |
|
1.9996646 |
2.0000007 |
1.9999999 |
1.9999999 |
0.0001 |
Table.1: The results for
f(x)=sin(x) on the interval [0, p]
f(x)=(100/x2)sin(10/x)
on the interval [1,3]
The following figure is the plot of this function on the provided
interval.
Figure.3: The plot of f(x)=(100/x2)sin(10/x)
on the interval [1,3]
The following table gives
the results of the integration of this function on the provided interval with
different methods.
|
Trapezoidal |
Simpson’s |
Gaussian Quad. |
Romberg |
Epsilon |
|
-1.3907865 |
-1.4294466 |
-1.4258159 |
-1.4260247 |
0.01 |
|
-1.4226971 |
-1.4258122 |
-1.4226024 |
-1.4260247 |
0.001 |
|
-1.4260349 |
-1.4226021 |
-1.4260241 |
-1.4260247 |
0.0001 |
Table.2: The results for
f(x)=(100/x2)sin(10/x) on the interval [1, 3]
f(x)=x2ex on the interval [0,1]
The following figure is plot of this
function on the provided interval.
Figure.4: The plot of f(x)=x2ex
on the interval [0,1]
The following table gives the results of the
integration of this function on the provided interval with different methods.
|
Trapezoidal |
Simpson’s |
Gaussian Quad. |
Romberg |
Epsilon |
|
0.1624884 |
0.1607224 |
0.1606027 |
0.1606280 |
0.01 |
|
0.1614894 |
0.1607224 |
0.1606027 |
0.1606280 |
0.001 |
|
0.1606369 |
0.1606103 |
0.1606027 |
0.1606279 |
0.0001 |
Table.3: The results for f(x)=x2ex
on the interval [0,1]
f(x)=cos(1/x) on the interval [0.1,2]
The following figure is plot of this
function on the provided interval.
Figure.5: The plot of
f(x)=cos(1/x) on the provided interval.
The following table gives the results of the
integration of this function on the provided interval with different methods.
|
Trapezoidal |
Simpson’s |
Gaussian Quad. |
Romberg |
Epsilon |
|
0.5814277 |
0.6739107 |
0.6528763 |
0.6738292 |
0.01 |
|
0.5814277 |
0.6737815 |
0.6738920 |
0.6738321 |
0.001 |
|
0.6736501 |
0.6738119 |
0.6738800 |
0.6738321 |
0.0001 |
Table.4: The results for
f(x)=cos(1/x) on the interval [0.1, 2]
f(x)=e2xsin(3x) on the interval [1,3]
The following figure is plot of this
function on the provided interval.
Figure.6: The plot of f(x)=e2xsin(3x)
on the interval [1,3]
The following table gives the results of the
integration of this function on the provided interval with different methods.
|
Trapezoidal |
Simpson’s |
Gaussian Quad. |
Romberg |
Epsilon |
|
108.5886632 |
108.5560839 |
108.5552994 |
108.5552994 |
0.01 |
|
108.5555266 |
108.5553434 |
108.5552994 |
108.5552994 |
0.001 |
|
108.5553235 |
108.5552971 |
108.5552812 |
108.5552812 |
0.0001 |
Table.5: The results for f(x)=e2xsin(3x)
on the interval [1,3]
By
looking at the above results, it can be said that all of the four methods are
good at calculating the integrals of oscillating functions. But the effect of
epsilon on different methods is different. The Trapezoidal and the Simpson’s
Rules are not giving sufficiently accurate results for epsilon bigger than
0.001. But the methods of Romberg Integration and Quassian Quadratures are
giving nearly same results and these results are better than the ones obtained
from the other two methods especially for bigger values of epsilon. The source
of this difference is the composite rules. The approximations done by
Trapezoidal and Simpson’s Methods for one subinterval are not as accurate as
the ones that are done by the other two methods. This also causes more levels
of halving for the Trapezoidal and Simpson’s Rules. Then although we didn’t
provide the exact numbers here, the time that is required for the calculation
of the integrals of the functions given above can be ordered as follows:
3. References
[1] Numerical Analysis, R.L. Burden, J.D.
Faires, PWS Publishing Company, Boston 1993.