Demostración de Integración de Romberg

Demostración de Integración de Romberg 

https://personal.math.ubc.ca/~CLP/CLP2/clp_2_ic/ap_Romberg.html

Romberg Integration

🔗

🔗

The formulae (E4a,b) for $𝐾$ and $\mathcal{𝐴}$ are, of course, only¹ approximate since they are based on (E2), which is an approximation to (E1). Let's repeat the derivation that leads to (E4), but using the full (E1),

“Only” is a bit strong. Don't underestimate the power of a good approximation (pun intended).

$$\mathcal{𝐴}=𝐴(ℎ)+𝐾{ℎ}^{𝑘}+{𝐾}_1{ℎ}^{𝑘+1}+{𝐾}_2{ℎ}^{𝑘+2}+\cdots$$

🔗

Once again, suppose that we have chosen some $ℎ$ and that we have evaluated $𝐴(ℎ)$ and $𝐴(ℎ/2)\text{.}$ They obey

$$\begin{array}{}\text{(E5a)}& \mathcal{𝐴}& =𝐴(ℎ)+𝐾{ℎ}^{𝑘}+{𝐾}_1{ℎ}^{𝑘+1}+{𝐾}_2{ℎ}^{𝑘+2}+\cdots \text{(E5b)}& \mathcal{𝐴}& =𝐴(ℎ/2)+𝐾(\frac{ℎ}{2}{)}^{𝑘}+{𝐾}_1(\frac{ℎ}{2}{)}^{𝑘+1}+{𝐾}_2(\frac{ℎ}{2}{)}^{𝑘+2}+\cdots \end{array}$$

🔗

Now, as we did in the derivation of (E4b), multiply (E5b) by $2^{𝑘}$ and then subtract (E5a). This gives

$$(2^{𝑘}-1)\mathcal{𝐴}=2^{𝑘}𝐴(ℎ/2)-𝐴(ℎ)-\frac{1}{2}{𝐾}_1{ℎ}^{𝑘+1}-\frac{3}{4}{𝐾}_2{ℎ}^{𝑘+1}+\cdots$$

🔗

and then, dividing across by $(2^{𝑘}-1)\text{,}$

$$\mathcal{𝐴}=\frac{2^{𝑘}𝐴(ℎ/2)-𝐴(ℎ)}{2^{𝑘}-1}-\frac{1/2}{2^{𝑘}-1}{𝐾}_1{ℎ}^{𝑘+1}-\frac{3/4}{2^{𝑘}-1}{𝐾}_2{ℎ}^{𝑘+2}+\cdots$$

🔗

Hence if we define our “new improved approximation”

$$\begin{array}{}\text{(E6)}& 𝐵(ℎ)=\frac{2^{𝑘}𝐴(ℎ/2)-𝐴(ℎ)}{2^{𝑘}-1}\text{and}\widetilde{𝐾}=-\frac{1/2}{2^{𝑘}-1}{𝐾}_1\text{and}{\widetilde{𝐾}}_1=-\frac{3/4}{2^{𝑘}-1}{𝐾}_2\end{array}$$

🔗

we have

$$\mathcal{𝐴}=𝐵(ℎ)+\widetilde{𝐾}{ℎ}^{𝑘+1}+{\widetilde{𝐾}}_1{ℎ}^{𝑘+2}+\cdots$$

🔗

which says that $𝐵(ℎ)$ is an approximation to $\mathcal{𝐴}$ whose error is of order $𝑘+1\text{,}$ one better² than $𝐴(ℎ)$'s.

🔗

🔗

If $𝐴(ℎ)$ has been computed for three values of $ℎ\text{,}$ we can generate $𝐵(ℎ)$ for two values of $ℎ$ and repeat the above procedure with a new value of $𝑘\text{.}$ And so on. One widely used numerical integration algorithm, called Romberg integration³, applies this procedure repeatedly to the trapezoidal rule. It is known that the trapezoidal rule approximation $𝑇(ℎ)$ to an integral $𝐼$ has error behaviour (assuming that the integrand $𝑓(𝑥)$ is smooth)

$$𝐼=𝑇(ℎ)+{𝐾}_1{ℎ}^2+{𝐾}_2{ℎ}^4+{𝐾}_3{ℎ}^6+\cdots$$

🔗

Only even powers of $ℎ$ appear. Hence

$$\begin{array}{rlrl}𝑇(ℎ)& & & \text{ has error of order 2}\end{array}$$

so that, using (E6) with $𝑘=2\text{,}$

$$\begin{array}{rlrl}{𝑇}_1(ℎ)& =\frac{4𝑇(ℎ/2)-𝑇(ℎ)}{3}& & \text{ has error of order 4}\end{array}$$

so that, using (E6) with $𝑘=4\text{,}$

$$\begin{array}{rlrl}{𝑇}_2(ℎ)& =\frac{16{𝑇}_1(ℎ/2)-{𝑇}_1(ℎ)}{15}& & \text{ has error of order 6}\end{array}$$

so that, using (E6) with $𝑘=6\text{,}$

$$\begin{array}{rlrl}{𝑇}_3(ℎ)& =\frac{64{𝑇}_2(ℎ/2)-{𝑇}_2(ℎ)}{63}& & \text{ has error of order 8}\end{array}$$

🔗

and so on. We know another method which produces an error of order 4 — Simpson's rule. In fact, ${𝑇}_1(ℎ)$ is exactly Simpson's rule (for step size $\frac{ℎ}{2}$).

🔗

Equation C.2.1. Romberg integration. 

 

Let $𝑇(ℎ)$ be the trapezoidal rule approximation, with step size $ℎ\text{,}$ to an integral $𝐼\text{.}$ The Romberg integration algorithm is

$$\begin{array}{rl}{𝑇}_1(ℎ)& =\frac{4𝑇(ℎ/2)-𝑇(ℎ)}{3}\\ {𝑇}_2(ℎ)& =\frac{16{𝑇}_1(ℎ/2)-{𝑇}_1(ℎ)}{15}\\ {𝑇}_3(ℎ)& =\frac{64{𝑇}_2(ℎ/2)-{𝑇}_2(ℎ)}{63}\\ & ⋮\\ {𝑇}_{𝑘}(ℎ)& =\frac{2^{2𝑘}{𝑇}_{𝑘-1}(ℎ/2)-{𝑇}_{𝑘-1}(ℎ)}{2^{2𝑘}-1}\\ & ⋮\end{array}$$

🔗

Example C.2.2. Finding $𝜋$ by Romberg integration. 

The following table⁴ illustrates Romberg integration by applying it to the area $𝐴$ of the integral $𝐴={\int }_0^1\frac{4}{1+{𝑥}^2}\mathrm{d}𝑥\text{.}$ The exact value of this integral is $𝜋$ which is $3.14159265358979\text{,}$ to fourteen decimal places.

$ℎ$	$𝑇(ℎ)$	${𝑇}_1(ℎ)$	${𝑇}_2(ℎ)$	${𝑇}_3(ℎ)$
1/4	3.13118	3.14159250246	3.14159266114	3.14159265359003
1/8	3.13899	3.141592651225	3.141592653708
1/16	3.14094	3.141592653553
1/32	3.14143

This computation required the evaluation of $𝑓(𝑥)=\frac{4}{1+{𝑥}^2}$ only for $𝑥=\frac{𝑛}{32}$ with $0\le 𝑛\le 32$ — that is, a total of $33$ evaluations of $𝑓\text{.}$ Those $33$ evaluations gave us $12$ correct decimal places. By way of comparison, $𝑇(\frac{1}{32})$ used the same $33$ evaluations of $𝑓\text{,}$ but only gave us $3$ correct decimal places.

🔗

As we have seen, Richardson extrapolation can be used to choose the step size so as to achieve some desired degree of accuracy. We are next going to consider a family of algorithms that extend this idea to use small step sizes in the part of the domain of integration where it is hard to get good accuracy and large step sizes in the part of the domain of integration where it is easy to get good accuracy. We will illustrate the ideas by applying them to the integral ${\int }_0^1\sqrt{𝑥}\mathrm{d}𝑥\text{.}$ The integrand $\sqrt{𝑥}$ changes very quickly when $𝑥$ is small and changes slowly when $𝑥$ is large. So we will make the step size small near $𝑥=0$ and make the step size large near $𝑥=1\text{.}$