Método de Runge-Kutta de cuarto orden para resolver ecuaciones diferenciales

Posted on by Rudeus Greyrat

Dadas las siguientes entradas, 

  • Una ecuación diferencial ordinaria que define el valor de dy/dx en la forma x e y.
  • Valor inicial de y, es decir, y(0)

Así nos lo dan a continuación.
\frac{\mathrm{dy} }{\mathrm{d} x} = f(x, y),y(0)= y_o
La tarea es encontrar el valor de la función desconocida y en un punto dado x.
El método de Runge-Kutta encuentra el valor aproximado de y para una x dada. Solo las ecuaciones diferenciales ordinarias de primer orden se pueden resolver utilizando el método de cuarto orden de Runge Kutta.
A continuación se muestra la fórmula utilizada para calcular el siguiente valor y n+1 a partir del valor anterior y n . Los valores de n son 0, 1, 2, 3, ….(x – x0)/h. Aquí h es la altura del escalón y x n+1 = x 0 + h
. Un tamaño de paso más bajo significa más precisión.
K_1 = hf(x_n, y_n)\\ K_2 = hf(x_n+\frac{h}{2}, y_n+\frac{k_1}{2})\\ K_3 = hf(x_n+\frac{h}{2}, y_n+\frac{k_2}{2})\\ K_4 = hf(x_n+h, y_n+k_3)\\ y_{n+1} = y_n + k_1/6 + k_2/3 + k_3/3 + k_4/6 + O(h^{5})

La fórmula básicamente calcula el siguiente valor y n+1 usando el y n actual más el promedio ponderado de cuatro incrementos. 

  • k 1 es el incremento basado en la pendiente al comienzo del intervalo, usando y
  • k 2 es el incremento basado en la pendiente en el punto medio del intervalo, usando y + hk 1/2 .
  • k 3 es nuevamente el incremento basado en la pendiente en el punto medio, usando y + hk 2 /2.
  • k 4 es el incremento basado en la pendiente al final del intervalo, usando y + hk 3 .

El método es un método de cuarto orden, lo que significa que el error de truncamiento local es del orden de O(h 5 ), mientras que el error total acumulado es del orden de O(h 4 ).
Fuente: https://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods

A continuación se muestra la implementación de la fórmula anterior. 

C++

// C++ program of th above approach
#include <bits/stdc++.h>
using namespace std;
 
// A sample differential equation "dy/dx = (x - y)/2"
float dydx(float x, float y)
{
    return((x - y)/2);
}
  
// Finds value of y for a given x using step size h
// and initial value y0 at x0.
float rungeKutta(float x0, float y0, float x, float h)
{
    // Count number of iterations using step size or
    // step height h
    int n = (int)((x - x0) / h);
  
    float k1, k2, k3, k4, k5;
  
    // Iterate for number of iterations
    float y = y0;
    for (int i=1; i<=n; i++)
    {
        // Apply Runge Kutta Formulas to find
        // next value of y
        k1 = h*dydx(x0, y);
        k2 = h*dydx(x0 + 0.5*h, y + 0.5*k1);
        k3 = h*dydx(x0 + 0.5*h, y + 0.5*k2);
        k4 = h*dydx(x0 + h, y + k3);
  
        // Update next value of y
        y = y + (1.0/6.0)*(k1 + 2*k2 + 2*k3 + k4);;
  
        // Update next value of x
        x0 = x0 + h;
    }
  
    return y;
}
 
// Driver Code
int main()
{
    float x0 = 0, y = 1, x = 2, h = 0.2;
    cout << "The value of y at x is : " <<
            rungeKutta(x0, y, x, h);
 
    return 0;
}
 
// This code is contributed by code_hunt.

C

// C program to implement Runge Kutta method
#include<stdio.h>
 
// A sample differential equation "dy/dx = (x - y)/2"
float dydx(float x, float y)
{
    return((x - y)/2);
}
 
// Finds value of y for a given x using step size h
// and initial value y0 at x0.
float rungeKutta(float x0, float y0, float x, float h)
{
    // Count number of iterations using step size or
    // step height h
    int n = (int)((x - x0) / h);
 
    float k1, k2, k3, k4, k5;
 
    // Iterate for number of iterations
    float y = y0;
    for (int i=1; i<=n; i++)
    {
        // Apply Runge Kutta Formulas to find
        // next value of y
        k1 = h*dydx(x0, y);
        k2 = h*dydx(x0 + 0.5*h, y + 0.5*k1);
        k3 = h*dydx(x0 + 0.5*h, y + 0.5*k2);
        k4 = h*dydx(x0 + h, y + k3);
 
        // Update next value of y
        y = y + (1.0/6.0)*(k1 + 2*k2 + 2*k3 + k4);;
 
        // Update next value of x
        x0 = x0 + h;
    }
 
    return y;
}
 
// Driver method
int main()
{
    float x0 = 0, y = 1, x = 2, h = 0.2;
    printf("\nThe value of y at x is : %f",
            rungeKutta(x0, y, x, h));
    return 0;
}

Java

// Java program to implement Runge Kutta method
import java.io.*;
class differential
{
    double dydx(double x, double y)
    {
        return ((x - y) / 2);
    }
     
    // Finds value of y for a given x using step size h
    // and initial value y0 at x0.
    double rungeKutta(double x0, double y0, double x, double h)
    {
        differential d1 = new differential();
        // Count number of iterations using step size or
        // step height h
        int n = (int)((x - x0) / h);
 
        double k1, k2, k3, k4, k5;
 
        // Iterate for number of iterations
        double y = y0;
        for (int i = 1; i <= n; i++)
        {
            // Apply Runge Kutta Formulas to find
            // next value of y
            k1 = h * (d1.dydx(x0, y));
            k2 = h * (d1.dydx(x0 + 0.5 * h, y + 0.5 * k1));
            k3 = h * (d1.dydx(x0 + 0.5 * h, y + 0.5 * k2));
            k4 = h * (d1.dydx(x0 + h, y + k3));
 
            // Update next value of y
            y = y + (1.0 / 6.0) * (k1 + 2 * k2 + 2 * k3 + k4);
             
            // Update next value of x
            x0 = x0 + h;
        }
        return y;
    }
     
    public static void main(String args[])
    {
        differential d2 = new differential();
        double x0 = 0, y = 1, x = 2, h = 0.2;
         
        System.out.println("\nThe value of y at x is : "
                    + d2.rungeKutta(x0, y, x, h));
    }
}
 
// This code is contributed by Prateek Bhindwar

Python3

# Python program to implement Runge Kutta method
# A sample differential equation "dy / dx = (x - y)/2"
def dydx(x, y):
    return ((x - y)/2)
 
# Finds value of y for a given x using step size h
# and initial value y0 at x0.
def rungeKutta(x0, y0, x, h):
    # Count number of iterations using step size or
    # step height h
    n = (int)((x - x0)/h)
    # Iterate for number of iterations
    y = y0
    for i in range(1, n + 1):
        "Apply Runge Kutta Formulas to find next value of y"
        k1 = h * dydx(x0, y)
        k2 = h * dydx(x0 + 0.5 * h, y + 0.5 * k1)
        k3 = h * dydx(x0 + 0.5 * h, y + 0.5 * k2)
        k4 = h * dydx(x0 + h, y + k3)
 
        # Update next value of y
        y = y + (1.0 / 6.0)*(k1 + 2 * k2 + 2 * k3 + k4)
 
        # Update next value of x
        x0 = x0 + h
    return y
 
# Driver method
x0 = 0
y = 1
x = 2
h = 0.2
print ('The value of y at x is:', rungeKutta(x0, y, x, h))
 
# This code is contributed by Prateek Bhindwar

C#

// C# program to implement Runge
// Kutta method
using System;
 
class GFG {
     
    static double dydx(double x, double y)
    {
        return ((x - y) / 2);
    }
     
    // Finds value of y for a given x
    // using step size h and initial
    // value y0 at x0.
    static double rungeKutta(double x0,
                double y0, double x, double h)
    {
     
        // Count number of iterations using
        // step size or step height h
        int n = (int)((x - x0) / h);
 
        double k1, k2, k3, k4;
 
        // Iterate for number of iterations
        double y = y0;
         
        for (int i = 1; i <= n; i++)
        {
             
            // Apply Runge Kutta Formulas
            // to find next value of y
            k1 = h * (dydx(x0, y));
             
            k2 = h * (dydx(x0 + 0.5 * h,
                             y + 0.5 * k1));
                              
            k3 = h * (dydx(x0 + 0.5 * h,
                            y + 0.5 * k2));
                             
            k4 = h * (dydx(x0 + h, y + k3));
 
            // Update next value of y
            y = y + (1.0 / 6.0) * (k1 + 2
                       * k2 + 2 * k3 + k4);
             
            // Update next value of x
            x0 = x0 + h;
        }
         
        return y;
    }
     
    // Driver code
    public static void Main()
    {
         
        double x0 = 0, y = 1, x = 2, h = 0.2;
         
        Console.WriteLine("\nThe value of y"
                             + " at x is : "
                 + rungeKutta(x0, y, x, h));
    }
}
 
// This code is contributed by Sam007.

PHP

<?php
// PHP program to implement
// Runge Kutta method
 
// A sample differential equation
// "dy/dx = (x - y)/2"
function dydx($x, $y)
{
    return(($x - $y) / 2);
}
 
// Finds value of y for a
// given x using step size h
// and initial value y0 at x0.
function rungeKutta($x0, $y0, $x, $h)
{
     
    // Count number of iterations
    // using step size or step
    // height h
    $n = (($x - $x0) / $h);
 
    $k1; $k2; $k3; $k4; $k5;
 
    // Iterate for number
    // of iterations
    $y = $y0;
    for($i = 1; $i <= $n; $i++)
    {
         
        // Apply Runge Kutta
        // Formulas to find
        // next value of y
        $k1 = $h * dydx($x0, $y);
        $k2 = $h * dydx($x0 + 0.5 * $h,
                        $y + 0.5 * $k1);
        $k3 = $h * dydx($x0 + 0.5 * $h,
                        $y + 0.5 * $k2);
        $k4 = $h * dydx($x0 + $h, $y + $k3);
 
        // Update next value of y
        $y = $y + (1.0 / 6.0) * ($k1 + 2 *
                    $k2 + 2 * $k3 + $k4);;
 
        // Update next value of x
        $x0 = $x0 + $h;
    }
 
    return $y;
}
 
    // Driver method
    $x0 = 0;
    $y = 1;
    $x = 2;
    $h = 0.2;
    echo "The value of y at x is : ",
          rungeKutta($x0, $y, $x, $h);
 
// This code is contributed by anuj_67.
?>

Javascript

<script>
 
// Javascript program to implement Runge Kutta method
 
// A sample differential equation "dy/dx = (x - y)/2"
function dydx(x, y)
{
    return((x - y) / 2);
}
 
// Finds value of y for a given x using step size h
// and initial value y0 at x0.
function rungeKutta(x0, y0, x, h)
{
     
    // Count number of iterations using
    // step size or step height h
    let n = parseInt((x - x0) / h, 10);
 
    let k1, k2, k3, k4, k5;
 
    // Iterate for number of iterations
    let y = y0;
    for(let i = 1; i <= n; i++)
    {
         
        // Apply Runge Kutta Formulas to find
        // next value of y
        k1 = h * dydx(x0, y);
        k2 = h * dydx(x0 + 0.5 * h, y + 0.5 * k1);
        k3 = h * dydx(x0 + 0.5 * h, y + 0.5 * k2);
        k4 = h * dydx(x0 + h, y + k3);
 
        // Update next value of y
        y = y + (1 / 6) * (k1 + 2 * k2 +
                            2 * k3 + k4);;
 
        // Update next value of x
        x0 = x0 + h;
    }
    return y.toFixed(6);
}
 
// Driver code
let x0 = 0, y = 1, x = 2, h = 0.2;
 
document.write("The value of y at x is : " +
               rungeKutta(x0, y, x, h));   
                
// This code is contributed by divyesh072019
 
</script>

Producción: 

The value of y at x is : 1.103639

La complejidad temporal de la solución anterior es O(n) donde n es (x-x0)/h.
Algunos recursos útiles para ejemplos detallados y más explicaciones. 
http://w3.gazi.edu.tr/~balbasi/mws_gen_ode_txt_runge4th.pdf  
https://www.youtube.com/watch?v=kUcc8vAgoQ0 
 

Este artículo es una contribución de Arpit Agarwal . Si le gusta GeeksforGeeks y le gustaría contribuir, también puede escribir un artículo y enviarlo por correo a review-team@geeksforgeeks.org. Vea su artículo que aparece en la página principal de GeeksforGeeks y ayude a otros Geeks.

Escriba comentarios si encuentra algo incorrecto o si desea compartir más información sobre el tema tratado anteriormente.
 

Publicación traducida automáticamente

Artículo escrito por GeeksforGeeks-1 y traducido por Barcelona Geeks. The original can be accessed here. Licence: CCBY-SA

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *