C.2: Integración Romberg
- Page ID
- 119266
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\( \newcommand{\id}{\mathrm{id}}\) \( \newcommand{\Span}{\mathrm{span}}\)
( \newcommand{\kernel}{\mathrm{null}\,}\) \( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\) \( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\) \( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\id}{\mathrm{id}}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\kernel}{\mathrm{null}\,}\)
\( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\)
\( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\)
\( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\) \( \newcommand{\AA}{\unicode[.8,0]{x212B}}\)
\( \newcommand{\vectorA}[1]{\vec{#1}} % arrow\)
\( \newcommand{\vectorAt}[1]{\vec{\text{#1}}} % arrow\)
\( \newcommand{\vectorB}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vectorC}[1]{\textbf{#1}} \)
\( \newcommand{\vectorD}[1]{\overrightarrow{#1}} \)
\( \newcommand{\vectorDt}[1]{\overrightarrow{\text{#1}}} \)
\( \newcommand{\vectE}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash{\mathbf {#1}}}} \)
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\(\newcommand{\avec}{\mathbf a}\) \(\newcommand{\bvec}{\mathbf b}\) \(\newcommand{\cvec}{\mathbf c}\) \(\newcommand{\dvec}{\mathbf d}\) \(\newcommand{\dtil}{\widetilde{\mathbf d}}\) \(\newcommand{\evec}{\mathbf e}\) \(\newcommand{\fvec}{\mathbf f}\) \(\newcommand{\nvec}{\mathbf n}\) \(\newcommand{\pvec}{\mathbf p}\) \(\newcommand{\qvec}{\mathbf q}\) \(\newcommand{\svec}{\mathbf s}\) \(\newcommand{\tvec}{\mathbf t}\) \(\newcommand{\uvec}{\mathbf u}\) \(\newcommand{\vvec}{\mathbf v}\) \(\newcommand{\wvec}{\mathbf w}\) \(\newcommand{\xvec}{\mathbf x}\) \(\newcommand{\yvec}{\mathbf y}\) \(\newcommand{\zvec}{\mathbf z}\) \(\newcommand{\rvec}{\mathbf r}\) \(\newcommand{\mvec}{\mathbf m}\) \(\newcommand{\zerovec}{\mathbf 0}\) \(\newcommand{\onevec}{\mathbf 1}\) \(\newcommand{\real}{\mathbb R}\) \(\newcommand{\twovec}[2]{\left[\begin{array}{r}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\ctwovec}[2]{\left[\begin{array}{c}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\threevec}[3]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\cthreevec}[3]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\fourvec}[4]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\cfourvec}[4]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\fivevec}[5]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\cfivevec}[5]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\mattwo}[4]{\left[\begin{array}{rr}#1 \amp #2 \\ #3 \amp #4 \\ \end{array}\right]}\) \(\newcommand{\laspan}[1]{\text{Span}\{#1\}}\) \(\newcommand{\bcal}{\cal B}\) \(\newcommand{\ccal}{\cal C}\) \(\newcommand{\scal}{\cal S}\) \(\newcommand{\wcal}{\cal W}\) \(\newcommand{\ecal}{\cal E}\) \(\newcommand{\coords}[2]{\left\{#1\right\}_{#2}}\) \(\newcommand{\gray}[1]{\color{gray}{#1}}\) \(\newcommand{\lgray}[1]{\color{lightgray}{#1}}\) \(\newcommand{\rank}{\operatorname{rank}}\) \(\newcommand{\row}{\text{Row}}\) \(\newcommand{\col}{\text{Col}}\) \(\renewcommand{\row}{\text{Row}}\) \(\newcommand{\nul}{\text{Nul}}\) \(\newcommand{\var}{\text{Var}}\) \(\newcommand{\corr}{\text{corr}}\) \(\newcommand{\len}[1]{\left|#1\right|}\) \(\newcommand{\bbar}{\overline{\bvec}}\) \(\newcommand{\bhat}{\widehat{\bvec}}\) \(\newcommand{\bperp}{\bvec^\perp}\) \(\newcommand{\xhat}{\widehat{\xvec}}\) \(\newcommand{\vhat}{\widehat{\vvec}}\) \(\newcommand{\uhat}{\widehat{\uvec}}\) \(\newcommand{\what}{\widehat{\wvec}}\) \(\newcommand{\Sighat}{\widehat{\Sigma}}\) \(\newcommand{\lt}{<}\) \(\newcommand{\gt}{>}\) \(\newcommand{\amp}{&}\) \(\definecolor{fillinmathshade}{gray}{0.9}\)Las fórmulas (E4a, b) para\(K\) y\(\cA\) son, por supuesto, solo 1 “Only” es un poco fuerte. No subestimes el poder de una buena aproximación (juego de palabras intencionado). aproximadas ya que se basan en (E2), que es una aproximación a (E1). Repetimos la derivación que lleva a (E4), pero usando el completo (E1),
\[ \cA=A(h)+Kh^k+K_1h^{k+1}+K_2h^{k+2}+\cdots \nonumber \]
Una vez más, supongamos que hemos elegido algunos\(h\) y que hemos evaluado\(A(h)\) y\(A(h/2)\text{.}\) Ellos obedecen
\ begin {align*}\ ca&=A (h) +kh^k+k_1h^ {k+1} +k_2h^ {k+2} +\ cdots\ tag {e5a}\\\ ca&=A (h/2) +K\ grande (\ tfrac {h} {2}\ grande) ^k +K_1\ grande (\ tfrac {h} 2}\ grande) ^ {k+1} +K_2\ grande (\ tfrac {h} {2}\ grande) ^ {k+2} +\ cdots\ tag {e5b}\ end {align*}
Ahora, como hicimos en la derivación de (E4b), multiplicar (E5b) por\(2^k\) y luego restar (E5a). Esto da
\[ \left(2^k-1\right)\cA=2^kA(h/2)-A(h)-\half K_1 h^{k+1}-\tfrac{3}{4} K_2 h^{k+1}+\cdots \nonumber \]
y luego, dividiendo entre\(\left(2^k-1\right)\text{,}\)
\[ \cA=\frac{2^kA(h/2)-A(h)}{2^k-1}-\frac{1/2}{2^k-1} K_1 h^{k+1}-\frac{3/4}{2^k-1} K_2 h^{k+2}+\cdots \nonumber \]
De ahí que si definimos nuestra “nueva aproximación mejorada”
\[ B(h)=\frac{2^kA(h/2)-A(h)}{2^k-1}\ \text{and}\ \tilde K = -\frac{1/2}{2^k-1}K_1\ \text{and}\ \tilde K_1 = -\frac{3/4}{2^k-1}K_2 \tag{E6} \nonumber \]
tenemos
\[ \cA=B(h)+\tilde K h^{k+1}+\tilde K_1 h^{k+2}+\cdots \nonumber \]
que dice que\(B(h)\) es una aproximación a\(\cA\) cuyo error es de orden\(k+1\text{,}\) uno mejor 2 Es decir, el error decae como\(h^{k+1}\) as opposed to \(h^k\) — so, as \(h\) decreases, it gets smaller faster. que\(A(h)\)'s.
Si se\(A(h)\) ha calculado para tres valores de\(h\text{,}\) podemos generar\(B(h)\) para dos valores de\(h\) y repetir el procedimiento anterior con un nuevo valor de\(k\text{.}\) Y así sucesivamente. Un algoritmo de integración numérica ampliamente utilizado, llamado Romberg integration 3 Romberg Integration fue introducido por el alemán Werner Romberg (1909—2003) en 1955. , aplica este procedimiento reiteradamente a la regla trapezoidal. Se sabe que la aproximación de la regla trapezoidal\(T(h)\) a una integral\(I\) tiene un comportamiento de error (asumiendo que el integrando\(f(x)\) es suave)
\[ I=T(h)+K_1h^2+K_2h^4+K_3h^6+\cdots \nonumber \]
Sólo incluso poderes de\(h\) aparecer. De ahí
\[\begin{align*} T(h)& &&\text{ has error of order 2}\\ \end{align*}\]
de modo que, usando (E6) con\(k=2\text{,}\)
\ begin {align*} T_1 (h) &=\ tfrac {4T (h/2) -T (h)} {3} &&\ text {tiene error de orden 4}\\\ end {align*}
de modo que, usando (E6) con\(k=4\text{,}\)
\ begin {align*} T_2 (h) &=\ tfrac {16T_1 (h/2) -T_1 (h)} {15} &&\ text {tiene error de orden 6}\\\ end {alinear*}
de modo que, usando (E6) con\(k=6\text{,}\)
\ begin {align*} T_3 (h) &=\ tfrac {64T_2 (h/2) -T_2 (h)} {63} &&\ text {tiene error de orden 8}\ end {align*}
y así sucesivamente. Conocemos otro método que produce un error de orden 4: la regla de Simpson. De hecho,\(T_1(h)\) es exactamente la regla de Simpson (para el tamaño del paso\(\frac{h}{2}\)).
\(T(h)\)Sea la aproximación de regla trapezoidal, con tamaño de paso\(h\text{,}\) a una integral\(I\text{.}\) El algoritmo de integración Romberg es
\ begin {align*} T_1 (h) &=\ tfrac {4T (h/2) -T (h)} {3}\\ T_2 (h) &=\ tfrac {16T_1 (h/2) -T_1 (h)} {15}\\ T_2 (h) &=\ tfrac {64T_2 (h/2) -T_2 (h)} {63}\\ &\\ vdots\\ t_k (h) &=\ tfrac {2^ {2k} T_ {k-1} (h/2) -T_ {k-1} (h)} {2^ {2k} -1}\\ &\\ vdots\ end {align*}
El siguiente cuadro 4 La segunda columna, por ejemplo, de la tabla solo informes\(5\) decimal places for \(T(h)\text{.}\) But many more decimal places of \(T(h)\) were used in the computations of \(T_1(h)\) etc. ilustra la integración Romberg aplicándola al área\(A\) de la integral\(A=\int_0^1\frac{4}{1+x^2}\, d{x}\text{.}\) El valor exacto de esta integral es\(\pi\) que es\(3.14159265358979\text{,}\) a catorce decimales.
\(h\) | \(T(h)\) | \(T_1(h)\) | \(T_2(h)\) | \(T_3(h)\) |
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 |
Este cómputo requirió la evaluación de\(f(x)=\frac{4}{1+x^2}\) solo para\(x=\frac{n}{32}\) con\(0\le n\le 32\) — es decir, un total de\(33\) evaluaciones de\(f\text{.}\) Esas\(33\) evaluaciones nos dieron cifras decimales\(12\) correctas. A modo de comparación,\(T\big(\frac{1}{32}\big)\) se utilizaron las mismas\(33\) evaluaciones de\(f\text{,}\) pero sólo nos dieron cifras decimales\(3\) correctas.
Como hemos visto, la extrapolación de Richardson se puede utilizar para elegir el tamaño de paso para lograr algún grado de precisión deseado. A continuación vamos a considerar una familia de algoritmos que amplían esta idea para utilizar pequeños tamaños de paso en la parte del dominio de la integración donde es difícil obtener una buena precisión y grandes tamaños de paso en la parte del dominio de integración donde es fácil obtener una buena precisión. Ilustraremos las ideas aplicándolas a la integral\(\int_0^1 \sqrt{x}\ \, d{x}\text{.}\) El integrando\(\sqrt{x}\) cambia muy rápidamente cuando\(x\) es pequeño y cambia lentamente cuando\(x\) es grande. Entonces haremos que el tamaño del paso sea pequeño cerca\(x=0\) y haremos que el tamaño del paso sea grande cerca\(x=1\text{.}\)