Métodos de Runge-Kutta
Los métodos de Runge-Kutta son métodos de aproximación de soluciones de análisis numérico de ecuaciones diferenciales . Fueron nombrados en honor a los matemáticos Carl Runge y Martin Wilhelm Kutta , quienes desarrollaron el método en 1901.
Estos métodos se basan en el principio de iteración , es decir, se utiliza una primera estimación de la solución para calcular una segunda estimación más precisa, etc.
Principio general
Considere el siguiente problema:
y′=F(t,y),y(t0)=y0{\ Displaystyle y '= f (t, y), \ quad y (t_ {0}) = y_ {0}}
vamos a ver a resolver en un conjunto discreto t 0 < t 1 <... < t N . Más que buscar un método directo, los métodos de Runge-Kutta proponen introducir los puntos intermedios para calcular por inducción los valores ( t n , y n ) con
{(tno,I,yno,I)}1⩽I⩽q{\ Displaystyle \ {(t_ {n, i}, y_ {n, i}) \} _ {1 \ leqslant i \ leqslant q}}
tno,I=tno+vsI⋅hno{\ Displaystyle t_ {n, i} = t_ {n} + c_ {i} \ cdot h_ {n}}
donde h n = t n + 1 - t n es el paso de tiempo y c i está en el intervalo [0; 1] . Para cada punto intermedio, notamos la pendiente correspondiente
pagno,I=F(tno,I,yno,I).{\ Displaystyle p_ {n, i} = f (t_ {n, i}, y_ {n, i}).}
Por lo tanto, para una solución exacta y del problema, tenemos
y(tno,I)=y(tno)+∫tnotno,IF(t,y(t))Dt=y(tno)+hno∫0vsIF(tno+tuhno,y(tno+tuhno))Dtu,∀I=1,...,q,{\ Displaystyle y (t_ {n, i}) = y (t_ {n}) + \ int _ {t_ {n}} ^ {t_ {n, i}} f (t, y (t)) \ mathrm {d} t = y (t_ {n}) + h_ {n} \ int _ {0} ^ {c_ {i}} f (t_ {n} + uh_ {n}, y (t_ {n} + uh_ {n})) \ mathrm {d} u, \, \ forall i = 1, \ dotsc, q,}
y(tno+1)=y(tno)+∫tnotno+1F(t,y(t))Dt=y(tno)+hno∫01F(tno+tuhno,y(tno+tuhno))Dtu.{\ Displaystyle y (t_ {n + 1}) = y (t_ {n}) + \ int _ {t_ {n}} ^ {t_ {n + 1}} f (t, y (t)) \ mathrm {d} t = y (t_ {n}) + h_ {n} \ int _ {0} ^ {1} f (t_ {n} + uh_ {n}, y (t_ {n} + uh_ {n}) )) \ mathrm {d} u.}
Calcularemos estas integrales mediante un método de cuadratura, que podemos elegir para que sea diferente para dos valores distintos de i :
∫0vsIgramo(tu)Dtu≈∑k=1I-1aIkgramo(vsk),∫01gramo(tu)Dtu≈∑k=1qBkgramo(vsk),{\ Displaystyle \ int _ {0} ^ {c_ {i}} g (u) \ mathrm {d} u \ approx \ sum _ {k = 1} ^ {i-1} a_ {ik} g (c_ { k}), \, \ int _ {0} ^ {1} g (u) \ mathrm {d} u \ approx \ sum _ {k = 1} ^ {q} b_ {k} g (c_ {k} ),}
calculado aquí para g ( u ) = f ( t n + uh n , y ( t n + uh n )) .
Por tanto, el método de Runge-Kutta de orden q vendrá dado por:
∀I=1,...,q,{tno,I=tno+vsIhno,yno,I=yno+hno∑k=1I-1aIkpagno,kpagno,I=F(tno,I,yno,I){\ Displaystyle \ forall i = 1, \ dotsc, q, {\ begin {cases} t_ {n, i} & = t_ {n} + c_ {i} h_ {n}, \\ y_ {n, i} & = y_ {n} + h_ {n} \ sum _ {k = 1} ^ {i-1} a_ {ik} p_ {n, k} \\ p_ {n, i} & = f (t_ {n , i}, y_ {n, i}) \ end {cases}}}
yno+1=yno+hno∑k=1qBkpagno,k.{\ Displaystyle y_ {n + 1} = y_ {n} + h_ {n} \ sum _ {k = 1} ^ {q} b_ {k} p_ {n, k}.}
El método a menudo se resume en la tabla de los diferentes pesos en cuadratura, llamada tabla de Butcher :
|
vs1{\ Displaystyle c_ {1}}
|
|
vs2{\ Displaystyle c_ {2}} |
a21{\ Displaystyle a_ {21}}
|
|
vs3{\ Displaystyle c_ {3}} |
a31{\ Displaystyle a_ {31}} |
a32{\ Displaystyle a_ {32}}
|
|
⋮{\ Displaystyle \ vdots} |
⋮{\ Displaystyle \ vdots} |
|
⋱{\ Displaystyle \ ddots}
|
|
vsq{\ Displaystyle c_ {q}}
|
aq1{\ Displaystyle a_ {q1}}
|
aq2{\ Displaystyle a_ {q2}}
|
⋯{\ Displaystyle \ cdots}
|
aq,q-1{\ Displaystyle a_ {q, q-1}} |
|
|
|
B1{\ Displaystyle b_ {1}} |
B2{\ Displaystyle b_ {2}} |
⋯{\ Displaystyle \ cdots} |
Bq-1{\ Displaystyle b_ {q-1}} |
Bq{\ Displaystyle b_ {q}}
|
Sin embargo, el método es consistente .
∀I=2,...,q,∑k=1I-1aIk=vsI{\ Displaystyle \ forall i = 2, \ dotsc, q, \ sum _ {k = 1} ^ {i-1} a_ {ik} = c_ {i}}
Ejemplos de
El método Runge-Kutta de primer orden (RK1)
Este método es equivalente al método de Euler , un método simple de la solución de ecuaciones diferenciales de 1 st grado.
Considere el siguiente problema:
y′=F(t,y),y(t0)=y0{\ Displaystyle y '= f (t, y), \ quad y (t_ {0}) = y_ {0}}El método RK1 usa la ecuación
yno+1=yno+hF(tno,yno){\ Displaystyle y_ {n + 1} = y_ {n} + hf \ left (t_ {n}, y_ {n} \ right)}
donde h es el paso de iteración. Por tanto, el problema está escrito:
|
0{\ displaystyle 0}
|
0{\ displaystyle 0}
|
|
|
1{\ Displaystyle 1}
|
El método de Runge-Kutta de segundo orden (RK2)
El método RK2 del punto medio es una composición del método de Euler :
yno+1=yno+hF(tno+h2,yno+h2F(tno,yno)){\ Displaystyle y_ {n + 1} = y_ {n} + hf \ left (t_ {n} + {\ frac {h} {2}}, y_ {n} + {\ frac {h} {2}} f \ left (t_ {n}, y_ {n} \ right) \ right)}
donde h es el paso de iteración.
Consiste en estimar la derivada en medio del paso de integración:
yno+12=yno+h2F(tno,yno){\ Displaystyle y_ {n + {\ frac {1} {2}}} = y_ {n} + {\ frac {h} {2}} f \ left (t_ {n}, y_ {n} \ right) }
yno+12′=F(tno+h2,yno+12){\ Displaystyle y '_ {n + {\ frac {1} {2}}} = f \ left (t_ {n} + {\ frac {h} {2}}, y_ {n + {\ frac {1 } {2}}} \ right)}
y para rehacer el paso de integración completa a partir de esta estimación:
yno+1=yno+hyno+12′.{\ Displaystyle y_ {n + 1} = y_ {n} + hy '_ {n + {\ frac {1} {2}}}.}
Este esquema se conoce comúnmente como un esquema predictor-correctivo explícito .
Este es el caso particular de α = 1/2 del método más general:
|
0{\ displaystyle 0} |
0{\ displaystyle 0} |
0{\ displaystyle 0}
|
|
α{\ Displaystyle \ alpha}
|
α{\ Displaystyle \ alpha}
|
0{\ displaystyle 0}
|
|
|
1-12α{\ Displaystyle 1 - {\ tfrac {1} {2 \ alpha}}} |
12α{\ displaystyle {\ tfrac {1} {2 \ alpha}}}
|
Por tanto, se reconoce que el método de cuadratura utilizado para los tiempos intermedios es el del punto medio .
Es un método de orden 2 porque el error es del orden de h 3 .
Otro caso común es el método de Heun, correspondiente al caso α = 1 . El método de cuadratura se basa en el método trapezoidal .
Método clásico de cuarto orden Runge-Kutta (RK4)
Es un caso particular de uso muy frecuente, señaló RK4.
Considere el siguiente problema:
y′=F(t,y),y(t0)=y0{\ Displaystyle y '= f (t, y), \ quad y (t_ {0}) = y_ {0}}El método RK4 viene dado por la ecuación:
yno+1=yno+h6(k1+2k2+2k3+k4){\ Displaystyle y_ {n + 1} = y_ {n} + {\ frac {h} {6}} \ left (k_ {1} + 2k_ {2} + 2k_ {3} + k_ {4} \ right) }
o
k1=F(tno,yno){\ Displaystyle k_ {1} = f \ left (t_ {n}, y_ {n} \ right)}
k2=F(tno+h2,yno+h2k1){\ Displaystyle k_ {2} = f \ left (t_ {n} + {\ frac {h} {2}}, y_ {n} + {\ frac {h} {2}} k_ {1} \ right) }
k3=F(tno+h2,yno+h2k2){\ Displaystyle k_ {3} = f \ left (t_ {n} + {\ frac {h} {2}}, y_ {n} + {\ frac {h} {2}} k_ {2} \ right) }
k4=F(tno+h,yno+hk3){\ Displaystyle k_ {4} = f \ left (t_ {n} + h, y_ {n} + hk_ {3} \ right)}
La idea es que el siguiente valor ( y n +1 ) se aproxima por la suma del valor actual ( y n ) y el producto del tamaño del intervalo ( h ) por la pendiente estimada. La pendiente se obtiene mediante un promedio ponderado de pendientes:
-
k 1 es la pendiente al inicio del intervalo;
-
k 2 es la pendiente en el medio del rango, usando la pendiente k 1 para calcular el valor de y en el punto t n + h / 2 mediante el método de Euler ;
-
k 3 es nuevamente la pendiente en el medio del intervalo, pero esta vez se obtiene usando la pendiente k 2 para calcular y ;
-
k 4 es la pendiente al final del intervalo, con el valor de y calculado usando k 3 .
En la media de las cuatro pistas, se da un mayor peso a las pendientes en el punto medio.
|
0{\ displaystyle 0} |
0{\ displaystyle 0} |
0{\ displaystyle 0} |
0{\ displaystyle 0} |
0{\ displaystyle 0}
|
|
12{\ displaystyle {\ tfrac {1} {2}}} |
12{\ displaystyle {\ tfrac {1} {2}}} |
0{\ displaystyle 0} |
0{\ displaystyle 0} |
0{\ displaystyle 0}
|
|
12{\ displaystyle {\ tfrac {1} {2}}} |
0{\ displaystyle 0} |
12{\ displaystyle {\ tfrac {1} {2}}} |
0{\ displaystyle 0} |
0{\ displaystyle 0}
|
|
1{\ Displaystyle 1}
|
0{\ displaystyle 0}
|
0{\ displaystyle 0}
|
1{\ Displaystyle 1}
|
0{\ displaystyle 0}
|
|
|
16{\ displaystyle {\ tfrac {1} {6}}} |
13{\ Displaystyle {\ tfrac {1} {3}}} |
13{\ Displaystyle {\ tfrac {1} {3}}} |
16{\ displaystyle {\ tfrac {1} {6}}}
|
El método RK4 es un método de orden 4, lo que significa que el error cometido en cada paso es del orden de h 5 , mientras que el error total acumulado es del orden de h 4 .
Estas fórmulas también son válidas para funciones con valores vectoriales.
Método de Runge-Kutta de cuarto orden con segunda derivada
En el caso , podemos descomponer el problema en un sistema de ecuaciones:
y″=F(t,y,y′),y(t0)=y0,y′(t0)=y0′{\ Displaystyle y '' = f (t, y, y '), \ quad y (t_ {0}) = y_ {0}, \ quad y' (t_ {0}) = y '_ {0}}
{y′=gramo(t,y,z)z′=s(t,y,z){\ Displaystyle \ left \ {{\ begin {matrix} y '= g (t, y, z) \\ z' = s (t, y, z) \ end {matrix}} \ right.}
con aquí y ' = g ( t , y , z ) = z y s = f .
Aplicando el método RK4 a cada una de estas ecuaciones, simplificando obtenemos:
k1=F(tno,yno,yno′){\ Displaystyle k_ {1} = f \ left (t_ {n}, y_ {n}, y '_ {n} \ right)}
k2=F(tno+h2,yno+h2yno′,yno′+h2k1){\ Displaystyle k_ {2} = f \ left (t_ {n} + {\ frac {h} {2}}, y_ {n} + {\ frac {h} {2}} y '_ {n}, y '_ {n} + {\ frac {h} {2}} k_ {1} \ right)}
k3=F(tno+h2,yno+h2yno′+h24k1,yno′+h2k2){\ Displaystyle k_ {3} = f \ left (t_ {n} + {\ frac {h} {2}}, y_ {n} + {\ frac {h} {2}} y '_ {n} + {\ frac {h ^ {2}} {4}} k_ {1}, y '_ {n} + {\ frac {h} {2}} k_ {2} \ right)}
k4=F(tno+h,yno+hyno′+h22k2,yno′+hk3){\ Displaystyle k_ {4} = f \ left (t_ {n} + h, y_ {n} + hy '_ {n} + {\ frac {h ^ {2}} {2}} k_ {2}, y '_ {n} + hk_ {3} \ right)}
Deducimos y n + 1 y y ' n + 1 gracias a:
yno+1=yno+hyno′+h26(k1+k2+k3){\ Displaystyle y_ {n + 1} = y_ {n} + hy '_ {n} + {\ frac {h ^ {2}} {6}} \ left (k_ {1} + k_ {2} + k_ {3} \ right)}
yno+1′=yno′+h6(k1+2k2+2k3+k4){\ Displaystyle y '_ {n + 1} = y' _ {n} + {\ frac {h} {6}} \ left (k_ {1} + 2k_ {2} + 2k_ {3} + k_ {4 } \ derecho)}
Estabilidad del método
Como métodos de un solo paso y autoexplicativos, uno tiene que cuestionar su estabilidad. Podemos mostrar el siguiente resultado:
Suponga que f es k -lipschitziano en y . Deje
Los métodos de Runge-Kutta son estables en [0, T ] , con la constante de estabilidad e Λ T , con
α=maxI=1,...,q∑k=1I|aIk|.{\ Displaystyle \ alpha = \ max _ {i = 1, \ dotsc, q} \ sum _ {k = 1} ^ {i} | a_ {ik} |.}
Λ=k∑I=1q|BI|(1+αkhmax+...+(αkhmax)q-1).{\ Displaystyle \ Lambda = k \ sum _ {i = 1} ^ {q} | b_ {i} | \ left (1+ \ alpha kh _ {\ max} + \ ldots + (\ alpha kh _ {\ max }) ^ {q-1} \ derecha).}
Referencias
-
[1] Cálculos detallados para RK4 con segunda derivada
-
Jean-Pierre Demailly , Análisis numérico y ecuaciones diferenciales [ detalle de ediciones ]
Artículo relacionado
<img src="https://fr.wikipedia.org/wiki/Special:CentralAutoLogin/start?type=1x1" alt="" title="" width="1" height="1" style="border: none; position: absolute;">