Saltar al contenido principal
LibreTexts Español

D.3: Métodos de Tamaño de Paso Variable

  • Page ID
    119319
  • \( \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}}\)

    Ahora introducimos una familia de procedimientos que deciden por sí mismos qué tamaño de paso usar. En todos estos procedimientos el usuario especifica una tasa de error aceptable y el procedimiento intenta ajustar el tamaño del paso para que cada paso introduzca error a no más de esa tasa. De esa manera, el procedimiento utiliza un tamaño de paso pequeño cuando es difícil obtener una aproximación precisa, y un tamaño de paso grande cuando es fácil obtener una buena aproximación.

    Supongamos que deseamos generar una aproximación al problema de valor inicial

    \[ y'=f(t,y),\qquad y(t_0)=y_0 \nonumber \]

    para algún rango de\(t\)'s y queremos que el error introducido por unidad aumente 1 Sabemos que el error va\(t\) a hacerse más grande cuanto más nos adentremos\(t\text{.}\) So it makes sense to try to limit the error per unit increase in \(t\text{.}\) de ser no más que sobre algún pequeño número fijo\(\varepsilon\text{.}\) Esto quiere decir que si\(y_n\approx y(t_0+nh)\) y\(y_{n+1}\approx y(t+(n+1)h)\text{,}\) entonces nosotros quiero que el error de truncamiento local en el paso de\(y_n\)\(y_{n+1}\) a no sea más que sobre\(\varepsilon h\text{.}\) Supongamos además que ya hemos producido la solución aproximada hasta donde\(t_n\text{.}\) La estrategia aproximada es la siguiente. Hacemos el paso de\(t_n\) a dos\(t_n+h\) veces usando dos algoritmos diferentes, dando dos aproximaciones diferentes a lo\(y(t_{n+1})\text{,}\) que llamamos\(A_{1,n+1}\) y\(A_{2,n+1}\text{.}\) Los dos algoritmos se eligen para que

    1. podemos usar\(A_{1,n+1}-A_{2,n+1}\) para calcular un error de truncamiento local aproximado y
    2. para la eficiencia, los dos algoritmos utilizan casi las mismas evaluaciones de\(f\text{.}\) Recuerde que evaluar la función\(f\) suele ser la parte que más tiempo consume de nuestro cálculo.

    En el caso de que el error de truncamiento local, dividido por\(h\text{,}\) (es decir, el error por unidad de incremento de\(t\)) sea menor de lo que\(\varepsilon\text{,}\) establecemos\(t_{n+1}=t_n+h\text{,}\) aceptar\(A_{2,n+1}\) como el valor aproximado 2 Mejor aún, aceptar\(A_{2,n+1}\) minus the computed approximate error in \(A_{2,n+1}\) as the approximate value for \(y(t_{n+1})\text{.}\) para\(y(t_{n+1})\text{,}\) y pasar a la siguiente paso. De lo contrario, elegimos, usando lo que hemos aprendido de\(A_{1,n+1}-A_{2,n+1}\text{,}\) un nuevo tamaño de paso de prueba\(h\) y comenzamos de nuevo en\(t_n\text{.}\)

    Ahora para los detalles. Empezamos con un procedimiento muy sencillo. Posteriormente lo haremos una sopa para conseguir un procedimiento mucho más eficiente.

    Euler y Euler-2step (versión preliminar)

    Denotar por\(\phi(t)\) la solución exacta a\(y'=f(t,y)\) que satisface la condición inicial\(\phi(t_n)=y_n\text{.}\) Si aplicamos un paso de Euler con tamaño de paso\(h\text{,}\) dando

    \[ A_{1,n+1}=y_n+hf(t_n,y_n) \nonumber \]

    sabemos, por (D.2.4), que

    \[ A_{1,n+1}=\phi(t_n+h)+Kh^2+O(h^3) \nonumber \]

    El problema, claro, es que no sabemos cuál es el error, ni siquiera aproximadamente, porque no sabemos cuál\(K\) es la constante. Pero podemos estimar\(K\) simplemente rehaciendo el paso de\(t_n\) a\(t_n+h\) usando un segundo algoritmo elegido juiciosamente. Hay una serie de segundos algoritmos diferentes que funcionarán. Llamamos al algoritmo simple que usamos en esta subsección EULER-2Step 3 Este nombre está pidiendo una nota al pie de página relacionada con la danza e invitamos al lector a que suministre la suya propia. . Un paso de Euler-2step con tamaño de paso\(h\) solo consiste en hacer dos pasos de Euler con tamaño de paso\(h/2\text{:}\)

    \[ A_{2,n+1} = y_n+\tfrac{h}{2}f(t_n,y_n) +\tfrac{h}{2}f\big(t_n+\tfrac{h}{2},y_n+\tfrac{h}{2}f(t_n,y_n)\big) \nonumber \]

    Aquí, el primer medio paso nos llevó de\(y_n\) a\(y_{\rm mid}=y_n+\frac{h}{2}f(t_n,y_n)\) y el segundo medio paso nos llevó de\(y_{\rm mid}\) a\(y_{\rm mid}+\frac{h}{2}f\big(t_n+\frac{h}{2},y_{\rm mid}\big)\text{.}\) El error de truncamiento local introducido en el primer medio paso es\(K(h/2)^2+O(h^3)\text{.}\) Eso para el segundo medio paso es\(K(h/2)^2+O(h^3)\) con el mismo 4 Porque los dos medios -los pasos comienzan en valores de\(t\) only \(h/2\) apart, and we are thinking of \(h\) as being very small, it should not be surprising that we can use the same value of \(K\) in both. In case you don't believe us, we have included a derivation of the local truncation error for Euler-2step later in this appendix.\(K\text{,}\) aunque con un\(O(h^3)\text{.}\) All together diferente

    \ begin {alinear*} A_ {2, n+1} &=\ phi (t_n+h) +\ grande [K\ grande (\ tfrac {h} {2}\ grande) ^2+O (h^3)\ grande] +\ grande [K\ grande (\ tfrac {h} {2}\ grande) ^2+O (h^3)\ grande]\ &= phi (t_n+h) +\ mitad Kh^2+O (h^3)\ final {alinear*}

    La diferencia es 5

    \ begin {align*} A_ {1, n+1} -A_ {2, n+1} &=\ grande [\ phi (t_n+h) +Kh^2+O (h^3)\ grande] -\ grande [\ phi (t_n+h) -\ mitad Kh^2-O (h^3)\ grande]\\ &=\ mitad Kh^2+O (h^3)\ fin {alinear*}

    Entonces, si hacemos un paso tanto de Euler como de Euler-2step, podemos estimar

    \[ \half Kh^2=A_{1,n+1}-A_{2,n+1}+O(h^3) \nonumber \]

    Ahora sabemos que en el paso recién completado EULER-2step introdujo un error de aproximadamente Es\(\half Kh^2\approx A_{1,n+1}-A_{2,n+1}\text{.}\) decir, la tasa de error actual es aproximadamente\(r=\frac{|A_{1,n+1}-A_{2,n+1}|}{h}\approx\half |K| h\) por unidad de aumento de\(t\text{.}\)

    • Si esto\(r=\frac{|A_{1,n+1}-A_{2,n+1}|}{h}>\varepsilon\text{,}\) rechazamos 6 La tasa de error medida,\(r\text{,}\) is bigger than the desired error rate \(\varepsilon\text{.}\) That means that it is harder to get the accuracy we want than we thought. So we have to take smaller steps.\(A_{2,n+1}\) y repetimos el paso actual con un nuevo tamaño de paso de prueba elegido de manera que\(\half |K|(\text{new }h)\lt\varepsilon\text{,}\) es decir,\(\frac{r}{h}(\text{new }h)\lt\varepsilon\text{.}\) Para darnos un pequeño margen de seguridad, podríamos usar 7 No queremos hacer el nuevo\(h\) too close to \(\frac{\varepsilon}{r}{h}\) since we are only estimating things and we might end up with an error rate bigger that \(\varepsilon\text{.}\) On the other hand, we don't want to make the new \(h\) too small because that means too much work — so we choose it to be just a little smaller than \(\frac{\varepsilon}{r}{h}\) \(\ldots\) say \(0.9\frac{\varepsilon}{r}{h}\).

      \[ \text{new }h=0.9\,\frac{\varepsilon}{r}\,h \nonumber \]

    • Si\(r=\frac{|A_{1,n+1}-A_{2,n+1}|}{h}\lt\varepsilon\) podemos aceptar 8 La tasa de error medida,\(r\text{,}\) is smaller than the desired error rate \(\varepsilon\text{.}\) That means that it is easier to get the accuracy we want than we thought. So we can make the next step larger.\(A_{2,n+1}\) como valor aproximado para\(y(t_{n+1})\text{,}\) con\(t_{n+1}=t_n+h\text{,}\) y pasar al siguiente paso, comenzando con el nuevo tamaño de paso de prueba 9 Tenga en cuenta que en este caso\(\frac{\varepsilon}{r}>1\text{.}\) So the new \(h\) can be bigger than the last \(h\text{.}\)

      \[ \text{new } h=0.9\,\frac{\varepsilon}{r}\,h \nonumber \]

    Esa es nuestra versión preliminar del método de tamaño de paso variable Euler/Euler-2step. Lo llamamos la versión preliminar, porque en breve la vamos a ajustar para obtener un procedimiento mucho más eficiente.

    Ejemplo D.3.1

    Como ejemplo concreto, supongamos que nuestro problema es

    \[ y(0)=e^{-2},\ y'=8(1-2t)y,\ \varepsilon=0.1 \nonumber \]

    y que hemos llegado hasta

    \[ t_n=0.33,\ y_n=0.75,\ \ \text{trial }h=0.094 \nonumber \]

    Luego, utilizando\(E=|A_{1,n+1}-A_{2,n+1}|\) para denotar la magnitud del error de truncamiento local estimado en\(A_{2,n+1}\) y\(r\) la tasa de error correspondiente

    \ begin {align*} f (t_n, y_n) &=8 (1-2\ times 0.33) 0.75=2.04\\ A_ {1, n+1} &=y_n+hf (t_n, y_n) =0.75+0.094\ times 2.04=0.942\ y_ {\ rm mid} &=y_n+\ tfrac {h} {2} f (t_n mid} &=y_n+\ tfrac {h} {2} f (t_n, y_n) =0.75+\ tfrac {0.094} {2}\ times 2.04=0.846\\ f\ grande (t_n+\ tfrac {h} {2}, y_ {\ rm mid}\ grande) &=8\ grande [1-2\ grande (0.33+\ tfrac {0.094} {2}\ grande)\ grande] 0.846=1.66\ \ A_ {2, n+1} &=y_ {\ rm mid} +\ tfrac {h} {2} f (t_n+\ tfrac {h} {2}, y_ {\ rm mid}) =0.846+\ tfrac {0.094} {2} 1.66=0.924\ E&=|A_ {1, n+1} -A_ {2, n+1} |=|0.942-0.924|=0.018\\ r&=\ frac {|E|} {h} =\ frac {0.018} {0.094} =0.19\ end {align*}

    Dado que\(r=0.19>\varepsilon=0.1\,\text{,}\) el tamaño de paso actual es inaceptable y tenemos que volver a calcular con el nuevo tamaño de paso

    \[ \text{new } h=0.9\frac{\varepsilon}{r}(\text{old }h) =0.9\ \frac{0.1}{0.19}\ 0.094 =0.045 \nonumber \]

    para dar

    \ begin {align*} f (t_n, y_n) &=8 (1-2\ tiempo0.33) 0.75=2.04\\ A_ {1, n+1} &=y_n+hf (t_n, y_n) =0.75+0.045\ times 2.04=0.842\ y_ {\ rm mid} &=y_n+\ tfrac {h} {2} f (t_n, y_n) =0.75+\ tfrac {0.045} {2}\ times 2.04=0.796\\ f\ big (t_n+\ tfrac {h} {2}, y_ {\ rm mid}\ big) &=8\ Grande [1-2\ grande (0.33+\ tfrac {0.045} {2}\ big)\ Grande] 0.796=1.88\ A_ {2, n+1} &=y_ {\ text {mid}} +\ tfrac {h} {2} f (t_n +\ tfrac {h} {2}, y_ {\ rm mid}) =0.796+\ tfrac {0.045} {2} 1.88 =0.838\\ E&=A_ {1.n+1} -A_ {2.n+1} =0.842-0.838=0.004\\ r&=\ frac {|E|} {h} =\ frac {0.004} {0.045} =0.09\ end {align*}

    Esta vez\(\,r=0.09\lt \varepsilon=0.1\,\text{,}\) es aceptable así que establecemos\(t_{n+1}=0.33+0.045=0.375\) y

    \[ y_{n+1}=A_{2,n+1}=0.838 \nonumber \]

    El tamaño inicial del paso de prueba de\(t_{n+1}\) a\(t_{n+2}\) es

    \[ \text{new }h = 0.9\frac{\varepsilon}{r}(\text{old }h) =0.9\,\frac{0.1}{0.09}\,.045=.045 \nonumber \]

    Por casualidad, ha resultado que lo nuevo\(h\) es lo mismo que el viejo\(h\) (a tres decimales). Si\(r\) hubiera sido significativamente más pequeño que\(\varepsilon\text{,}\) entonces el nuevo\(h\) habría sido significativamente más grande que el antiguo\(h\), lo que indica que es (relativamente) fácil estimar cosas en esta región, haciendo suficiente un tamaño de paso más grande.

    Como dijimos anteriormente, en breve actualizaremos el método de tamaño de paso variable anterior, que estamos llamando a la versión preliminar del método Euler/EUER-2step, para obtener un procedimiento mucho más eficiente. Antes de hacerlo, hagamos una pausa para investigar un poco qué tan bien funciona nuestro procedimiento preliminar para controlar la tasa de producción de errores.

    Nos hemos venido refiriendo, vagamente,\(\varepsilon\) como la tasa deseada para la introducción del error, por nuestro método de tamaño de paso variable, a medida que\(t\) avanza. Si la tasa de aumento del error fuera exactamente\(\varepsilon\text{,}\) entonces en\(t_f\) el momento final el error acumulado sería exactamente\(\varepsilon(t_f-t_0)\text{.}\) Pero nuestro algoritmo en realidad elige el tamaño del paso\(h\) para cada paso para que el error de truncamiento local estimado en\(A_{2,n+1}\) para ese paso sea sobre\(\varepsilon h\text{.}\) Nosotros han visto que, una vez introducido algún error de truncamiento local, su contribución al error de truncamiento global puede crecer exponencialmente con\(t_f\text{.}\)

    Aquí están los resultados de un experimento numérico que ilustra este efecto. En este experimento, el método Euler/Euler-2step preliminar anterior se aplica al problema de valor inicial\(\ y'=t-2y,\ y(0)=3 \ \) para\(\varepsilon=\frac{1}{16},\frac{1}{32},\cdots\) (diez valores diferentes) y para\(t_f=0.2,\ 0.4,\ \cdots,\ 3.8\text{.}\) Aquí hay una gráfica del resultante\(\frac{\text{actual error at }t=t_f}{\varepsilon t_f}\) contra\(t_f\text{.}\)

    vble_error.svg

    Si la tasa de introducción de error fuera exactamente\(\varepsilon\text{,}\) tendríamos\(\frac{\text{actual error at }t=t_f}{\varepsilon t_f}=1\text{.}\) Hay un pequeño cuadrado en la gráfica para cada par diferente\(\varepsilon,t_f\text{.}\) Entonces para cada valor de\(t_f\) hay diez cuadrados (posiblemente superpuestos) en la línea\(x=t_f\text{.}\) Este experimento numérico sugiere que\(\frac{\text{actual error at }t=t_f}{\varepsilon t_f}\) es relativamente independiente\(\varepsilon\) y comienza, cuando\(t_f\) es pequeño, en aproximadamente uno, como queremos, pero crece (quizás exponencialmente) con\(t_f\text{.}\)

    Euler y Euler-2step (versión final)

    Ahora estamos listos para usar un poco astuta de arithemtic para potenciar nuestro método Euler/Euler-2step. Al igual que en nuestro desarrollo de la versión preliminar del método, denotar por\(\phi(t)\) la solución exacta a\(y'=f(t,y)\) que satisface la condición inicial\(\phi(t_n)=y_n\text{.}\) Hemos visto, al inicio del § D.3.1, que aplicar un paso de Euler con tamaño de paso\(h\text{,}\) da

    \ begin {align*} A_ {1, n+1} &=y_n+hf (t_n, y_n)\ noetiqueta\\ &=\ phi (t_n+h) +Kh^2+O (h^3)\ tag {E5}\ end {align*}

    y aplicar un paso de Euler-2step con tamaño de paso\(h\) (es decir, aplicar dos pasos de Euler con tamaño de paso\(h/2\)) da

    \ begin {alinear*} A_ {2, n+1} &= y_n+\ tfrac {h} {2} f (t_n, y_n) +\ tfrac {h} {2} f\ grande (t_n+\ tfrac {h} {2}\,,\, y_n+\ tfrac {h} {2} f (t_n, y_n)\ grande)\ noetiqueta\\ &=\ phi (t_n+h) +\ tfrac {1} {2} Kh^2+O (h^3)\ tag {E6}\ end {align*}

    porque el error de truncamiento local introducido en el primer medio paso fue\(K(h/2)^2+O(h^3)\) y el introducido en el segundo medio paso fue\(K(h/2)^2+O(h^3)\text{.}\) Ahora aquí está el bocado furtivo. Las ecuaciones (E5) y (E6) son muy similares y podemos eliminar todas las\(Kh^2\) restando (E5) de los\(2\) tiempos (E6). Esto da

    \[ 2\text{(E6)} - \text{(E5):}\qquad 2A_{2,n+1}-A_{1,n+1} = \phi(t_n+h) +O(h^3) \nonumber \]

    (¡no más\(h^2\) término!) o

    \[ \phi(t_n+h)= 2A_{2,n+1}-A_{1,n+1}+O(h^3) \tag{E7} \nonumber \]

    lo que nos dice que elegir

    \[ y_{n+1}=2A_{2,n+1}-A_{1,n+1} \tag{E8} \nonumber \]

    daría un error de orden de truncamiento local\(h^3\text{,}\) en lugar del orden\(h^2\) del método preliminar Euler/Euler-2step. Para convertir la versión preliminar a la versión final, solo reemplazamos\(y_{n+1}=A_{2,n+1}\) por\(y_{n+1} = 2A_{2,n+1}-A_{1,n+1}\text{:}\)

    Ecuación D.3.2 Método Euler/Euler-2step

    Dado\(\varepsilon>0\text{,}\)\(t_n\text{,}\)\(y_n\) y el tamaño de paso actual\(h\)

    • computar

      \ begin {align*} A_ {1, n+1} &=y_n+hf (t_n, y_n)\\ A_ {2, n+1} &= y_n+\ tfrac {h} {2} f (t_n, y_n) +\ tfrac {h} {2} f\ big (t_n+\ tfrac {h} {2}\,, y_n n+\ tfrac {h} {2} f (t_n, y_n)\ grande)\\ r&=\ frac {|A_ {1, n+1} -A_ {2, n+1} |} {h}\ final {alinear*}

    • Si\(r>\varepsilon\text{,}\) repite la primera bala pero con el nuevo tamaño de paso

      \[ (\text{new }h)=0.9\,\frac{\varepsilon}{r}\,(\text{old }h) \nonumber \]

    • Si se\(r\lt\varepsilon\) establece

      \ begin {align*} t_ {n+1} &=t_n+h\\ y_ {n+1} &=2A_ {2, n+1} -A_ {1, n+1}\ quad\ text {y el nuevo tamaño del paso de prueba}\\ (\ text {new} h) &=0.9\,\ frac {\ varepsilon} {r}\, (\ text {old} h)\ final {alinear*}

      y pasar al siguiente paso. Tenga en cuenta\(r\lt\varepsilon\text{,}\)\(\frac{r}{\varepsilon}h\gt h\) que ya que lo que indica que el nuevo\(h\) puede ser más grande que el\(0.9\) viejo\(h\text{.}\) Incluimos el de tener cuidado de no hacer demasiado grande el error del siguiente paso.

    Pensemos un poco sobre cómo debería funcionar nuestro método final Euler/Euler-2step.

    • El tamaño del paso aquí, como en la versión preliminar, se elige de manera que el error de truncamiento local en el aumento\(A_{2,n+1}\) por unidad de\(t\text{,}\) a saber\(r=\frac{|A_{1,n+1}-A_{2,n+1}|}{h}\approx\frac{Kh^2/2}{h}=\frac{K}{2}h\text{,}\) sea aproximadamente\(\varepsilon\text{.}\) So\(h\) es aproximadamente proporcional a\(\varepsilon\text{.}\)
    • Por otro lado, (E7) muestra que, en el método completo, se está sumando el error de truncamiento local\(y_{n+1}\) a una tasa de incremento\(\frac{O(h^3)}{h}=O(h^2)\) por unidad en\(t\text{.}\)
    • Entonces uno esperaría que el truncamiento local incremente el error a una tasa proporcional al incremento\(\varepsilon^2\) por unidad en\(t\text{.}\)
    • Si la tasa de incremento del error fuera exactamente un tiempo constante\(\varepsilon^2\text{,}\) entonces el error acumulado entre el tiempo inicial\(t=0\) y el tiempo final\(t=t_f\) sería exactamente un tiempo constante\(\varepsilon^2\,t_f\text{.}\)
    • Sin embargo hemos visto que, una vez que se ha introducido algún error de truncamiento local, su contribución al error global puede crecer exponencialmente con\(t_f\text{.}\) Así que esperaríamos que, bajo el método completo de Euler/Euler-2step,\(\frac{{\rm actual\ error\ at\ }t=t_f}{\varepsilon^2 t_f}\) sea más o menos independiente de\(\varepsilon\text{,}\) pero aún creciendo exponencialmente en\(t_f\text{.}\)

    Aquí están los resultados de un experimento numérico que ilustra esto. En este experimento, el método final Euler/Euler-2step anterior, (D.3.2), se aplica al problema de valor inicial\(\ y'=t-2y,\ y(0)=3 \ \) para\(\varepsilon=\frac{1}{16},\frac{1}{32},\cdots\) (diez valores diferentes) y para\(t_f=0.2,\ 0.4,\ \cdots,\ 3.8\text{.}\) En la siguiente gráfica, hay un pequeño cuadrado para el resultante\(\frac{\text{actual error at }t=t_f}{\varepsilon^2 t_f}\) para cada par diferente\(\varepsilon,t_f\text{.}\)

    vble4_error.svg

    De hecho, parece que\(\frac{{\rm actual\ error\ at\ }t=t_f}{\varepsilon^2 t_f}\) es relativamente independiente\(\varepsilon\) pero crece (quizás exponencialmente) con\(t_f\text{.}\) Note que\(\frac{{\rm actual\ error\ at\ }t=t_f}{\varepsilon^2 t_f}\) contiene un factor de\(\varepsilon^2\) en el denominador. La tasa de error real\(\frac{{\rm actual\ error\ at\ }t=t_f}{t_f}\) es mucho menor de lo que sugiere la gráfica.

    Método de Fehlberg

    Por supuesto, en la práctica métodos más eficientes y más precisos 10 Hay un número muy grande de tales métodos. Solo veremos brevemente un par de los más simples. El lector interesado puede encontrar más mediante el motor de búsqueda de palabras clave como “métodos Runge-Kutta” y “tamaño de paso adaptativo”. que Euler y Euler-2step se utilizan. El método de Fehlberg 11 E. Fehlberg, el Informe Técnico de la NASA R315 (1969) y el Informe Técnico de la NASA R287 (1968). utiliza Euler mejorado y un segundo método más preciso. Cada paso implica tres cálculos de\(f\text{:}\)

    \ begin {align*} f_ {1, n} &=f (t_n, y_n)\\ f_ {2, n} &=f (t_n+h, y_n+hf_ {1, n})\\ f_ {3, n} &=f\ left (t_n+\ tfrac {h} {2}, y_n+\ tfrac {h} {4} [f_ {1, n} +f_ {2, n}]\ derecha)\ final {alinear*}

    Una vez realizadas estas tres evaluaciones, el método genera dos aproximaciones para\(y(t_n+h)\text{:}\)

    \ begin {alinear*} A_ {1, n+1} &=y_n+\ tfrac {h} {2}\ izquierda [f_ {1, n} +f_ {2, n}\ derecha]\\ A_ {2, n+1} &=y_n+\ tfrac {h} {6}\ izquierda [f_ {1, n} +f_ {2, n} 4f_ {3, n}\ derecha]\ final {alinear*}

    Denotar por\(\phi(t)\) la solución exacta a\(y'=f(t,y)\) que satisface la condición inicial\(\phi(t_n)=y_n\text{.}\) Ahora\(A_{1,n+1}\) es solo el\(y_{n+1}\) producido por el método mejorado de Euler. El error de truncamiento local para el método mejorado de Euler es de orden\(h^3\text{,}\) uno de potencia\(h\) menor que el del método de Euler. Entonces

    \[ A_{1,n+1} = \phi(t_n+h) + Kh^3+O(h^4) \nonumber \]

    y resulta 12 El lector interesado puede encontrar el papel original de Fehlberg en línea (¡en la NASA!) y seguir la derivación. Requiere cuidadosas expansiones de Taylor y luego álgebra inteligente para cancelar los términos de error más grandes. que

    \[ A_{2,n+1} = \phi(t_n+h) + O(h^4) \nonumber \]

    Entonces el error en\(A_{1,n+1}\) es

    \ begin {align*} E&=\ big|kh^3+o (h^4)\ big| =\ big|a_ {1, n+1} -\ phi (t_n+h)\ big|+O (h^4)\\ &=\ big|a_ {1, n+1} -A_ {2, n+1}\ big|+O (h^4)\ end {align*}

    y nuestra estimación para la tasa a la que se introduce el error\(A_{1,n+1}\) es

    \[ r=\frac{|A_{1,n+1}-A_{2,n+1}|}{h}\approx |K|h^2 \nonumber \]

    por unidad de aumento de\(t\text{.}\)

    • Si\(r>\varepsilon\) rehacemos este paso con un nuevo tamaño de paso de prueba elegido para que,\(|K|(\text{new }h)^2\lt\varepsilon\text{,}\) es decir,\(\frac{r}{h^2}(\text{new }h)^2\lt\varepsilon\text{.}\) Con nuestro factor de seguridad tradicional, tomamos

      \[ \text{new }h=0.9\sqrt{\frac{\varepsilon}{r}}\,h\qquad \text{(the new }h\text{ is smaller)} \nonumber \]

    • Si\(r\le \varepsilon\) configuramos\(t_{n+1}=t_n+h\) y\(y_{n+1}=A_{2,n+1}\) (ya que\(A_{2,n+1}\) debería ser considerablemente más preciso que\(A_{1,n+1}\)) y pasamos al siguiente paso con el tamaño del paso de prueba

      \[ \text{new }h=0.9\sqrt{\frac{\varepsilon}{r}}\,h\qquad \text{(the new }h\text{ is usually bigger)} \nonumber \]

    El proceso Kutta-Merson

    El proceso Kutta-Merson 13 R.H. Merson, ``Un método operativo para el estudio de los procesos de integración”, Proc. Symp. Procesamiento de Datos, Armas Res. Establo. Salisbury, Salisbury (1957) pp. 110—125. utiliza dos variaciones del método Runge-Kutta. Cada paso implica cinco cálculos 14 Al igual que los otros métodos descritos anteriormente, los coeficientes\(1/3\text{,}\) \(1/6\text{,}\) \(1/8\) etc. are chosen so as to cancel larger error terms. While determining the correct choice of coefficients is not conceptually difficult, it does take some work and is beyond the scope of this appendix. The interested reader should search-engine their way to a discussion of adaptive Runge-Kutta methods. de\(f\text{:}\)

    \ begin {alinear*} k_ {1, n} &=f (t_n, y_n)\\ k_ {2, n} &=f\ grande (t_n+\ tfrac {1} {3} h, y_n+\ tfrac {1} {3} hk_ {1, n}\ grande)\ k_ {3, n} &= f\ grande (t_n+\ tfrac {1} {3} h, y_n+\ tfrac {1} {6} hk_ {1, n} +\ tfrac {1} {6} hk_ {2, n}\ grande)\\ k_ {4, n} &= f\ grande (t_n+\ tfrac {1} {2} h, y_n+\ tfrac {1} 8} hk_ {1, n} +\ tfrac {3} {8} hk_ {3, n}\ grande)\\ k_ {5, n} &= f\ grande (t_n+h, y_n+\ tfrac {1} {2} hk_ {1, n} -\ tfrac {3} {2} hk_ {3, n} +2hk_ {4, n}\ grande)\ final {alinear*}

    Una vez realizadas estas cinco evaluaciones, el proceso genera dos aproximaciones para\(y(t_n+h)\text{:}\)

    \ begin {align*} A_ {1, n+1} &=y_n+h\ izquierda [\ tfrac {1} {2} k_ {1, n} -\ tfrac {3} {2} k_ {3, n} +2k_ {4, n}\ derecha]\\ A_ {2, n+1} &=y_n+h\ izquierda [\ tfrac {1} {6} k_ {1, n} +\ tfrac {2} {3} k_ {4, n} +\ tfrac {1} {6} k_ {5, n}\ derecha]\ final {alinear*}

    El error (firmado) en\(A_{1,n+1}\) es\(\frac{1}{120}h^5K+O(h^6)\) mientras que en\(A_{2,n+1}\) está\(\frac{1}{720}h^5K+O(h^6)\) con la misma constante\(K\text{.}\) Así\(A_{1,n+1}-A_{2,n+1} = \frac{5}{720}Kh^5+O(h^6)\) y la constante desconocida se\(K\) puede determinar, dentro de un error\(O(h)\text{,}\) por

    \[ K=\frac{720}{5\,h^5}(A_{1,n+1}-A_{2,n+1}) \nonumber \]

    y el error aproximado (firmado) en\(A_{2,n+1}\) y su tasa correspondiente por incremento unitario de\(t\) son

    \ begin {alinear*} E&=\ frac {1} {720} K h^5=\ frac {1} {5} (A_ {1, n+1} -A_ {2, n+1})\\ r=\ frac {|E|} {h} &=\ frac {1} {720} |K| h^4=\ frac {1} {5\, h}\ big|a_ {1, n+1} -A_ {2, n+1}\ big|\ end {alinear*}

    • Si\(r>\varepsilon\) rehacemos este paso con un nuevo tamaño de paso de prueba elegido para que,\(\frac{1}{720}|K|(\text{new h})^4\lt\varepsilon\text{,}\) es decir,\(\frac{r}{h^4}(\text{new }h)^4\lt\varepsilon\text{.}\) Con nuestro factor de seguridad tradicional, tomamos

      \[ \text{new }h=0.9\left(\frac{\varepsilon}{r}\right)^{1/4}\,h \nonumber \]

    • Si\(r\le \varepsilon\) establecemos\(t_{n+1}=t_n+h\) y\(y_{n+1}=A_{2,n+1}-E\) (ya que\(E\) es nuestra estimación del error firmado en\(A_{2,n+1}\)) y pasamos al siguiente paso con tamaño de paso de prueba

      \[ \text{new }h=0.9\left(\frac{\varepsilon}{r}\right)^{1/4}\,h \nonumber \]

    El error de truncamiento local para EULER-2step

    En nuestra descripción de Euler/euler-2step anterior simplemente declaramos el error de truncamiento local sin una explicación. En esta sección, mostramos cómo se puede derivar. Observamos que cálculos muy similares sustentan los otros métodos que hemos descrito.

    En esta sección, estaremos usando derivadas parciales y, en particular, la regla de cadena para funciones de dos variables. Ese material está cubierto en el Capítulo 2 del texto CLP-3. Si aún no te sientes cómodo con él, puedes tomar nuestra palabra por esos bits, o puedes retrasar la lectura de esta sección hasta que hayas aprendido un poco de cálculo multivariable.

    Recordemos que, por definición, el error de truncamiento local para un algoritmo es el error (firmado) generado por un solo paso del algoritmo, bajo los supuestos de que iniciamos el paso con la solución exacta y que no hay error de redondoff 15 Debemos señalar que en serio grande cálculos numéricos, uno realmente tiene que tomar en cuenta los errores de redondeo porque pueden causar serios problemas. El lector interesado debe buscar su camino hacia la historia de las simulaciones numéricas de Edward Lorenz y los inicios de la teoría del caos. Desafortunadamente simplemente no tenemos espacio en este texto para discutir todos los aspectos de las matemáticas. . Denote por\(\phi(t)\) la solución exacta a

    \ begin {align*} y' (t) &=f (t, y)\\ y (t_n) &=y_n\ end {align*}

    En otras palabras,\(\phi(t)\) obedece

    \ begin {align*}\ phi' (t) &= f\ big (t,\ phi (t)\ big)\ qquad\ text {para todos} t\\ phi (t_n) &=y_n\ end {align*}

    En particular\(\phi'(t_n)=f\big(t_n,\phi(t_n)\big)=f(t_n,y_n)\) y, utilizando cuidadosamente la regla de la cadena, que es (2.4.2) en el texto CLP-3,

    \ begin {align*}\ phi "(t_n) &=\ dfrac {d} {dt} f\ big (t,\ phi (t)\ big)\ big|_ {t=t_n} =\ grande [f_t\ big (t,\ phi (t)\ big) +f_y\ big (t,\ phi (t)\ big)\ phi' (t) Big\] _ {t=t_n}\\ &=f_t (t_n, y_n) +f_y (t_n, y_n)\, f (t_n, y_n)\ tag {E9}\ final {alinear*}

    Recuerde que\(f_t\) es la derivada parcial de\(f\) con respecto a\(t\text{,}\) y esa\(f_y\) es la derivada parcial de\(f\) con respecto a\(y\text{.}\) Vamos a necesitar (E9) a continuación.

    Por definición, el error de truncamiento local para Euler es

    \[ E_1(h)=\phi(t_n+h)-y_n-hf\big(t_n,y_n\big) \nonumber \]

    mientras que para EULER-2step es

    \[ E_2(h)=\phi(t_n+h)-y_n-\tfrac{h}{2}f(t_n,y_n) -\tfrac{h}{2}f\big(t_n+\tfrac{h}{2},y_n+\tfrac{h}{2}f(t_n,y_n)\big) \nonumber \]

    Para entender cómo\(E_1(h)\) y\(E_2(h)\) comportarnos para los pequeños\(h\) podemos usar las expansiones de Taylor ((3.4.10) en el texto CLP-1) para escribirlas como series de poder en\(h\text{.}\) Para ser precisos, usamos

    \[ g(h)=g(0)+g'(0)\,h+\half g''(0)\,h^2+O(h^3) \nonumber \]

    para expandir ambos\(E_1(h)\) y\(E_2(h)\) en poderes de\(h\) ordenar\(h^2\text{.}\) Tenga en cuenta que, en la expresión para\(E_1(h)\text{,}\)\(t_n\) y\(y_n\) son constantes — no varían con\(h\text{.}\) Así que calcular derivados de\(E_1(h)\) con respecto a\(h\) es en realidad bastante simple.

    \ begin {alignat*} {2} E-1 (h) &=\ phi (t_n+h) -y_n-hf\ grande (t_n, y_n\ grande)\ qquad y E-1 (0) &=\ phi (t_n) -y_n=0\\ E'_1 (h) &=\ phi' (t_n+h) -f\ big (t_n, y_n grande\) & E'_1 (0) &=\ phi' (t_n) -f\ grande (t_n, y_n\ grande) =0\\ E"_1 (h) &=\ phi "(t_n+h) & E"_1 (0) &=\ phi" (t_n)\ end {alignat*}

    Por Taylor, el error de truncamiento local para Euler obedece

    Ecuación D.3.3

    \[ E_1(h)=\half\phi''(t_n)h^2+O(h^3)=Kh^2+O(h^3)\qquad\text{with } K=\half\phi''(t_n) \nonumber \]

    Computar argumentos de\(E_2(h)\) con respecto a\(h\) es un poco más difícil, ya que\(h\) ahora aparece en los argumentos de la función\(f\text{.}\) Como consecuencia, tenemos que incluir algunas derivadas parciales.

    \ begin {align*} E_2 (h) &=\ phi (t_n+h) -y_n-\ frac {h} {2} f (t_n, y_n) -\ frac {h} {2} f\ grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ grande)\ E'_2 (h) &=\ phi' (t_n+h) -\ frac {1} {2} f (t_n, y_n) -\ frac {1} {2} f\ grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ grande)\\ &\ hskip2.5in-\ frac {h} {2}\ underbrackets {\ dfrac {d} {dh} f\ Grande (t_n+ \ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ Grande)} _ {\ text {deja esta expresión como está por ahora}}\\ E"_2 (h) &=\ phi "(t_n+h) -2\ veces\ frac {1} {2}\ underbrackets {\ dfrac {d} {dh} f\ Grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ Grande)} _ {\ text {deja este también}}\\ &\ hskip2.5in-\ frac {h} {2}\ underbrackets {\ frac {\ mathrm {d^2}} {\ mathrm {d} h^2} f\ Grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ grande)} _ {\ text {y deja este también}}\ end {alinear*}

    Ya que solo necesitamos\(E_2(h)\) y sus derivados en no\(h=0\text{,}\) tenemos que computar el\(\frac{\mathrm{d^2 f}}{\mathrm{d}h^2}\) término (por suerte) y tampoco necesitamos computar el\(\dfrac{df}{dh}\) término en\(E_2'\text{.}\) Hacemos, sin embargo,\(\dfrac{df}{dh}\Big|_{h=0}\) necesidad de\(E_2''(0)\text{.}\)

    \ begin {align*} E_2 (0) &=\ phi (t_n) -y_n=0\\ E'_2 (0) &=\ phi' (t_n) -\ frac {1} {2} f (t_n, y_n) -\ frac {1} {2} f (t_n, y_n) =0\\ E"_2 (0) &=\ phi "(t_n) -\ dfrac {d} {dh} f\ Grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ Grande)\ Big|_ {h=0}\\ &=\ phi" (t_n) -\ frac {1} {2} f_t\ Grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ Grande)\ Gran|_ {h=0}\\ &\ hskip1.5in -\ frac {1} {2} f (t_n, y_n)\, f_y\ Grande (t_n+\ frac {h} {2}, y_n+\ frac {h} {2} f (t_n, y_n)\ Grande)\ Big|_ {h=0}\\ &= phi\ "(t_n) -\ frac {1} {2} f_t (t_n, y_n) -\ frac {1} {2} f_y (t_n, y_n)\, f (t_n, y_n)\\ &=\ mitad\ phi" (t_n)\ qquad\ text {por (E9)}\ end {align*}

    Por Taylor, el error de truncamiento local para EULER-2step obedece

    Ecuación D.3.4

    \[ E_2(h)=\frac{1}{4}\phi''(t_n)\,h^2+O(h^3) =\half Kh^2+O(h^3)\qquad{\rm with\ }K=\half\phi''(t_n) \nonumber \]

    Observe que el\(K\) in (D.3.4) es idéntico al\(K\) in (D.3.3). Esto es exactamente lo que necesitábamos en nuestro análisis de las Secciones D.3.1 y D.3.2.


    This page titled D.3: Métodos de Tamaño de Paso Variable is shared under a CC BY-NC-SA 4.0 license and was authored, remixed, and/or curated by Joel Feldman, Andrew Rechnitzer and Elyse Yeager via source content that was edited to the style and standards of the LibreTexts platform; a detailed edit history is available upon request.