【C++】微分方程式の計算とグラフ表示
はじめに
微分方程式は物理や経済学、工学などの分野で広く使われています。本記事では、C++を用いて基本的な微分方程式を数値的に解く方法を解説し、グラフで結果を可視化する方法も紹介します。
微分方程式の基本概念
微分方程式は、関数とその導関数の関係を記述する方程式です。最も基本的な形は以下のようになります。
\[ \frac{dy}{dx} = f(x, y) \]
例えば、指数関数の微分方程式は \( \frac{dy}{dx} = y \) です。
オイラー法
オイラー法は最も単純な数値解法の一つです。次のように近似計算を行います。
\[ y_{n+1} = y_n + h f(x_n, y_n) \]
以下はオイラー法を使ったC++の実装例です。
#include <iostream>
#include <cmath>
using namespace std;
double f(double x, double y) {
return y; // dy/dx = y
}
void euler(double x0, double y0, double h, int n) {
double x = x0, y = y0;
for (int i = 0; i < n; i++) {
y = y + h * f(x, y);
x = x + h;
cout << "x: " << x << ", y: " << y << endl;
}
}
int main() {
euler(0, 1, 0.1, 10);
return 0;
}
ルンゲ・クッタ法
オイラー法よりも精度の高い方法としてルンゲ・クッタ法(特に4次のもの)がよく使われます。
その基本的な形式は次のようになります。
\[ k_1 = h f(x_n, y_n) \]
\[ k_2 = h f(x_n + h/2, y_n + k_1/2) \]
\[ k_3 = h f(x_n + h/2, y_n + k_2/2) \]
\[ k_4 = h f(x_n + h, y_n + k_3) \]
\[ y_{n+1} = y_n + \frac{k_1 + 2k_2 + 2k_3 + k_4}{6} \]
以下はルンゲ・クッタ法を用いたC++の実装例です。
void rungeKutta(double x0, double y0, double h, int n) {
double x = x0, y = y0;
for (int i = 0; i < n; i++) {
double k1 = h * f(x, y);
double k2 = h * f(x + h/2, y + k1/2);
double k3 = h * f(x + h/2, y + k2/2);
double k4 = h * f(x + h, y + k3);
y = y + (k1 + 2*k2 + 2*k3 + k4) / 6;
x = x + h;
cout << "x: " << x << ", y: " << y << endl;
}
}
グラフ表示
数値計算の結果を視覚的に確認するために、Pythonのmatplotlibを使ってC++の出力をプロットする方法を紹介します。
import matplotlib.pyplot as plt
x_values = []
y_values = []
with open("data.txt", "r") as file:
for line in file:
x, y = map(float, line.split())
x_values.append(x)
y_values.append(y)
plt.plot(x_values, y_values, marker="o")
plt.xlabel("x")
plt.ylabel("y")
plt.title("微分方程式の解")
plt.grid()
plt.show()
まとめ
本記事では、C++での微分方程式の基本的な解法としてオイラー法とルンゲ・クッタ法を紹介しました。また、Pythonを使って結果を可視化する方法についても説明しました。これらの手法を応用すれば、より複雑な問題にも対応できるようになります。