1階常微分方程式の解法まとめ

自然現象の多くは「変化の速度」で記述されます。人口の増加率は現在の人口に比例する。放射性物質の崩壊率は現在の量に比例する。コンデンサの充電速度は現在の電圧差に比例する。これらは全て1階常微分方程式（first-order ordinary differential equation, ODE）として表されます。

$$ \frac{dy}{dx} = f(x, y) $$

微分方程式を解くとは、「瞬間の変化率の情報から全体の振る舞いを復元する」ことです。身近なアナロジーで言えば、速度計の読みだけから自動車の位置を求めるようなものです。速度は微分（位置の変化率）なので、位置を復元するには積分（微分の逆操作）が必要です。1階ODEの解法は、まさにこの「逆操作」を体系的に行う技術です。

1階ODEは微分方程式論の出発点であり、その解法を体系的に理解することは、より高度な微分方程式（2階以上、偏微分方程式）への足がかりとなります。2階以上のODEも、変数変換によって1階の連立ODEに帰着できるため、1階ODEの理解は全ての微分方程式解法の基礎です。

1階ODEを理解すると、以下のような応用が開けます。

• 物理学: 放射性崩壊、RC回路の過渡応答、ニュートンの冷却法則。半減期や時定数の計算に直接使われます
• 生物学: ロジスティック増殖モデル、薬物動態学。個体数の時間発展や血中薬物濃度の予測に必要です
• 工学: タンクの混合問題、化学反応速度論。プロセス制御や環境工学の基礎計算に使われます
• 経済学: 連続複利、ソロー成長モデル。金融工学や経済成長理論の出発点です

本記事の内容

• 変数分離法
• 同次形微分方程式
• 1階線形微分方程式と積分因子法
• 完全微分方程式
• ベルヌーイの方程式
• Pythonによる解析解と数値解の比較

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

• 微積分の基礎 — 微分と積分の基本
• 偏微分の基礎 — 完全微分方程式に必要

1階ODEの分類

1階ODEの右辺 $f(x, y)$ の形によって、適用できる解法が異なります。

実際の微分方程式を解く際には、まず式がどの型に属するかを見極め、対応する解法を適用します。ある方程式が複数の型に該当する場合もあります（例えば変数分離形はつねに1階線形でもある場合がある）。その場合はより簡単な手法を選ぶのが実用的です。それでは各解法を順に見ていきましょう。

変数分離法

基本的な考え方

$y’ = g(x)h(y)$ の形の方程式は、$x$ だけの関数と $y$ だけの関数に分離できます。

$$ \frac{dy}{h(y)} = g(x)\, dx $$

両辺を積分すると解が得られます。

$$ \int \frac{dy}{h(y)} = \int g(x)\, dx + C $$

これは最も基本的で広く適用できる解法です。多くの微分方程式は、適切な変数変換により変数分離形に帰着できます。

具体例1: 指数的成長・減衰

$$ \frac{dy}{dx} = ky \quad (k \text{ は定数}) $$

変数を分離すると

$$ \frac{dy}{y} = k\, dx $$

両辺を積分すると

$$ \ln|y| = kx + C_1 $$

$$ y = Ce^{kx} \quad (C = \pm e^{C_1}) $$

初期条件 $y(0) = y_0$ を適用すると $C = y_0$ であり、$y = y_0 e^{kx}$ が得られます。$k > 0$ なら指数的成長（人口増加、細菌の増殖など）、$k < 0$ なら指数的減衰（放射性崩壊、ニュートンの冷却法則など）を表します。放射性物質の半減期 $T_{1/2}$ は $y_0/2 = y_0 e^{kT_{1/2}}$ から $T_{1/2} = \frac{\ln 2}{|k|}$ で求められます。

具体例2: ロジスティック方程式

$$ \frac{dy}{dx} = ry\left(1 – \frac{y}{K}\right) $$

この方程式は生態学で個体数のモデリングに使われます。$y$ が小さいときは $y’ \approx ry$（指数成長）ですが、$y$ が環境収容力 $K$ に近づくと $(1 – y/K) \to 0$ となり成長が鈍化します。

変数を分離します。

$$ \frac{dy}{y(1 – y/K)} = r\, dx $$

左辺を部分分数分解します。$\frac{1}{y(1 – y/K)} = \frac{K}{y(K – y)}$ とした上で

$$ \frac{K}{y(K-y)} = \frac{A}{y} + \frac{B}{K-y} $$

$y = 0$ のとき $A = 1$、$y = K$ のとき $B = 1$ なので $\frac{1}{y} + \frac{1}{K – y}$ を $K$ で割ったものです。整理すると

$$ \frac{1}{y} + \frac{1}{K – y} = \frac{K}{y(K-y)} $$

両辺を積分すると

$$ \ln|y| – \ln|K – y| = rx + C_1 $$

$$ \ln\left|\frac{y}{K-y}\right| = rx + C_1 $$

指数をとり、初期条件 $y(0) = y_0$ を適用すると $C_1 = \ln\frac{y_0}{K – y_0}$ となり

$$ y(x) = \frac{K}{1 + \left(\frac{K}{y_0} – 1\right)e^{-rx}} $$

$y_0 < K$ のとき、$y$ は $K$（環境収容力）に向かって S 字型に増加します。増加が最も速いのは $y = K/2$ のときで、そこでの増加率は $rK/4$ です。

変数分離法は直感的でわかりやすいですが、全ての方程式に適用できるわけではありません。右辺が $g(x)h(y)$ の積の形に分離できない場合には使えません。次に、比率 $y/x$ の関数として書ける同次形の解法を見ましょう。

同次形微分方程式

定義と解法

$y’ = F(y/x)$ の形を同次形と呼びます。「同次」とは、右辺 $f(x, y)$ が $f(tx, ty) = f(x, y)$ を満たすこと（0次同次関数）を意味します。直感的には、$x$ と $y$ を同じ比率で拡大しても方程式の形が変わらないということです。

$v = y/x$（$y = vx$）と置換します。$y = vx$ を $x$ で微分すると、積の微分法則により $y’ = v + xv’$ です。元の方程式に代入すると

$$ v + xv’ = F(v) $$

$xv’$ について整理すると $xv’ = F(v) – v$ で、$v$ と $x$ を分離できます。

$$ \frac{dv}{F(v) – v} = \frac{dx}{x} $$

変数分離形に帰着しました。

具体例

$$ y’ = \frac{x + y}{x} = 1 + \frac{y}{x} $$

右辺が $y/x$ の関数なので同次形です。$v = y/x$ とおくと $y’ = v + xv’$ で、元の方程式は

$$ v + xv’ = 1 + v $$

両辺から $v$ を引くと $xv’ = 1$ となります。これは変数分離形で

$$ dv = \frac{dx}{x} $$

両辺を積分すると

$$ v = \ln|x| + C $$

$v = y/x$ を代入して $y$ について解くと

$$ y = x(\ln|x| + C) $$

検算: $y’ = \ln|x| + C + x \cdot \frac{1}{x} = \ln|x| + C + 1$ です。一方 $1 + y/x = 1 + \ln|x| + C$ なので一致します。

同次形は変数分離形に帰着できましたが、$y’ = F(y/x)$ の形でない方程式も多くあります。次に、より広い範囲の線形方程式に適用できる積分因子法を見てみましょう。

1階線形微分方程式と積分因子法

標準形

$$ \begin{equation} \frac{dy}{dx} + P(x)y = Q(x) \end{equation} $$

積分因子の導出

積分因子法のアイデアは、「方程式の両辺に適切な関数を掛けて、左辺を何かの微分の形にする」ことです。積の微分法則 $\frac{d}{dx}(\mu y) = \mu’ y + \mu y’$ を思い出してください。これと $\mu y’ + \mu P y$ を見比べると、$\mu’ = \mu P$ であれば $\mu y’ + \mu P y = \mu y’ + \mu’ y = \frac{d}{dx}(\mu y)$ が成立します。

$\mu’ = P\mu$ は変数分離形の方程式であり、$\frac{d\mu}{\mu} = P\, dx$ と書けます。両辺を積分すると $\ln\mu = \int P\, dx$、すなわち

$$ \mu(x) = e^{\int P(x)\,dx} $$

この $\mu(x)$ を積分因子（integrating factor）と呼びます。

元の方程式 $y’ + Py = Q$ の両辺に $\mu$ を掛けると

$$ \mu y’ + \mu P y = \mu Q $$

左辺は $\frac{d}{dx}(\mu y)$ と書けるので

$$ \frac{d}{dx}(\mu y) = \mu Q $$

両辺を $x$ で積分すると

$$ \mu(x) y = \int \mu(x) Q(x)\, dx + C $$

$y$ について解くと

$$ \begin{equation} y = \frac{1}{\mu(x)}\left[\int \mu(x) Q(x)\, dx + C\right] \end{equation} $$

この公式を暗記するよりも、「$\mu$ を掛けて左辺を微分の形にする→両辺積分」という手順を理解する方が実用的です。

具体例: RC回路

RC直列回路にステップ電圧 $E$ を印加したときの電荷 $q(t)$ は、キルヒホッフの電圧則から

$$ R\frac{dq}{dt} + \frac{q}{C} = E $$

と表されます。$R$ で割って標準形にすると $q’ + \frac{1}{RC}q = \frac{E}{R}$ です。

ここで $P(t) = \frac{1}{RC}$ なので、積分因子は

$$ \mu(t) = e^{\int \frac{1}{RC}\, dt} = e^{t/(RC)} $$

両辺に $\mu$ を掛けると $\frac{d}{dt}(e^{t/(RC)} q) = \frac{E}{R} e^{t/(RC)}$ です。両辺を積分すると

$$ e^{t/(RC)} q = \frac{E}{R} \cdot RC \cdot e^{t/(RC)} + C_1 = CE \cdot e^{t/(RC)} + C_1 $$

$q_0 = q(0)$ として初期条件を適用すると $C_1 = q_0 – CE$ です。整理すると

$$ q(t) = CE(1 – e^{-t/(RC)}) + q_0 e^{-t/(RC)} $$

時定数 $\tau = RC$ が充電の速さを決定します。$t = \tau$ で、初期値 $q_0 = 0$ の場合、$q(\tau) = CE(1 – e^{-1}) \approx 0.632 CE$ となり、最終値の約63.2%まで充電されます。$t = 5\tau$ で $q \approx 0.993 CE$ となり、ほぼ完全に充電されます。

具体例: 混合問題

容量 100L のタンクに初期濃度 0 の溶液が入っています。濃度 0.5 kg/L の溶液が 2 L/min で流入し、よく混ぜられた溶液が 2 L/min で流出します。流入量と流出量が等しいため、タンク内の液体量は一定に保たれます。

タンク内の溶質量 $y(t)$ (kg) の変化率は「流入量 – 流出量」で決まります。

• 流入: $0.5 \text{ kg/L} \times 2 \text{ L/min} = 1 \text{ kg/min}$
• 流出: タンク内の濃度 $\frac{y}{100}$ kg/L の溶液が 2 L/min で流出するので $\frac{y}{100} \times 2 = \frac{y}{50}$ kg/min

$$ \frac{dy}{dt} = 1 – \frac{y}{50} $$

標準形に書き直すと $y’ + \frac{1}{50}y = 1$ です。

積分因子は $\mu = e^{\int \frac{1}{50}dt} = e^{t/50}$ です。$\frac{d}{dt}(e^{t/50}y) = e^{t/50}$ を積分すると $e^{t/50}y = 50e^{t/50} + C$ です。初期条件 $y(0) = 0$ から $C = -50$ で

$$ y(t) = 50(1 – e^{-t/50}) $$

平衡値は $y = 50$ kg（タンク内の濃度が流入濃度 $0.5$ kg/L に等しくなるとき）で、時定数は $\tau = 50$ 分です。

積分因子法はすべての1階線形ODEに適用できる万能の手法です。しかし、世の中には $y’$ に関して線形でない方程式も多く存在します。次に見る完全微分方程式は、線形・非線形を問わず適用できる別のアプローチです。

完全微分方程式

定義

$$ M(x,y)\, dx + N(x,y)\, dy = 0 $$

が完全微分方程式であるとは

$$ \frac{\partial M}{\partial y} = \frac{\partial N}{\partial x} $$

が成り立つことです。

この条件はどこから来るのでしょうか。もしポテンシャル関数 $F(x,y)$ が存在して $dF = M\, dx + N\, dy$ と書けるなら、$F$ の全微分の定義 $dF = \frac{\partial F}{\partial x}dx + \frac{\partial F}{\partial y}dy$ と比較して $M = \frac{\partial F}{\partial x}$, $N = \frac{\partial F}{\partial y}$ です。このとき、2階偏微分の対称性 $\frac{\partial^2 F}{\partial y \partial x} = \frac{\partial^2 F}{\partial x \partial y}$ から

$$ \frac{\partial M}{\partial y} = \frac{\partial^2 F}{\partial y \partial x} = \frac{\partial^2 F}{\partial x \partial y} = \frac{\partial N}{\partial x} $$

が成り立ちます。つまり、$M_y = N_x$ は「ポテンシャル関数が存在するための必要十分条件」です。

ポテンシャル関数 $F(x,y)$ が見つかれば、元の方程式は $dF = 0$、すなわち $F(x,y) = C$ で解が与えられます。

求解手順

1. $F = \int M\, dx + g(y)$（$M$ を $x$ について積分し、$y$ のみの未知関数 $g(y)$ を残す）
2. $\frac{\partial F}{\partial y} = N$ の条件から $g'(y)$ を求める
3. $g(y)$ を積分して $F$ を完成させる

具体例

$$ (2xy + 3)\, dx + (x^2 + 4y)\, dy = 0 $$

完全性の確認: $M = 2xy + 3$, $N = x^2 + 4y$ として

$$ \frac{\partial M}{\partial y} = 2x, \quad \frac{\partial N}{\partial x} = 2x $$

$M_y = N_x = 2x$ なので完全微分方程式です。

ステップ1: $F = \int M\, dx = \int (2xy + 3)\, dx = x^2 y + 3x + g(y)$

ステップ2: $\frac{\partial F}{\partial y} = x^2 + g'(y)$ これが $N = x^2 + 4y$ に等しいので $g'(y) = 4y$

ステップ3: $g(y) = 2y^2$

したがって、解は $F(x,y) = x^2 y + 3x + 2y^2 = C$ です。

検算: $dF = (2xy + 3)dx + (x^2 + 4y)dy = 0$ で元の方程式に一致します。

完全微分方程式はポテンシャル関数の存在を利用する美しい解法ですが、$M_y = N_x$ が成り立たない場合は直接適用できません。最後に、非線形な方程式を線形に変換する巧妙な手法を見ましょう。

ベルヌーイの方程式

定義と解法

$$ y’ + P(x)y = Q(x)y^n \quad (n \neq 0, 1) $$

$n = 0$ なら普通の1階線形方程式、$n = 1$ なら変数分離形なので、$n \neq 0, 1$ の場合が本質的に非線形です。右辺に $y^n$ があるため、積分因子法がそのままでは使えません。

しかし、巧妙な置換 $v = y^{1-n}$ を行うと線形方程式に変換できます。その手順を丁寧に追ってみましょう。

まず、元の方程式の両辺を $y^n$ で割ります。

$$ y^{-n} y’ + P(x) y^{1-n} = Q(x) $$

ここで $v = y^{1-n}$ とおくと、$v$ を $x$ で微分すると連鎖律から

$$ v’ = (1-n) y^{-n} y’ $$

よって $y^{-n} y’ = \frac{v’}{1-n}$ です。これを先ほどの式に代入すると

$$ \frac{v’}{1-n} + P(x) v = Q(x) $$

両辺に $(1-n)$ を掛けて整理すると

$$ v’ + (1-n)P(x)v = (1-n)Q(x) $$

$v$ に関する1階線形方程式に帰着しました。これは積分因子法で解けます。

具体例: $y’ + y = y^2$, $y(0) = 0.5$

$P(x) = 1$, $Q(x) = 1$, $n = 2$ です。$v = y^{1-2} = y^{-1} = 1/y$ とおくと、$v’ + (1-2) \cdot 1 \cdot v = (1-2) \cdot 1$、すなわち

$$ v’ – v = -1 $$

これは積分因子 $\mu = e^{-x}$ の1階線形方程式です。$\frac{d}{dx}(e^{-x}v) = -e^{-x}$ を積分すると $e^{-x}v = e^{-x} + C$、よって $v = 1 + Ce^{x}$ です。

$v(0) = 1/y(0) = 2$ より $C = 1$ で、$v = 1 + e^x$ です。$y = 1/v$ に戻すと

$$ y = \frac{1}{1 + e^x} $$

$x \to \infty$ で $y \to 0$（$e^x$ が支配的）、$x \to -\infty$ で $y \to 1$ となるシグモイド型の曲線です。

Pythonでの実装

各解法の解析解と数値解の比較

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp

fig, axes = plt.subplots(2, 3, figsize=(16, 10))

# (1) 変数分離法: y' = ky (指数成長)
ax = axes[0, 0]
k = 0.5
y0 = 1.0
t = np.linspace(0, 5, 200)
y_exact = y0 * np.exp(k * t)
sol = solve_ivp(lambda t, y: k*y, [0, 5], [y0], t_eval=t)
ax.plot(t, y_exact, "r-", linewidth=2, label="Exact: $y_0 e^{kt}$")
ax.plot(sol.t, sol.y[0], "b--", linewidth=2, label="Numerical (RK45)")
ax.set_title("Exponential Growth: $y' = 0.5y$", fontsize=11)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
ax.set_xlabel("x", fontsize=10)

# (2) ロジスティック方程式
ax = axes[0, 1]
r, K, y0_log = 1.0, 10.0, 1.0
y_logistic = K / (1 + (K/y0_log - 1)*np.exp(-r*t))
sol = solve_ivp(lambda t, y: r*y*(1-y/K), [0, 5], [y0_log], t_eval=t)
ax.plot(t, y_logistic, "r-", linewidth=2, label="Exact")
ax.plot(sol.t, sol.y[0], "b--", linewidth=2, label="Numerical")
ax.axhline(K, color="gray", linestyle=":", alpha=0.5, label=f"K = {K}")
ax.set_title("Logistic: $y' = y(1-y/10)$", fontsize=11)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
ax.set_xlabel("x", fontsize=10)

# (3) 1階線形: y' + y = sin(x)
ax = axes[0, 2]
y_exact_lin = np.exp(-t) * (1 + 0.5) + 0.5*(np.sin(t) - np.cos(t))
# 正確な解析解
y0_lin = 1.0
# y' + y = sin(x), μ = e^x
# y = e^{-x}[∫e^x sin(x) dx + C] = e^{-x}[e^x(sin x - cos x)/2 + C]
# = (sin x - cos x)/2 + Ce^{-x}
# y(0) = (-1)/2 + C = 1 → C = 3/2
y_exact_lin = (np.sin(t) - np.cos(t))/2 + 1.5*np.exp(-t)
sol = solve_ivp(lambda t, y: -y + np.sin(t), [0, 5], [y0_lin], t_eval=t)
ax.plot(t, y_exact_lin, "r-", linewidth=2, label="Exact (integrating factor)")
ax.plot(sol.t, sol.y[0], "b--", linewidth=2, label="Numerical")
ax.set_title("Linear: $y' + y = \\sin(x)$", fontsize=11)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
ax.set_xlabel("x", fontsize=10)

# (4) RC回路
ax = axes[1, 0]
R, C_cap, E = 100, 0.001, 5.0
tau = R * C_cap
t_rc = np.linspace(0, 0.5, 200)
q_exact = C_cap * E * (1 - np.exp(-t_rc/tau))
sol = solve_ivp(lambda t, q: (E - q/C_cap)/R, [0, 0.5], [0], t_eval=t_rc)
ax.plot(t_rc*1000, q_exact*1000, "r-", linewidth=2, label="Exact")
ax.plot(sol.t*1000, sol.y[0]*1000, "b--", linewidth=2, label="Numerical")
ax.set_title(f"RC Circuit: $\\tau$ = {tau*1000:.0f} ms", fontsize=11)
ax.set_xlabel("Time (ms)", fontsize=10)
ax.set_ylabel("Charge (mC)", fontsize=10)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)

# (5) ベルヌーイ: y' + y = y^2
ax = axes[1, 1]
y0_bern = 0.5
# v = 1/y → v' - v = -1 → v = Ce^x + 1
# y(0) = 0.5 → v(0) = 2 → C = 1
# v = e^x + 1 → y = 1/(e^x + 1)
t_bern = np.linspace(0, 4, 200)
y_bern = 1.0 / (np.exp(t_bern) + 1)
sol = solve_ivp(lambda t, y: -y + y**2, [0, 4], [y0_bern], t_eval=t_bern)
ax.plot(t_bern, y_bern, "r-", linewidth=2, label="Exact: $1/(e^x + 1)$")
ax.plot(sol.t, sol.y[0], "b--", linewidth=2, label="Numerical")
ax.set_title("Bernoulli: $y' + y = y^2$", fontsize=11)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
ax.set_xlabel("x", fontsize=10)

# (6) 混合問題
ax = axes[1, 2]
t_tank = np.linspace(0, 200, 200)
y_tank = 50 * (1 - np.exp(-t_tank/50))
sol = solve_ivp(lambda t, y: 1 - y/50, [0, 200], [0], t_eval=t_tank)
ax.plot(t_tank, y_tank, "r-", linewidth=2, label="Exact")
ax.plot(sol.t, sol.y[0], "b--", linewidth=2, label="Numerical")
ax.axhline(50, color="gray", linestyle=":", alpha=0.5, label="Equilibrium = 50 kg")
ax.set_title("Mixing Problem", fontsize=11)
ax.set_xlabel("Time (min)", fontsize=10)
ax.set_ylabel("Solute (kg)", fontsize=10)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig("ode_first_order.png", dpi=150, bbox_inches="tight")
plt.show()

このグラフから、各型の1階ODEの解析解（赤い実線）と数値解（青い破線）が完全に一致していることが確認できます。これは解析解の正しさの検証であると同時に、scipy の solve_ivp（RK45法）の精度の高さも示しています。

1. 指数成長（左上）: $y’ = 0.5y$ の解 $y = e^{0.5x}$ は指数的に増大します。$x = 5$ で $y \approx 12.2$ です。成長率 $k$ が正なので解は発散しますが、$k < 0$ なら指数減衰（放射性崩壊など）になります。変数分離法の最も基本的な例です。

2. ロジスティック（中央上）: S字型の成長曲線が環境収容力 $K = 10$（灰色の点線）に漸近します。初期に近い段階では指数的に増加し、$y = K/2 = 5$ 付近で増加率が最大になり、その後 $K$ に向かって増加が緩やかになります。変数分離法で部分分数分解を使って解きました。

3. 1階線形（右上）: $y’ + y = \sin(x)$ を積分因子法で解いた結果です。解は $y = \frac{\sin x – \cos x}{2} + \frac{3}{2}e^{-x}$ で、過渡項 $\frac{3}{2}e^{-x}$ が減衰した後、定常状態の $\frac{\sin x – \cos x}{2}$ に落ち着きます。入力 $\sin x$ と同じ周波数で振動しますが、位相がずれている点に注目してください。

4. RC回路（左下）: 時定数 $\tau = RC = 100$ ms に支配される指数的な充電曲線です。$t = \tau$ で最終値の約63%、$t = 3\tau$ で約95%に達します。横軸はミリ秒で表示されており、RC回路の過渡応答が数百ミリ秒のオーダーであることがわかります。

5. ベルヌーイ（中央下）: $y’ + y = y^2$ を $v = 1/y$ の置換で線形化して解きました。解 $y = 1/(e^x + 1)$ はシグモイド関数の左右反転した形であり、$y(0) = 0.5$ からスタートして $y \to 0$ に単調減少します。非線形ODEでも適切な変数変換で解析解が得られる好例です。

6. 混合問題（右下）: 溶質量が平衡値 50 kg に指数的に近づきます。時定数 $\tau = 50$ 分で、$t = 150$ 分（$= 3\tau$）で約95%に達しています。タンクの容量と流量の比が時定数を決定するという物理的な対応が見て取れます。

解法選択のガイドライン

5つの解法を学んできましたが、実際の問題に直面したとき「どの解法を使えばよいか」を判断する手順を整理しておきます。

ステップ1: 変数分離形かどうか確認する

右辺が $f(x,y) = g(x)h(y)$ の形に書けるなら、変数分離法を使います。これが最もシンプルで計算が楽です。例えば $y’ = xy^2$ は $g(x) = x$, $h(y) = y^2$ と分離できます。

ステップ2: 1階線形かどうか確認する

$y’ + P(x)y = Q(x)$ の形なら積分因子法を使います。$y$ について1次式であることが重要です。$y^2$ や $\sin(y)$ が含まれていたら線形ではありません。

ステップ3: 同次形かどうか確認する

右辺が $y/x$ の関数として書ける場合（$f(tx, ty) = f(x, y)$ が成り立つ場合）、$v = y/x$ の置換で変数分離形に帰着します。例えば $y’ = \frac{x^2 + y^2}{xy}$ は $\frac{1 + (y/x)^2}{y/x}$ と書けるので同次形です。

ステップ4: 完全微分方程式かどうか確認する

$M\,dx + N\,dy = 0$ の形に書いたとき $M_y = N_x$ が成り立てば完全微分方程式です。ポテンシャル関数を求めます。

ステップ5: ベルヌーイ形かどうか確認する

$y’ + P(x)y = Q(x)y^n$ の形なら、$v = y^{1-n}$ の置換で線形に帰着します。

ステップ6: 上記のいずれにも該当しない場合

変数変換を工夫するか、数値解法（オイラー法、ルンゲ-クッタ法など）を使います。解析的に解けない微分方程式は実は非常に多く、数値解法の重要性はそこにあります。

実際には、1つの方程式が複数の型に該当する場合もあります。例えば $y’ = -2y$ は変数分離形でもあり、1階線形でもあります。そのような場合は、より計算の簡単な方法を選べばよいのです。

まとめ

本記事では、1階常微分方程式の主要な解法を体系的に解説しました。

• 変数分離法: $y’ = g(x)h(y)$ の形で、$x$ と $y$ を分離して両辺を積分する最も基本的な手法。指数成長・減衰やロジスティック方程式に適用
• 同次形: $v = y/x$ の置換で変数分離形に帰着させる。右辺が $y/x$ の関数として書ける場合に適用
• 積分因子法: $y’ + Py = Q$ に $\mu = e^{\int P\,dx}$ を掛けて積の微分法則の形にする。すべての1階線形ODEに適用可能な万能の手法
• 完全微分方程式: $M_y = N_x$ の条件のもとでポテンシャル関数を構成する。2階偏微分の対称性が完全性の条件の由来
• ベルヌーイの方程式: $v = y^{1-n}$ の置換で線形方程式に帰着させる。非線形ODEを線形化する巧妙なテクニック
• 解法選択: まず変数分離→線形→同次→完全→ベルヌーイの順に確認し、該当する型の手法を適用する

次のステップとして、以下の記事も参考にしてください。

• 2階線形常微分方程式の解法 — 振動・減衰系の解析
• 級数解法とべき級数展開による微分方程式の解法 — 特殊関数への入口
• 常微分方程式の数値解法まとめ — 数値計算の手法