四阶龙格库塔（Runge-Kutta）求解微分方程-多种编程语言

Posted 2022-12-03 studyer_domi

前期是分享了matlab下面实现四阶龙格库塔（Runge-Kutta）求解微分方程，这期分享一下C++、C、Java、Python下面的四阶龙格库塔（Runge-Kutta）求解微分方程。

前文传送门：matlab代码实现四阶龙格库塔求解微分方程

C++方法

#include <bits/stdc++.h>using namespace std;
//  "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 Codeint 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;

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 methodint 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 methodimport 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));  

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 methodx0 = 0y = 1x = 2h = 0.2print ('The value of y at x is:', rungeKutta(x0, y, x, h))

The value of y at x is: 1.1036393232374955