1階線形常微分方程式は、工学や物理学で最も頻繁に遭遇する微分方程式の一つです。RC回路の過渡応答、放射性崩壊、ニュートンの冷却法則など、多くの現象がこの形式で記述されます。
本記事では、1階線形常微分方程式を積分因子法(integrating factor method)によって体系的に解く方法を解説します。この方法は、どのような係数を持つ1階線形ODEにも適用できる万能な手法です。
本記事の内容
- 1階線形常微分方程式の標準形
- 積分因子 $\mu(x)$ の導出
- 一般解の公式
- 具体的な計算例
- Pythonでの数値解法と解析解の比較
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
1階線形常微分方程式の標準形
1階線形常微分方程式は、次の標準形に書けます。
$$ \frac{dy}{dx} + P(x)y = Q(x) $$
ここで $P(x)$ と $Q(x)$ は $x$ の既知関数です。
- $Q(x) = 0$ のとき斉次(homogeneous)
- $Q(x) \neq 0$ のとき非斉次(nonhomogeneous)
斉次方程式の解法
まず、$Q(x) = 0$ の場合を考えます。
$$ \frac{dy}{dx} + P(x)y = 0 $$
これは変数分離形に帰着します。
$$ \begin{align} \frac{dy}{y} &= -P(x) \, dx \\ \int \frac{dy}{y} &= -\int P(x) \, dx \\ \ln|y| &= -\int P(x) \, dx + C_0 \\ y &= Ce^{-\int P(x) \, dx} \end{align} $$
積分因子法の導出
非斉次方程式 $y’ + P(x)y = Q(x)$ を解くために、積分因子 $\mu(x)$ を導入します。
両辺に $\mu(x)$ を掛けます。
$$ \mu(x) \frac{dy}{dx} + \mu(x) P(x) y = \mu(x) Q(x) $$
ここで、左辺が $\frac{d}{dx}[\mu(x) y]$ と書ければ直接積分できます。積の微分法則より、
$$ \frac{d}{dx}[\mu(x) y] = \mu(x) \frac{dy}{dx} + \mu'(x) y $$
これが先の式と一致するためには、
$$ \mu'(x) = \mu(x) P(x) $$
が必要です。この微分方程式を解くと、
$$ \begin{align} \frac{d\mu}{\mu} &= P(x) \, dx \\ \ln|\mu| &= \int P(x) \, dx \\ \mu(x) &= e^{\int P(x) \, dx} \end{align} $$
これが積分因子です。
一般解の導出
積分因子 $\mu(x) = e^{\int P(x) \, dx}$ を用いると、元の方程式は
$$ \frac{d}{dx}\left[\mu(x) y\right] = \mu(x) Q(x) $$
と書けます。両辺を積分すると、
$$ \begin{align} \mu(x) y &= \int \mu(x) Q(x) \, dx + C \\ y &= \frac{1}{\mu(x)} \left[\int \mu(x) Q(x) \, dx + C\right] \end{align} $$
これが1階線形ODEの一般解の公式です。
$$ \boxed{y = e^{-\int P \, dx} \left[\int e^{\int P \, dx} Q(x) \, dx + C\right]} $$
具体例
例1: $y’ + 2y = 6e^x$
Step 1: $P(x) = 2$, $Q(x) = 6e^x$ を特定
Step 2: 積分因子を計算
$$ \mu(x) = e^{\int 2 \, dx} = e^{2x} $$
Step 3: 両辺に $\mu(x)$ を掛ける
$$ \begin{align} e^{2x} y’ + 2e^{2x} y &= 6e^{3x} \\ \frac{d}{dx}\left[e^{2x} y\right] &= 6e^{3x} \end{align} $$
Step 4: 両辺を積分
$$ \begin{align} e^{2x} y &= \int 6e^{3x} \, dx \\ &= 2e^{3x} + C \end{align} $$
Step 5: $y$ について解く
$$ y = 2e^x + Ce^{-2x} $$
例2: $y’ + \frac{y}{x} = x^2$, $y(1) = 2$
$P(x) = 1/x$, $Q(x) = x^2$ です。
$$ \mu(x) = e^{\int \frac{1}{x} dx} = e^{\ln x} = x $$
$$ \begin{align} \frac{d}{dx}[xy] &= x^3 \\ xy &= \frac{x^4}{4} + C \\ y &= \frac{x^3}{4} + \frac{C}{x} \end{align} $$
初期条件 $y(1) = 2$ より、
$$ 2 = \frac{1}{4} + C \implies C = \frac{7}{4} $$
$$ y = \frac{x^3}{4} + \frac{7}{4x} $$
Pythonでの実装
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
# 例1: y' + 2y = 6e^x, y(0) = 3
def ode1(x, y):
return -2 * y + 6 * np.exp(x)
sol1 = solve_ivp(ode1, [0, 3], [3], t_eval=np.linspace(0, 3, 200))
# 解析解: y = 2e^x + Ce^{-2x}, y(0)=3 → C=1 → y = 2e^x + e^{-2x}
x_exact = np.linspace(0, 3, 200)
y_exact1 = 2 * np.exp(x_exact) + np.exp(-2 * x_exact)
# 例2: y' + y/x = x^2, y(1) = 2
def ode2(x, y):
return -y / x + x**2
sol2 = solve_ivp(ode2, [1, 5], [2], t_eval=np.linspace(1, 5, 200))
# 解析解: y = x^3/4 + 7/(4x)
y_exact2 = sol2.t**3 / 4 + 7 / (4 * sol2.t)
# 可視化
fig, axes = plt.subplots(1, 2, figsize=(12, 5))
axes[0].plot(sol1.t, sol1.y[0], 'b-', linewidth=2, label='数値解')
axes[0].plot(x_exact, y_exact1, 'r--', linewidth=1.5, label='解析解')
axes[0].set_xlabel('x')
axes[0].set_ylabel('y')
axes[0].set_title("$y' + 2y = 6e^x, \\; y(0)=3$")
axes[0].legend()
axes[0].grid(True, alpha=0.3)
axes[1].plot(sol2.t, sol2.y[0], 'b-', linewidth=2, label='数値解')
axes[1].plot(sol2.t, y_exact2, 'r--', linewidth=1.5, label='解析解')
axes[1].set_xlabel('x')
axes[1].set_ylabel('y')
axes[1].set_title("$y' + y/x = x^2, \\; y(1)=2$")
axes[1].legend()
axes[1].grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('first_order_linear_ode.png', dpi=150, bbox_inches='tight')
plt.show()
数値解と解析解が完全に一致していることが確認できます。
まとめ
本記事では、1階線形常微分方程式の積分因子法による解法を解説しました。
- 標準形は $y’ + P(x)y = Q(x)$
- 積分因子は $\mu(x) = e^{\int P(x) \, dx}$
- $\mu$ を掛けると左辺が $\frac{d}{dx}[\mu y]$ の形になり、直接積分できる
- 一般解は $y = e^{-\int P \, dx}\left[\int e^{\int P \, dx} Q \, dx + C\right]$
- 初期条件で定数 $C$ を決定する
次のステップとして、以下の記事も参考にしてください。