روش رانگ کوتا مرتبه چهارم (RK4) — از صفر تا صد (+ دانلود فیلم آموزش گام به گام)

1# Python program to implement Runge Kutta method 
2# A sample differential equation "dy / dx = (x - y)/2" 
3def dydx(x, y): 
4    return ((x - y)/2) 
5  
6# Finds value of y for a given x using step size h 
7# and initial value y0 at x0. 
8def rungeKutta(x0, y0, x, h): 
9    # Count number of iterations using step size or 
10    # step height h 
11    n = (int)((x - x0)/h)  
12    # Iterate for number of iterations 
13    y = y0 
14    for i in range(1, n + 1): 
15        "Apply Runge Kutta Formulas to find next value of y"
16        k1 = h * dydx(x0, y) 
17        k2 = h * dydx(x0 + 0.5 * h, y + 0.5 * k1) 
18        k3 = h * dydx(x0 + 0.5 * h, y + 0.5 * k2) 
19        k4 = h * dydx(x0 + h, y + k3) 
20  
21        # Update next value of y 
22        y = y + (1.0 / 6.0)*(k1 + 2 * k2 + 2 * k3 + k4) 
23  
24        # Update next value of x 
25        x0 = x0 + h 
26    return y 
27  
28# Driver method 
29x0 = 0
30y = 1
31x = 2
32h = 0.2
33print 'The value of y at x is:', rungeKutta(x0, y, x, h) 
34  
35# This code is contributed by Prateek Bhindwar

1// Java program to implement Runge Kutta method 
2import java.io.*; 
3class differential 
4{ 
5    double dydx(double x, double y) 
6    { 
7        return ((x - y) / 2); 
8    } 
9      
10    // Finds value of y for a given x using step size h 
11    // and initial value y0 at x0. 
12    double rungeKutta(double x0, double y0, double x, double h) 
13    { 
14        differential d1 = new differential(); 
15        // Count number of iterations using step size or 
16        // step height h 
17        int n = (int)((x - x0) / h); 
18  
19        double k1, k2, k3, k4, k5; 
20  
21        // Iterate for number of iterations 
22        double y = y0; 
23        for (int i = 1; i <= n; i++)  
24        { 
25            // Apply Runge Kutta Formulas to find 
26            // next value of y 
27            k1 = h * (d1.dydx(x0, y)); 
28            k2 = h * (d1.dydx(x0 + 0.5 * h, y + 0.5 * k1)); 
29            k3 = h * (d1.dydx(x0 + 0.5 * h, y + 0.5 * k2)); 
30            k4 = h * (d1.dydx(x0 + h, y + k3)); 
31  
32            // Update next value of y 
33            y = y + (1.0 / 6.0) * (k1 + 2 * k2 + 2 * k3 + k4); 
34              
35            // Update next value of x 
36            x0 = x0 + h; 
37        } 
38        return y; 
39    } 
40      
41    public static void main(String args[]) 
42    { 
43        differential d2 = new differential(); 
44        double x0 = 0, y = 1, x = 2, h = 0.2; 
45          
46        System.out.println("\nThe value of y at x is : " 
47                    + d2.rungeKutta(x0, y, x, h)); 
48    } 
49} 
50  
51// This code is contributed by Prateek Bhindwar

1// C program to implement Runge Kutta method 
2#include<stdio.h> 
3  
4// A sample differential equation "dy/dx = (x - y)/2" 
5float dydx(float x, float y) 
6{ 
7    return((x - y)/2); 
8} 
9  
10// Finds value of y for a given x using step size h 
11// and initial value y0 at x0. 
12float rungeKutta(float x0, float y0, float x, float h) 
13{ 
14    // Count number of iterations using step size or 
15    // step height h 
16    int n = (int)((x - x0) / h); 
17  
18    float k1, k2, k3, k4, k5; 
19  
20    // Iterate for number of iterations 
21    float y = y0; 
22    for (int i=1; i<=n; i++) 
23    { 
24        // Apply Runge Kutta Formulas to find 
25        // next value of y 
26        k1 = h*dydx(x0, y); 
27        k2 = h*dydx(x0 + 0.5*h, y + 0.5*k1); 
28        k3 = h*dydx(x0 + 0.5*h, y + 0.5*k2); 
29        k4 = h*dydx(x0 + h, y + k3); 
30  
31        // Update next value of y 
32        y = y + (1.0/6.0)*(k1 + 2*k2 + 2*k3 + k4);; 
33  
34        // Update next value of x 
35        x0 = x0 + h; 
36    } 
37  
38    return y; 
39} 
40  
41// Driver method 
42int main() 
43{ 
44    float x0 = 0, y = 1, x = 2, h = 0.2; 
45    printf("\nThe value of y at x is : %f", 
46            rungeKutta(x0, y, x, h)); 
47    return 0; 
48}

1<?php 
2// PHP program to implement  
3// Runge Kutta method 
4  
5// A sample differential equation 
6// "dy/dx = (x - y)/2" 
7function dydx($x, $y) 
8{ 
9    return(($x - $y) / 2); 
10} 
11  
12// Finds value of y for a 
13// given x using step size h 
14// and initial value y0 at x0. 
15function rungeKutta($x0, $y0, $x, $h) 
16{ 
17      
18    // Count number of iterations 
19    // using step size or step  
20    // height h 
21    $n = (($x - $x0) / $h); 
22  
23    $k1; $k2; $k3; $k4; $k5; 
24  
25    // Iterate for number  
26    // of iterations 
27    $y = $y0; 
28    for($i = 1; $i <= $n; $i++) 
29    { 
30          
31        // Apply Runge Kutta 
32        // Formulas to find 
33        // next value of y 
34        $k1 = $h * dydx($x0, $y); 
35        $k2 = $h * dydx($x0 + 0.5 * $h, 
36                        $y + 0.5 * $k1); 
37        $k3 = $h * dydx($x0 + 0.5 * $h,  
38                        $y + 0.5 * $k2); 
39        $k4 = $h * dydx($x0 + $h, $y + $k3); 
40  
41        // Update next value of y 
42        $y = $y + (1.0 / 6.0) * ($k1 + 2 *  
43                    $k2 + 2 * $k3 + $k4);; 
44  
45        // Update next value of x 
46        $x0 = $x0 + $h; 
47    } 
48  
49    return $y; 
50} 
51  
52    // Driver method 
53    $x0 = 0; 
54    $y = 1;  
55    $x = 2;  
56    $h = 0.2; 
57    echo "The value of y at x is : ", 
58          rungeKutta($x0, $y, $x, $h); 
59  
60// This code is contributed by anuj_67. 
61?>

1// C# program to implement Runge 
2// Kutta method 
3using System; 
4  
5class GFG { 
6      
7    static double dydx(double x, double y) 
8    { 
9        return ((x - y) / 2); 
10    } 
11      
12    // Finds value of y for a given x 
13    // using step size h and initial 
14    // value y0 at x0. 
15    static double rungeKutta(double x0, 
16                double y0, double x, double h) 
17    { 
18      
19        // Count number of iterations using 
20        // step size or step height h 
21        int n = (int)((x - x0) / h); 
22  
23        double k1, k2, k3, k4; 
24  
25        // Iterate for number of iterations 
26        double y = y0; 
27          
28        for (int i = 1; i <= n; i++)  
29        { 
30              
31            // Apply Runge Kutta Formulas 
32            // to find next value of y 
33            k1 = h * (dydx(x0, y)); 
34              
35            k2 = h * (dydx(x0 + 0.5 * h, 
36                             y + 0.5 * k1)); 
37                               
38            k3 = h * (dydx(x0 + 0.5 * h, 
39                            y + 0.5 * k2)); 
40                              
41            k4 = h * (dydx(x0 + h, y + k3)); 
42  
43            // Update next value of y 
44            y = y + (1.0 / 6.0) * (k1 + 2  
45                       * k2 + 2 * k3 + k4); 
46              
47            // Update next value of x 
48            x0 = x0 + h; 
49        } 
50          
51        return y; 
52    } 
53      
54    // Driver code 
55    public static void Main() 
56    { 
57          
58        double x0 = 0, y = 1, x = 2, h = 0.2; 
59          
60        Console.WriteLine("\nThe value of y"
61                             + " at x is : "
62                 + rungeKutta(x0, y, x, h)); 
63    } 
64} 
65  
66// This code is contributed by Sam007.

1clc;
2clear;
3close all;
4%%
5h = 0.1 ; 
6t = 0:h:2 ; 
7x = t ; 
8Nt = numel(t) ; 
9y = zeros(1 , Nt) ; y(1) = 0.5 ; 
10f = @(t , y) y - t^2 +1 ;
11
12for i = 1:Nt-1
13    K1 = h * f(x(i) , y(i)) ; 
14    K2 = h * f(x(i)+h/2 , y(i)+K1/2) ; 
15    K3 = h * f(x(i)+h/2 , y(i)+K2/2) ; 
16    K4 = h * f(x(i)+h , y(i)+K3) ; 
17    y(i+1) = y(i) + 1/6*(K1 + 2*K2 + 2*K3+K4) ;
18end
19
20yr = (1+t).^2-0.5.*exp(t); 
21
22
23figure(1) ; 
24subplot(211);
25plot(t , yr , t , y , 'LineWidth' , 3) ;
26xlabel('t'); ylabel('y') ; 
27grid on
28legend('Real Y' , 'App Y');
29
30subplot(212);
31plot(t , yr-y , 'LineWidth' , 3) ;
32xlabel('t'); ylabel('error') ; 
33grid on
34legend('Error');

1clc;
2clear;
3close all;
4%%
5h = 0.1; 
6x = 0:h:5 ; 
7Nt = numel(x) ; 
8y = zeros(2 , Nt) ; y(: , 1) = [0.5;0] ; 
9f = @(x , y) [y(2) ; x*y(1)];
10
11for i = 1:Nt-1
12    K1 = h * f(x(i) , y(: , i)) ; 
13    K2 = h * f(x(i)+h/2 , y(: , i)+K1/2) ; 
14    K3 = h * f(x(i)+h/2 , y(: , i)+K2/2) ; 
15    K4 = h * f(x(i)+h , y(: , i)+K3) ; 
16    y(: , i+1) = y(:  , i) + 1/6*(K1 + 2*K2 + 2*K3+K4) ;
17end
18
19figure(1) ; 
20plot(x , y , 'LineWidth' , 3) ;
21xlabel('x'); ylabel('y') ; 
22grid on
23legend('Y1 = z(x)' , 'Y2 = dot Z(x)');