Método Runge-Kutta

La versión actual de la página aún no ha sido revisada por colaboradores experimentados y puede diferir significativamente de la versión revisada el 12 de mayo de 2021; las comprobaciones requieren 2 ediciones .

Los métodos de Runge-Kutta (hay un nombre de métodos de Runge-Kutta en la literatura ) es una gran clase de métodos numéricos para resolver el problema de Cauchy para ecuaciones diferenciales ordinarias y sus sistemas. Los primeros métodos de esta clase fueron propuestos alrededor de 1900 por los matemáticos alemanes K. Runge y M. V. Kutta .

La clase de métodos de Runge-Kutta incluye el método de Euler explícito y el método de Euler modificado con recálculo, que son métodos de primer y segundo orden de precisión, respectivamente. Existen métodos explícitos estándar de tercer orden de precisión que no se utilizan mucho. El más utilizado e implementado en varios paquetes matemáticos ( Maple , MathCAD , Maxima ) es el método clásico de Runge-Kutta , que tiene el cuarto orden de precisión. Cuando se realizan cálculos con mayor precisión, se utilizan cada vez más métodos de quinto y sexto orden de precisión [1] [2] . La construcción de circuitos de orden superior está asociada a grandes dificultades computacionales [3] .

Los métodos de séptimo orden deben tener al menos nueve etapas y los métodos de octavo orden deben tener al menos 11 etapas. Para los métodos del noveno y superiores órdenes (que, sin embargo, no tienen gran importancia práctica), no se sabe cuántas etapas son necesarias para lograr el orden de precisión correspondiente [3] .

El método clásico de Runge-Kutta de cuarto orden

El método de Runge-Kutta de cuarto orden para cálculos con un paso de integración constante está tan extendido que a menudo se lo denomina simplemente método de Runge-Kutta.

Considere el problema de Cauchy para un sistema de ecuaciones diferenciales ordinarias de primer orden. (Además , a ).

Luego, el valor aproximado en puntos posteriores se calcula mediante la fórmula iterativa:

El cálculo de un nuevo valor se realiza en cuatro etapas:

donde  está el valor del paso de cuadrícula en .

Este método tiene el cuarto orden de precisión. Esto significa que el error en un paso es de orden y el error total en el intervalo de integración final es de orden .

Métodos explícitos de Runge-Kutta

La familia de métodos explícitos de Runge-Kutta es una generalización tanto del método explícito de Euler como del método clásico de Runge-Kutta de cuarto orden. viene dada por las formulas

donde  está el valor del paso de cuadrícula y el cálculo del nuevo valor se lleva a cabo en los siguientes pasos:

El método específico está determinado por el número y los coeficientes y . Estos coeficientes a menudo se ordenan en una tabla (llamada tabla Butcher ):

Para los coeficientes del método de Runge-Kutta, se deben cumplir las condiciones para . Si desea que el método tenga orden , también debe proporcionar la condición

donde  es la aproximación obtenida por el método de Runge-Kutta. Después de la diferenciación múltiple, esta condición se transforma en un sistema de ecuaciones polinómicas con respecto a los coeficientes del método.

Métodos implícitos de Runge-Kutta

Todos los métodos de Runge-Kutta mencionados hasta ahora son métodos explícitos . Desafortunadamente, los métodos explícitos de Runge-Kutta, por regla general, no son adecuados para resolver ecuaciones rígidas debido a la pequeña región de su estabilidad absoluta [4] . La inestabilidad de los métodos explícitos de Runge-Kutta también crea problemas muy serios en la solución numérica de ecuaciones diferenciales parciales .

La inestabilidad de los métodos explícitos de Runge-Kutta motivó el desarrollo de métodos implícitos. El método implícito de Runge-Kutta tiene la forma [5] [6] :

dónde

El método explícito se caracteriza por el hecho de que la matriz de coeficientes para él tiene una forma triangular inferior (incluida la diagonal principal cero), en contraste con el método implícito, donde la matriz tiene una forma arbitraria. Esto también se ve en Butcher's table .

Una consecuencia de esta diferencia es la necesidad de resolver el sistema de ecuaciones para , donde  es el número de etapas, en cada paso. Esto aumenta el costo computacional, sin embargo, con un valor suficientemente pequeño, uno puede aplicar el principio de mapeos de contracción y resolver este sistema por iteración simple [7] . En el caso de una iteración, esto aumenta el costo computacional por solo un factor de dos.

Por otro lado, Jean Kunzman (1961) y John Butcher (1964) demostraron que para cualquier número de etapas existe un método implícito de Runge-Kutta con orden de precisión . Esto significa, por ejemplo, que para el método explícito de cuatro etapas de cuarto orden descrito anteriormente, hay un análogo implícito con el doble de precisión.

Método implícito de Runge-Kutta de segundo orden

El método de Runge-Kutta implícito más simple es el método "con recálculo" de Euler modificado. Está dada por la fórmula:

.

Para implementarlo en cada paso, se requieren al menos dos iteraciones (y dos evaluaciones de funciones).

Pronóstico:

.

Corrección:

.

La segunda fórmula es una iteración simple de la solución del sistema de ecuaciones con respecto a , escrita en forma de mapeo de contracción. Para mejorar la precisión, la iteración-corrección se puede hacer varias veces sustituyendo . El método de Euler modificado "con recálculo" tiene el segundo orden de precisión.

Sostenibilidad

La ventaja de los métodos implícitos de Runge-Kutta en comparación con los explícitos es su mayor estabilidad, lo que es especialmente importante cuando se resuelven ecuaciones rígidas . Considere como ejemplo la ecuación lineal y' = λ y . El método habitual de Runge-Kutta aplicado a esta ecuación se reduce a la iteración , con r igual a

donde denota un vector columna de unidades [8] . La función se llama función de estabilidad [9] . Se puede ver en la fórmula que es la razón de dos polinomios de grado si el método tiene etapas. Los métodos explícitos tienen una matriz triangular estrictamente inferior , lo que implica que y que la función de estabilidad es un polinomio [10] .

La solución numérica de este ejemplo converge a cero bajo la condición c . El conjunto de tales se denomina región de estabilidad absoluta . En particular, un método se llama A-estable si todos los c están en la región de estabilidad absoluta. La función de estabilidad del método explícito de Runge-Kutta es un polinomio, por lo que, en principio, los métodos explícitos de Runge-Kutta no pueden ser A-estables [10] .

Si el método tiene orden , entonces la función de estabilidad satisface la condición de . Por lo tanto, la relación de polinomios de un grado dado es de interés, aproximando la función exponencial de la mejor manera. Estas relaciones se conocen como aproximaciones de Padé . La aproximante de Padé con numerador de grado y denominador de grado es A-estable si y solo si [11]

El método Gauss-Legendre de etapas tiene orden , por lo que su función de estabilidad es la aproximación de Padé . Esto implica que el método es A-estable [12] . Esto muestra que los métodos de Runge-Kutta estables en A pueden ser de orden arbitrariamente alto. Por el contrario, el orden de estabilidad A de los métodos de Adams no puede exceder de dos.

Pronunciación

De acuerdo con las normas gramaticales del idioma ruso, el apellido Kutta se declina, por lo que dicen: "El método Runge-Kutta". Las reglas de la gramática rusa prescriben declinar todos los apellidos (incluidos los masculinos) que terminan en -а, -я, precedidos por una consonante. La única excepción son los apellidos de origen francés con acento en la última sílaba como Dumas, Zola [13] . Sin embargo, a veces hay una versión indeclinable del "Método Runge-Kutta" (por ejemplo, en el libro [14] ).

Un ejemplo de una solución en lenguajes de programación algorítmica

haciendo un reemplazo y transfiriendo al lado derecho, obtenemos el sistema:

Código Java para resolver un sistema de ecuaciones diferenciales utilizando el método de Runge-Kutta clase pública MainClass { public static void main ( String [] args ) { int k = 2 ; doble Xo , Yo , Y1 , Zo , Z1 ; doble k1 , k2 , k4 , k3 , h ; doble q1 , q2 , q4 , q3 ; /* *Condiciones iniciales */ Xo = 0 ; Yo = 0,8 ; Zo = 2 ; h = 0,1 ; // paso sistema _ fuera _ println ( "\tX\t\tY\t\tZ" ); para (; r ( Xo , 2 ) < 1.0 ; Xo += h ){ k1 = h * f ( Xo , Yo , Zo ); q1 = h * g ( Xo , Yo , Zo ); k2 = h * f ( Xo + h / 2.0 , Yo + q1 / 2.0 , Zo + k1 / 2.0 ); q2 = h * g ( Xo + h / 2.0 , Yo + q1 / 2.0 , Zo + k1 / 2.0 ); k3 = h * f ( Xo + h / 2.0 , Yo + q2 / 2.0 , Zo + k2 / 2.0 ); q3 = h * g ( Xo + h / 2.0 , Yo + q2 / 2.0 , Zo + k2 / 2.0 ); k4 = h * f ( Xo + h , Yo + q3 , Zo + k3 ); q4 = h * g ( Xo + h , Yo + q3 , Zo + k3 ); Z1 = Zo + ( k1 + 2,0 * k2 + 2,0 * k3 + k4 ) / 6,0 ; Y1 = Yo + ( q1 + 2.0 * q2 + 2.0 * q3 + q4 ) / 6.0 ; sistema _ fuera _ println ( "\t" + r ( Xo + h , k ) + "\t\t" + r ( Y1 , k ) + "\t\t" + r ( Z1 , k )); Yo = Y1 ; Zo = Z1 ; } } /** * función para redondear y descartar la "cola" */ public static double r ( double value , int k ){ return ( double ) Math . ronda (( Math . pow ( 10 , k ) * value )) / Math . pow ( 10 , k ); } /** * funciones que se obtienen del sistema */ public static double f ( double x , double y , double z ){ return ( Math . cos ( 3 * x ) - 4 * y ); } public static double g ( double x , double y , double z ){ return ( z ); } } Código en C# utilizando el sistema ; usando System.Collections.Generic ; espacio de nombres PRJ_RungeKutta { /// <summary> /// Implementación del método Runge-Kutta para una ecuación diferencial ordinaria /// </summary> public abstract class RungeKutta { /// <summary> /// Hora actual /// </resumen> doble t pública ; /// <summary> /// La solución deseada Y[0] es la solución misma, Y[i] es la i-ésima derivada de la solución /// </summary> public double [] Y ; /// <resumen> /// Variables internas /// </resumen> double [] YY , Y1 , Y2 , Y3 , Y4 ; doble protegido [] FY ; /// <summary> /// Constructor /// </summary> /// <param name="N">dimensión del sistema</param> public RungeKutta ( uint N ) { Init ( N ); } /// <summary> /// Constructor /// </summary> public RungeKutta (){} /// <summary> /// Asignar memoria para matrices de trabajo /// </summary> /// <param name="N">Dimensiones de arreglos</param> protected void Init ( uint N ) { Y = new double [ N ]; YY = nuevo doble [ N ]; Y1 = nuevo doble [ N ]; Y2 = nuevo doble [ N ]; Y3 = nuevo doble [ N ]; Y4 = nuevo doble [ N ]; FY = nuevo doble [ N ]; } /// <summary> /// Establecer condiciones iniciales /// </summary> /// <param name="t0">Hora de inicio</param> /// <param name="Y0">Condición inicial </param> public void SetInit ( doble t0 , doble [] Y0 ) { t = t0 ; if ( Y == null ) Init (( uint ) Y0 . Longitud ); for ( int i = 0 ; i < Y . Longitud ; i ++) Y [ i ] = Y0 [ i ]; } /// <summary> /// Cálculo de las partes correctas del sistema /// </summary> /// <param name="t">hora actual</param> /// <param name=" Y">soluciones vectoriales</param> /// <returns>lado derecho</returns> abstract public double [] F ( double t , double [] Y ); /// <summary> /// El siguiente paso del método Runge-Kutta /// </summary> /// <param name="dt">paso de tiempo actual (puede ser variable)</param> public void SiguientePaso ( doble dt ) { int i ; si ( dt < 0 ) devuelve ; // calcular Y1 Y1 = F ( t , Y ); for ( i = 0 ; i < Y . Longitud ; i ++) YY [ i ] = Y [ i ] + Y1 [ i ] * ( dt / 2.0 ); // calcular Y2 Y2 = F ( t + dt / 2.0 , YY ); for ( i = 0 ; i < Y . Longitud ; i ++) YY [ i ] = Y [ i ] + Y2 [ i ] * ( dt / 2.0 ); // calcular Y3 Y3 = F ( t + dt / 2.0 , YY ); for ( i = 0 ; i < Y . Longitud ; i ++) YY [ i ] = Y [ i ] + Y3 [ i ] * dt ; // calcular Y4 Y4 = F ( t + dt , YY ); // calcular la solución en el nuevo paso para ( i = 0 ; i < Y . Longitud ; i ++) Y [ i ] = Y [ i ] + dt / 6.0 * ( Y1 [ i ] + 2.0 * Y2 [ i ] + 2.0 * Y3 [ i ] + Y4 [ i ]); // calcular la hora actual t = t + dt ; } } clase TMyRK : RungeKutta { pública TMyRK ( uint N ) : base ( N ) { } /// <summary> /// ejemplo de péndulo matemático /// y''(t)+y(t)=0 /// </summary> /// <param name="t">Tiempo</param > /// <param name="Y">Solución</param> /// <returns>Lado derecho</returns> public override double [] F ( double t , double [] Y ) { FY [ 0 ] = Y [ 1 ]; AF [ 1 ] = - Y [ 0 ]; volver año fiscal ; } /// <summary> /// Ejemplo de uso /// </summary> static public void Test () { // Paso de tiempo double dt = 0.001 ; // Tarea del objeto del método TMyRK = new TMyRK ( 2 ); // Definir las condiciones iniciales y(0)=0, y'(0)=1 del problema double [] Y0 = { 0 , 1 }; // Establecer las condiciones iniciales para la tarea tarea . EstablecerInicio ( 0 , Y0 ); // resuelve hasta 15 segundos while ( tarea . t <= 15 ) { Console . WriteLine ( "Tiempo = {0:F5}; Func = {1:F8}; d Func / dx = {2:F8}" , tarea . t , tarea . Y [ 0 ], tarea . Y [ 1 ]); // salida t, y, y' // calcular el siguiente paso, tarea de paso de integración . SiguientePaso ( dt ); } Consola . Línea de lectura (); } } Programa de clase { static void Main ( cadena [] argumentos ) { TMyRK . prueba (); } } }

El programa C# usa la clase abstracta RungeKutta , que debe anular el método abstracto F que especifica los lados derechos de las ecuaciones.

Un ejemplo de solución en el entorno MATLAB

La resolución de sistemas de ecuaciones diferenciales por el método de Runge-Kutta es uno de los métodos de solución numérica más comunes en ingeniería. En el entorno MATLAB , se implementa una de sus variedades: el método Dorman-Prince . Para resolver un sistema de ecuaciones, primero debes escribir una función que calcule las derivadas, es decir, las funciones y = g ( x ,  y ,  z ) y z = cos(3 x ) − 4 y = f ( x ,  y ,  z ), como se describe anteriormente . En uno de los directorios a los que se puede acceder desde el sistema MATLAB , debe crear un archivo de texto llamado (por ejemplo) runge.m con el siguiente contenido (para la versión 5.3 de MATLAB):

MATLAB , runge.m función Dy = rango ( x, y ) Dy = y (:); Dy ( 1 ) = y ( 2 ); Dy ( 2 ) = cos ( 3 * x ) - 4 * y ( 1 );

El nombre del archivo y el nombre de la función deben coincidir, pero puede ser cualquier cosa que no se haya usado anteriormente.

Luego, debe crear un archivo principal llamado, por ejemplo, main.m , que hará los cálculos básicos. Este archivo principal contendrá el siguiente texto:

MATLAB , principal.m claro ; clc ; % Borrar memoria y pantalla h = 0,1 ; % Paso de integración aleta_x = 8 ; % Tiempo de integración final y0 = 0,8 ; % Valor inicial de la función Dy0 = 2 ; % Valor inicial de la derivada de la función [ x , y ] = ode45 ( 'runge' , [ 0 : h : x_fin ], [ y0 Dy0 ]); % método de Runge-Kutta parcela ( x , y , 'LineWidth' , 2 ); cuadrícula ; % Parcela y Cuadrícula leyenda ( 'y(x)' , 'y''(x)' , 0 ); % Leyenda en el gráfico

Dado que MATLAB está orientado a matrices, la solución de Runge-Kutta es muy fácil de realizar para un rango completo de x , como, por ejemplo, en el programa de ejemplo anterior. Aquí la solución es la gráfica de la función dentro de los tiempos de 0 a x_fin .

Las variables xey devueltas por la función ODE45 son vectores de valor. Obviamente, la solución al ejemplo anterior es el segundo elemento de x , ya que el primer valor es 0, el paso de integración es h = 0.1 y el valor de interés x = 0.1. La siguiente entrada en la ventana de comandos de MATLAB dará la solución deseada:

soluciónMATLAB _ y1 = y ( encontrar ( x == 0.1 ))

Respuesta: y1 = 0.98768

Notas

  1. Bakhvalov N. S. , Zhidkov N. P. , Kobelkov G. M.  . Métodos numéricos. - M. : Laboratorio de Conocimientos Básicos, 2001. - 630 p. — ISBN 5-93208-043-4 .  - S. 363-375.
  2. Ilyina V. A., Silaev P. K. . Métodos numéricos para físicos teóricos. Parte 2. - Moscú-Izhevsk: Instituto de Investigación Informática, 2004. - 118 p. — ISBN 5-93972-320-9 .  - S. 16-30.
  3. 12 Carnicero J. C.  . Métodos Numéricos para Ecuaciones Diferenciales Ordinarias. — Nueva York: John Wiley & Sons , 2008. — ISBN 978-0-470-72335-7 .
  4. Süli y Mayers, 2003 , pág. 349-351.
  5. Iserlés, 1996 , pág. 41.
  6. Süli y Mayers, 2003 , pág. 351-352.
  7. Métodos implícitos de Runge: Kutta . Archivado el 6 de marzo de 2019 en Wayback Machine .
  8. Hairer & Wanner, 1996 , pág. 40-41.
  9. Hairer & Wanner, 1996 , pág. 40
  10. 1 2 Iserles, 1996 , p. 60
  11. Iserlés, 1996 , pág. 62-63.
  12. Iserlés, 1996 , pág. 63.
  13. Cómo rechazar apellidos (casos difíciles) - "Gramota.ru" - portal de Internet de referencia e información "Idioma ruso" . Consultado el 5 de julio de 2016. Archivado desde el original el 23 de mayo de 2018.
  14. Demidovich B. P., Maron I. A., Shuvalova E. Z. . Métodos numéricos de análisis. 3ra ed. — M .: Nauka , 1967.

Literatura

  • Hairer E., Wanner G. Resolución de ecuaciones diferenciales ordinarias II: Problemas rígidos y algebraicos diferenciales. 2ª ed. - Berlín, Nueva York: Springer-Verlag , 1996. - ISBN 978-3-540-60452-5 .
  • Iserles A. . Un primer curso en el análisis numérico de ecuaciones diferenciales. - Cambridge: Cambridge University Press , 1996. - ISBN 978-0-521-55655-2 .
  • Suli E., Mayers D. . Una introducción al análisis numérico. -Cambridge: Cambridge University Press , 2003. -ISBN 0-521-00794-1 .