Subject : The following demonstrates the Romberg method of estimating integrals of continuous functions.
Authors : Nathan Collier, Autar Kaw, Loubna Guennoun
Date : 26 October 2005

Introduction

Romberg integration is based on the trapezoidal rule, where we use two estimates of an integral to compute a value that is more accurate than the previous estimates. This is called Richardson's extrapolation. Thus,
	I = I(h) + E(h)
	h = (b-a) / n
where I is the exact value of the integral, I(h) is the approximate integral using the trapezoidal rule with n segments, and E(h) is the truncation error. A general form of Romberg integration is:
	\(TraditionalForm\`I\_\(j, k\)\  = \(\(\ \)\(4\^\(k - 1\)\ I\_\(j + 1, k - 1\)\  - \ I\_\(j, k - 1\)\)\)\/\(\(\ \)\(4\^\(k - 1\)\  - \ 1\)\)\)
where the index j is the order of the estimate integral, and k is the level of integration. [click here for textbook notes] [click here for power point presentation].

Inputs

The following simulation will illustrate Romberg integration. This section is the only section where the user may interacts with the program. The user may enter any function f(x) and the lower and upper limit for the function. By entering this data, the program will calculate the exact value of the integral, followed by the results using the trapezoidal rule with n = 1, 2, 4, 8 segments, and the Romberg integration for each segments.

Integrand f[x] = 0

f[x_] := 300.0*x/(1.0 + Exp[x]);

Lower and upper limit of the integral, a and b respectively.

a = 0.0;
b = 10.0;

curve = 
        Plot[f[x], {x, a, b}, 
          PlotLabel -> "Entered function on given interval", 
          TextStyle -> {FontSize -> 
                11}];
General::spelll  Off

1 Segment

n = 1;
h1 = (b - a)/n

The integral of the function f(x) from a to b using the trapezoidal rule with one segment will be equal to:

(Subscript[I, 1, 1]) = h1*(f[a] + f[b])/2

NOTE: In the index 1,1, the first number "1" means we are integrating with n=1 segment, and the second number "1" is the first iteration, using the original trapezoidal rule, which corresponds to O(h2).

2 Segment

n = 2;
h2 = (b - a)/n

The integral of the function f(x) from a to b using the trapezoidal rule with two segments will be equal to:

(Subscript[I, 2, 1]) = h2*(f[a] + 2*f[a + h2] + f[b])/2

NOTE: In the index of I the number "2" means we are integrating with n=2 segments, and the second number "1" is the first iteration, using the original trapezoidal rule, which corresponds to O(h2).

j = 1;
k = 2;
(Subscript[I, 1, 2]) = (4^(k - 1)*(Subscript[I, j + 1, k - 1]) - 
            (Subscript[I, j, k - 1]))/(4^(k - 1) - 1)

NOTE: In the index of I the number "1" corresponds to the first result of the second iteration, and the second number "2" (2nd iteration) corresponds to O(h4).

The approximate error is:

Ea2 = (Subscript[I, 1, 2]) - (Subscript[I, 1, 1])

The absolute relative approximate error is

ea2 = Abs[Ea2/(Subscript[I, 1, 2])]*100

4 Segment

n = 4;
h4 = (b - a)/n

The integral of the function f(x) from a to b using the trapezoidal rule with four segments will be equal to:

summ = 
        f[a + h4] + f[a + 2*h4] + f[a + 3*h4];
(Subscript[I, 3, 1]) = h4*(f[a] + 2*summ + f[b])/2

NOTE: In the index of I the first number "3" corresponds to n=4 segments, and the second number "1" is the first iteration using the original trapezoidal rule, which corresponds to O(h2).

j = 2;
k = 2;
(Subscript[I, 2, 2]) = (4^(k - 1)*(Subscript[I, j + 1, k - 1]) - 
            (Subscript[I, j, k - 1]))/(4^(k - 1) - 1)

NOTE: In the index of I the number "2" corresponds to the second result of the second iteration, and the second number "2" (2nd iteration) corresponds to O(h4).

j = 1;
k = 3;
(Subscript[I, 1, 3]) = (4^(k - 1)*(Subscript[I, j + 1, k - 1]) - 
            (Subscript[I, j, k - 1]))/(4^(k - 1) - 1)

NOTE: In the index of I the number "1" corresponds to the first result of the second iteration, and the second number "3" (3rd iteration) corresponds to O(h6).

The approximate error is:

Ea4 = (Subscript[I, 1, 3]) - (Subscript[I, 1, 2])

The absolute relative approximate error is

ea4 = Abs[Ea4/(Subscript[I, 1, 3])]*100

8 Segment

n = 8;
h8 = (b - a)/n

The integral of the function f(x) from a to b using the trapezoidal rule with eight segments will be equal to:

summ = 
        f[a + h8] + f[a + 2*h8] + f[a + 3*h8] + f[a + 4*h8] + f[a + 5*h8] + 
          f[a + 6*h8] + f[a + 7*h8];
(Subscript[I, 4, 1]) = h8*(f[a] + 2*summ + f[b])/2

NOTE: In the index of I the first number "4" corresponds to n=8 segments, and the second number "1" is the first iteration using the original trapezoidal rule, which corresponds to O(h2).

j = 3;
k = 2;
(Subscript[I, 3, 2]) = (4^(k - 1)*(Subscript[I, j + 1, k - 1]) - 
            (Subscript[I, j, k - 1]))/(4^(k - 1) - 1)

NOTE: In the index of I the number "3" corresponds to the third result of the second iteration, and the second number "2" (2nd iteration) corresponds to O(h4).

j = 2;
k = 3;
(Subscript[I, 2, 3]) = (4^(k - 1)*(Subscript[I, j + 1, k - 1]) - 
            (Subscript[I, j, k - 1]))/(4^(k - 1) - 1)

NOTE: In the index of I the number "3" corresponds to the third result of the second iteration, and the second number "2" (2nd iteration) corresponds to O(h4).

j = 1;
k = 4;
(Subscript[I, 1, 4]) = (4^(k - 1)*(Subscript[I, j + 1, k - 1]) - 
            (Subscript[I, j, k - 1]))/(4^(k - 1) - 1)

NOTE: In the index of I the number "1" corresponds to the first result of the second iteration, and the second number "4" (4th iteration) corresponds to O(h8).

The approximate error is:

Ea8 = (Subscript[I, 1, 4]) - (Subscript[I, 1, 3])

The absolute relative approximate error is

ea8 = Abs[Ea8/(Subscript[I, 1, 4])]*100