定数係数の2階ODEは特性方程式で解けましたが、係数が変数関数の場合はどうでしょうか。たとえばベッセルの方程式 $x^2y” + xy’ + (x^2 – n^2)y = 0$ は、係数に $x$ が入っているため特性方程式のテクニックが使えません。
こうした変数係数のODEに対して威力を発揮するのが級数解法（power series method）です。解を $y = \sum a_n x^n$ とべき級数で仮定し、方程式に代入して係数 $a_n$ の漸化式を導出します。テイラー展開が「既知の関数をべき級数に分解する」作業だとすれば、級数解法は「未知の関数をべき級数で組み立てる」逆方向の作業です。方程式そのものが係数 $a_n$ の間の関係を教えてくれるので、漸化式を順に解いていけば解の全貌が見えてきます。
級数解法を学ぶ最大のモチベーションは、物理学に登場する多くの特殊関数（ルジャンドル多項式、ベッセル関数、エルミート多項式など）の理論的基盤を理解できることです。これらの関数は天から降ってくるものではなく、具体的な物理法則を記述する微分方程式をべき級数で解いた結果として自然に現れるものです。
- 量子力学: 水素原子のシュレーディンガー方程式（ルジャンドル多項式、ラゲール多項式）
- 電磁気学: 円筒座標・球座標でのラプラス方程式（ベッセル関数、球面調和関数）
- 振動工学: ドラムの振動モード（ベッセル関数）
- 統計力学: エルミート多項式（量子調和振動子）
本記事の内容
- べき級数解法の基本手順
- 正則点での級数解
- 確定特異点とフロベニウス法
- ルジャンドルの方程式とルジャンドル多項式
- ベッセルの方程式とベッセル関数
- Pythonによる実装
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
- 2階線形常微分方程式の解法 — 定数係数の場合
- テイラー展開の基礎 — べき級数の収束
べき級数解法の基本手順
級数解法の考え方
級数解法のアイデアを直感的に説明しましょう。関数 $y(x)$ は、$x = 0$ の周りで十分滑らかであれば
$$ y(x) = y(0) + y'(0)x + \frac{y”(0)}{2!}x^2 + \frac{y”'(0)}{3!}x^3 + \cdots $$
とテイラー展開できます。つまり $y(0), y'(0), y”(0), \dots$ さえわかれば関数全体が決まります。微分方程式はまさにこれらの導関数の間の関係式を与えてくれるので、初期値 $y(0)$ と $y'(0)$ から出発して $y”(0), y”'(0), \dots$ を芋づる式に決定できるのです。
級数解法ではテイラー係数を $a_n$ と書いて一般論を展開します。
手順
解を次のべき級数で仮定します（$x = 0$ の周りの展開）。
$$ y = \sum_{n=0}^{\infty} a_n x^n = a_0 + a_1 x + a_2 x^2 + \cdots $$
ここで $a_n = y^{(n)}(0)/n!$ ですが、この値を個別に求めるのではなく、微分方程式に代入して $a_n$ 同士の関係（漸化式）を導出するのがポイントです。
これを微分方程式に代入するため、導関数のべき級数表現も用意しておきます。
$$ y’ = \sum_{n=1}^{\infty} n a_n x^{n-1}, \quad y” = \sum_{n=2}^{\infty} n(n-1) a_n x^{n-2} $$
微分方程式に代入した後、すべての項を $x^k$ の形に整理し、各 $x^k$ の係数が独立にゼロになるという条件から $a_n$ の漸化式を得ます。この「各べきの係数がゼロ」という条件は、べき級数の一意性定理に基づいています。
正則点と特異点
べき級数解が存在するかどうかは、展開中心 $x_0$ が正則点か特異点かで決まります。この分類は級数解法の適用可否を判断するための最初のステップです。
$y” + P(x)y’ + Q(x)y = 0$ に対して、$x = x_0$ が正則点（ordinary point）であるとは、$P(x)$ と $Q(x)$ が $x_0$ で解析的（テイラー展開可能）であることです。直感的には、方程式の係数が $x_0$ で「おとなしい」（発散しない）ということです。正則点ではべき級数解が必ず存在し、その収束半径は $P(x)$ と $Q(x)$ の最も近い特異点までの距離以上であることが保証されます（Fuchsの定理）。
$x_0$ が特異点（singular point）の場合、つまり $P(x)$ または $Q(x)$ が $x_0$ で解析的でない場合は、通常のべき級数解が存在しない可能性があります。しかし特異点にも「程度」があります。
- 確定特異点（regular singular point）: $(x – x_0)P(x)$ と $(x – x_0)^2 Q(x)$ が $x_0$ で解析的である場合です。特異性が「穏やか」で、後述するフロベニウス法で対処できます
- 不確定特異点（irregular singular point）: 上記の条件を満たさない場合です。特異性が「激しく」、級数解法の標準的な手法では扱えません
たとえば $y” + y’/x + y/x^2 = 0$ では、$P(x) = 1/x$、$Q(x) = 1/x^2$ です。$xP(x) = 1$、$x^2Q(x) = 1$ はいずれも $x = 0$ で解析的なので、$x = 0$ は確定特異点です。
一方、$y” + y/x^3 = 0$ では $P(x) = 0$、$Q(x) = 1/x^3$ です。$x^2Q(x) = 1/x$ は $x = 0$ で解析的でないため、$x = 0$ は不確定特異点です。
この分類を理解しておくと、与えられた微分方程式に対してどの解法を適用すべきかがすぐに判断できます。正則点なら通常のべき級数解法、確定特異点ならフロベニウス法を使います。
具体例: エアリの方程式
正則点でのべき級数解法の手順を、エアリの方程式を使って丁寧に確認しましょう。エアリの方程式は量子力学のWKB近似や光学の回折理論に現れる重要な方程式です。
$$ y” – xy = 0 $$
まず $x = 0$ が正則点であることを確認します。標準形 $y” + P(x)y’ + Q(x)y = 0$ と比較すると $P(x) = 0$、$Q(x) = -x$ です。どちらも $x = 0$ で解析的なので、$x = 0$ は正則点です。
$y = \sum_{n=0}^{\infty} a_n x^n$ を代入します。$y” = \sum_{n=2}^{\infty} n(n-1)a_n x^{n-2}$ なので
$$ \sum_{n=2}^{\infty} n(n-1)a_n x^{n-2} – \sum_{n=0}^{\infty} a_n x^{n+1} = 0 $$
$x^k$ の係数をそろえるため、第1項で $m = n – 2$（$n = m + 2$）と置換し、第2項で $m = n + 1$（$n = m – 1$）と置換すると
$$ \sum_{m=0}^{\infty} (m+2)(m+1)a_{m+2} x^m – \sum_{m=1}^{\infty} a_{m-1} x^m = 0 $$
ここで重要なのは、第2項の和が $m = 1$ から始まることです。これにより $m = 0$ の項は第1項にしか存在しません。
$m = 0$ の係数: $(2)(1)a_2 = 0$ なので $a_2 = 0$
$m \geq 1$ の係数: $(m+2)(m+1)a_{m+2} – a_{m-1} = 0$
漸化式は
$$ a_{m+2} = \frac{a_{m-1}}{(m+2)(m+1)}, \quad m \geq 1 $$
$a_0$ と $a_1$ が自由パラメータ（2階ODEなので2つの初期条件に対応）です。$a_2 = 0$ なので、$a_0$ 系列と $a_1$ 系列が独立に計算できます。
$a_0$ 系列（$a_1 = 0$ としたとき）: $a_0, 0, 0, a_3, 0, 0, a_6, \dots$ という3つ飛びのパターンです。
$$ a_3 = \frac{a_0}{3 \cdot 2}, \quad a_6 = \frac{a_3}{6 \cdot 5} = \frac{a_0}{6 \cdot 5 \cdot 3 \cdot 2}, \quad a_9 = \frac{a_6}{9 \cdot 8} = \frac{a_0}{9 \cdot 8 \cdot 6 \cdot 5 \cdot 3 \cdot 2} $$
$a_1$ 系列（$a_0 = 0$ としたとき）: $0, a_1, 0, 0, a_4, 0, 0, a_7, \dots$
$$ a_4 = \frac{a_1}{4 \cdot 3}, \quad a_7 = \frac{a_4}{7 \cdot 6} = \frac{a_1}{7 \cdot 6 \cdot 4 \cdot 3}, \quad a_{10} = \frac{a_7}{10 \cdot 9} = \frac{a_1}{10 \cdot 9 \cdot 7 \cdot 6 \cdot 4 \cdot 3} $$
$a_0$ 系列からエアリ関数 $\text{Ai}(x)$ が、$a_1$ 系列からエアリ関数 $\text{Bi}(x)$ が構成されます（適切な正規化を施した後）。2つの自由パラメータから2つの独立な解が得られており、2階ODEの一般解の構造と整合しています。
ここまでは正則点での級数解法を見てきました。次に、特異点を持つ微分方程式に対応するフロベニウス法を学びましょう。
フロベニウス法
フロベニウス法の直感的理解
正則点でのべき級数 $y = \sum a_n x^n$ は、$x = 0$ で $y(0) = a_0$ となり、解は $x = 0$ で有限値を取ります。しかし確定特異点では、解が $x = 0$ で発散する可能性があります。たとえば $x^2 y” + xy’ – y = 0$ の解の一つは $y = x^{-1}$ であり、$x \to 0$ で発散します。
通常のべき級数 $\sum a_n x^n$ では $x^{-1}$ のような負のべきを表現できません。そこでフロベニウスは、べき級数の前に $x^r$ を掛けるという一般化を考えました。$r$ を実数（場合によっては複素数）とすることで、発散する解も正しく捕捉できるようになります。
確定特異点での解法
確定特異点 $x = 0$ では、一般化べき級数（フロベニウス級数）
$$ y = x^r \sum_{n=0}^{\infty} a_n x^n = \sum_{n=0}^{\infty} a_n x^{n+r} $$
を仮定します。$a_0 \neq 0$ とします（$a_0 = 0$ なら $r$ の値を変えて吸収できる）。
指数方程式の導出
一般的な確定特異点の方程式
$$ x^2 y” + x\,p(x)\,y’ + q(x)\,y = 0 $$
ここで $p(x) = p_0 + p_1 x + \cdots$、$q(x) = q_0 + q_1 x + \cdots$ は $x = 0$ で解析的です。フロベニウス級数を代入すると、最低次 $x^r$ の係数から
$$ r(r-1) + p_0 r + q_0 = 0 $$
が得られます。これが指数方程式（indicial equation）です。指数方程式は $r$ についての2次方程式であり、2つの根 $r_1, r_2$（$\text{Re}(r_1) \geq \text{Re}(r_2)$）が得られます。
指数方程式の2つの根の差 $r_1 – r_2$ の値に応じて、2つの独立解の求め方が変わります。
- $r_1 – r_2$ が整数でない場合: 2つの根のそれぞれに対してフロベニウス級数解が得られ、それらが2つの独立解になります
- $r_1 = r_2$ の場合: 第1解はフロベニウス級数で得られますが、第2解には対数項 $\ln x$ が含まれます
- $r_1 – r_2$ が正の整数の場合: $r_1$ からは必ずフロベニウス解が得られます。$r_2$ からの解は対数項を含む場合と含まない場合があります
この分類はODEの理論において重要な定理（Fuchsの定理の拡張）として知られています。
ベッセルの方程式
ベッセルの方程式は、円筒座標でのラプラス方程式やヘルムホルツ方程式を変数分離したときに現れます。たとえば円形膜（太鼓）の振動モードを求めると、動径方向の関数がベッセルの方程式を満たします。
$$ x^2 y” + xy’ + (x^2 – \nu^2)y = 0 $$
まず特異点の分類を行います。標準形に変換すると $y” + (1/x)y’ + (1 – \nu^2/x^2)y = 0$ なので $P(x) = 1/x$、$Q(x) = 1 – \nu^2/x^2$ です。$xP(x) = 1$、$x^2Q(x) = x^2 – \nu^2$ はいずれも $x = 0$ で解析的なので、$x = 0$ は確定特異点です。
フロベニウス級数 $y = \sum_{n=0}^{\infty} a_n x^{n+r}$ を代入します。各導関数は
$$ y’ = \sum_{n=0}^{\infty} (n+r) a_n x^{n+r-1}, \quad y” = \sum_{n=0}^{\infty} (n+r)(n+r-1) a_n x^{n+r-2} $$
元の方程式に代入すると
$$ \sum_{n=0}^{\infty} (n+r)(n+r-1) a_n x^{n+r} + \sum_{n=0}^{\infty} (n+r) a_n x^{n+r} + \sum_{n=0}^{\infty} a_n x^{n+r+2} – \nu^2 \sum_{n=0}^{\infty} a_n x^{n+r} = 0 $$
第1項・第2項・第4項をまとめて $[(n+r)^2 – \nu^2]a_n x^{n+r}$ とし、第3項の添字を $n \to n-2$ にシフトすると
$$ \sum_{n=0}^{\infty} [(n+r)^2 – \nu^2] a_n x^{n+r} + \sum_{n=2}^{\infty} a_{n-2} x^{n+r} = 0 $$
$n = 0$: $(r^2 – \nu^2)a_0 = 0$ であり、$a_0 \neq 0$ なので指数方程式
$$ r^2 – \nu^2 = 0 \implies r = \pm\nu $$
が得られます。
$n = 1$: $[(1+r)^2 – \nu^2]a_1 = 0$ です。$r = \nu$ のとき $(1+2\nu)a_1 = 0$ なので（$\nu \geq 0$ なら）$a_1 = 0$ です。
$n \geq 2$: 漸化式は
$$ [(n+r)^2 – \nu^2] a_n + a_{n-2} = 0 \quad \Rightarrow \quad a_n = -\frac{a_{n-2}}{(n+r)^2 – \nu^2} $$
$r = \nu$ を代入すると
$$ a_n = -\frac{a_{n-2}}{n(n + 2\nu)} $$
$a_1 = 0$ より奇数番目の係数はすべてゼロです。偶数番目の係数を $n = 2m$ として
$$ a_{2m} = -\frac{a_{2m-2}}{2m(2m + 2\nu)} = -\frac{a_{2m-2}}{4m(m + \nu)} $$
この漸化式を繰り返し適用すると
$$ a_{2m} = \frac{(-1)^m a_0}{4^m \cdot m! \cdot (\nu+1)(\nu+2)\cdots(\nu+m)} $$
ガンマ関数の性質 $\Gamma(\nu+m+1) = (\nu+m)(\nu+m-1)\cdots(\nu+1)\Gamma(\nu+1)$ を用いて $a_0 = 1/[2^\nu \Gamma(\nu+1)]$ と選ぶと
$$ a_{2m} = \frac{(-1)^m}{m! \cdot \Gamma(\nu + m + 1) \cdot 2^{2m+\nu}} $$
これが第1種ベッセル関数 $J_\nu(x)$ の級数表現です。
$$ J_\nu(x) = \sum_{m=0}^{\infty} \frac{(-1)^m}{m!\, \Gamma(\nu + m + 1)}\left(\frac{x}{2}\right)^{2m+\nu} $$
$\nu = 0$ の場合を具体的に書き下すと
$$ J_0(x) = 1 – \frac{x^2}{4} + \frac{x^4}{64} – \frac{x^6}{2304} + \cdots $$
級数の各項が交互に正負を取るため、ベッセル関数が振動的な性質を持つことが級数表現から読み取れます。
ベッセルの方程式の解を導出できたところで、次はもう一つの重要な特殊関数であるルジャンドル多項式を級数解法で導出しましょう。
ルジャンドルの方程式
ルジャンドルの方程式は、球座標でのラプラス方程式を変数分離したときに極角 $\theta$ に関する部分から現れます。$x = \cos\theta$ と変数変換すると、極角方向の関数がルジャンドルの方程式を満たします。
$$ (1-x^2)y” – 2xy’ + l(l+1)y = 0 $$
$x = 0$ は正則点です（$P(x) = -2x/(1-x^2)$、$Q(x) = l(l+1)/(1-x^2)$ はいずれも $x = 0$ で解析的）。
べき級数 $y = \sum_{n=0}^{\infty} a_n x^n$ を代入します。
$$ (1-x^2)\sum_{n=2}^{\infty} n(n-1)a_n x^{n-2} – 2x\sum_{n=1}^{\infty} na_n x^{n-1} + l(l+1)\sum_{n=0}^{\infty} a_n x^n = 0 $$
第1項を展開すると
$$ \sum_{n=2}^{\infty} n(n-1)a_n x^{n-2} – \sum_{n=2}^{\infty} n(n-1)a_n x^n – \sum_{n=1}^{\infty} 2na_n x^n + \sum_{n=0}^{\infty} l(l+1)a_n x^n = 0 $$
第1項で $m = n – 2$ と置換し、残りの項で $m = n$ とすると
$$ \sum_{m=0}^{\infty} (m+2)(m+1)a_{m+2} x^m + \sum_{m=0}^{\infty} [-m(m-1) – 2m + l(l+1)]a_m x^m = 0 $$
第2項の角括弧を整理すると $-m^2 – m + l(l+1) = l(l+1) – m(m+1) = (l-m)(l+m+1)$ となるので、漸化式は
$$ a_{m+2} = -\frac{(l-m)(l+m+1)}{(m+2)(m+1)} a_m $$
ここで重要な観察があります。$l$ が非負整数のとき、$m = l$ で分子が $(l – l)(l + l + 1) = 0$ となるため $a_{l+2} = 0$ となり、それ以降の係数もすべてゼロになります。つまり級数が有限項で打ち切られ、多項式になるのです。
これがルジャンドル多項式 $P_l(x)$ です。最初の数個を具体的に書くと
$$ P_0(x) = 1, \quad P_1(x) = x, \quad P_2(x) = \frac{3x^2 – 1}{2}, \quad P_3(x) = \frac{5x^3 – 3x}{2} $$
ルジャンドル多項式はロドリゲスの公式でも表現できます。
$$ P_l(x) = \frac{1}{2^l l!}\frac{d^l}{dx^l}(x^2 – 1)^l \quad \text{（ロドリゲスの公式）} $$
さらにルジャンドル多項式は $[-1, 1]$ 上の直交多項式系を構成し
$$ \int_{-1}^{1} P_m(x)P_n(x)\,dx = \frac{2}{2n+1}\delta_{mn} $$
この直交性は、球面調和関数の展開やルジャンドル級数展開の基盤となります。
級数解法と収束の問題
べき級数解法で得られた級数が実際に意味のある解を与えるためには、級数が収束しなければなりません。収束半径は展開中心から最も近い特異点までの距離で決まります。
たとえばルジャンドルの方程式では、$x = \pm 1$ が特異点（$1 – x^2 = 0$）なので、$x = 0$ の周りの展開の収束半径は1です。ベッセルの方程式では、$x = 0$ 以外に有限距離の特異点がないため、収束半径は無限大であり、任意の $x$ で収束します。
この収束半径の情報は実用上も重要です。収束が遅い領域では多くの項数が必要となり、数値計算のコストが増大します。Pythonの実装では、このような場合に漸近展開や数値的ODE解法との使い分けが必要になります。
それでは、Pythonで級数解法の結果を可視化し、特殊関数の性質を確認してみましょう。
Pythonでの実装
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import jv, airy, legendre
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# (a) エアリ関数（級数解から計算）
ax = axes[0, 0]
x = np.linspace(-15, 5, 1000)
ai, aip, bi, bip = airy(x)
ax.plot(x, ai, "b-", linewidth=2, label="Ai(x)")
ax.plot(x, bi, "r-", linewidth=2, label="Bi(x)")
ax.set_xlabel("x", fontsize=12)
ax.set_ylabel("y", fontsize=12)
ax.set_title("Airy Functions: $y'' - xy = 0$", fontsize=13)
ax.legend(fontsize=10)
ax.set_ylim(-1, 2)
ax.grid(True, alpha=0.3)
# (b) ベッセル関数
ax = axes[0, 1]
x_bes = np.linspace(0, 20, 500)
for nu in [0, 1, 2, 3]:
ax.plot(x_bes, jv(nu, x_bes), linewidth=2, label=f"$J_{{{nu}}}(x)$")
ax.set_xlabel("x", fontsize=12)
ax.set_ylabel("$J_\\nu(x)$", fontsize=12)
ax.set_title("Bessel Functions of the First Kind", fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
ax.set_ylim(-0.5, 1.1)
# (c) ルジャンドル多項式
ax = axes[1, 0]
x_leg = np.linspace(-1, 1, 500)
for l in range(6):
Pl = legendre(l)
ax.plot(x_leg, Pl(x_leg), linewidth=2, label=f"$P_{{{l}}}(x)$")
ax.set_xlabel("x", fontsize=12)
ax.set_ylabel("$P_l(x)$", fontsize=12)
ax.set_title("Legendre Polynomials", fontsize=13)
ax.legend(fontsize=9, ncol=2)
ax.grid(True, alpha=0.3)
# (d) 級数の部分和の収束（ベッセル関数）
ax = axes[1, 1]
x_conv = np.linspace(0, 10, 300)
nu = 0
for N in [1, 3, 5, 10, 20]:
y_partial = np.zeros_like(x_conv)
for m in range(N):
coeff = (-1)**m / (np.math.factorial(m) * np.math.factorial(m + int(nu)))
y_partial += coeff * (x_conv/2)**(2*m + nu)
ax.plot(x_conv, y_partial, linewidth=1.5, label=f"N = {N}")
ax.plot(x_conv, jv(0, x_conv), "k--", linewidth=2, label="Exact $J_0(x)$")
ax.set_xlabel("x", fontsize=12)
ax.set_ylabel("$J_0(x)$ partial sum", fontsize=12)
ax.set_title("Convergence of Power Series for $J_0(x)$", fontsize=13)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
ax.set_ylim(-1, 1.5)
plt.tight_layout()
plt.savefig("ode_series_solution.png", dpi=150, bbox_inches="tight")
plt.show()
このグラフから、級数解法で得られる特殊関数の性質が読み取れます。
-
左上（エアリ関数）: $y” – xy = 0$ の2つの独立な解 Ai(x) と Bi(x) です。$x > 0$ では指数的に減衰/増大し、$x < 0$ では振動します。これは方程式の形から直感的に理解できます。$x > 0$ のとき $y” = xy > 0$（$y > 0$ なら）で解は上に凸にならず指数的に変化しますが、$x < 0$ のとき $y'' = xy$ は $y$ と逆符号になり振動的になります。Ai(x) は $x \to +\infty$ で減衰するためWKB近似との接続関数として量子力学で重要な役割を果たします
-
右上（ベッセル関数）: 次数 $\nu$ が大きいほど原点近くでの立ち上がりが遅く、ゼロ点の位置もシフトしています。これは級数表現 $J_\nu(x) \sim (x/2)^\nu/\Gamma(\nu+1)$（$x \to 0$）から理解できます。$\nu$ が大きいと $x^\nu$ の立ち上がりが遅いためです。また振動しながら振幅が減衰していく様子は、遠方では $J_\nu(x) \approx \sqrt{2/(\pi x)}\cos(x – \nu\pi/2 – \pi/4)$ という漸近展開で記述されます
-
左下（ルジャンドル多項式）: $l$ 次のルジャンドル多項式は $[-1, 1]$ 上でちょうど $l$ 個のゼロ点を持ち、直交多項式系を構成しています。$P_0 = 1$（定数）、$P_1 = x$（直線）、$P_2 = (3x^2-1)/2$（放物線的）と次数が上がるにつれて振動が増えていく様子が見えます。この直交性はガウス求積法の標本点としてルジャンドル多項式のゼロ点が使われる理由でもあります
-
右下（級数の収束）: $J_0(x)$ のべき級数が項数 $N$ の増加とともに真の値に収束する様子が見えます。$N = 1$ や $N = 3$ では $x$ が大きい領域で大きなずれがありますが、$N = 20$ でほぼ完全に一致しています。ベッセルの方程式の唯一の特異点は $x = 0$（とそのフロベニウス級数の展開中心）であり、有限距離に他の特異点がないため、収束半径は無限大です。したがって十分な項数を取ればどの $x$ でも収束します
応用例 — 物理学・工学での活用
級数解法で得られる特殊関数は、物理学と工学の至る所に現れます。ここではいくつかの具体例を紹介します。
量子力学 — 水素原子
シュレーディンガー方程式を球座標で変数分離すると、極角方向にはルジャンドルの方程式（一般にはルジャンドル陪関数の方程式）、動径方向にはラゲールの陪方程式が現れます。級数解法により、$l$ や $n$（量子数）が整数の場合にのみ級数が有限で打ち切られて物理的に許容される解（波動関数）が得られ、これがエネルギーの量子化条件を与えます。
電磁気学 — 円筒導波管
マイクロ波やミリ波の伝送に使われる円筒導波管では、管内の電磁界分布がベッセル関数で記述されます。管壁での境界条件がベッセル関数のゼロ点で決まり、伝搬可能なモード（TE/TM モード）のカットオフ周波数が決定されます。
振動工学 — 円形膜の振動
太鼓のような円形膜の振動モードは、動径方向がベッセル関数、角度方向が三角関数で表されます。膜の固定端での境界条件 $J_n(k a) = 0$（$a$ は膜の半径）から固有振動数が決まります。
重力場 — 球対称ポテンシャル
地球や惑星の重力場を球面調和関数で展開する際、ルジャンドル多項式が展開基底になります。人工衛星の軌道摂動を計算する際にはルジャンドル多項式の具体的な値が必要になります。
まとめ
本記事では、級数解法による常微分方程式の解法について解説しました。
- べき級数解法: 解を $\sum a_n x^n$ とおいて微分方程式に代入し、各 $x^n$ の係数を比較して漸化式を導出する
- 正則点（係数関数が解析的な点）ではべき級数解が必ず存在し、収束半径は最も近い特異点までの距離以上
- 確定特異点では通常のべき級数解が存在しない可能性があるが、フロベニウス法 $y = x^r \sum a_n x^n$ を使えば、指数方程式から $r$ を決定して級数解を構成できる
- ベッセル関数 $J_\nu(x)$ はベッセルの方程式のフロベニウス解として、ルジャンドル多項式 $P_l(x)$ はルジャンドルの方程式のべき級数解が有限で打ち切られたものとして、エアリ関数はエアリの方程式のべき級数解として導出される
- 級数の収束半径は特異点の配置で決まり、項数を増やすと解が真の値に収束する。収束半径の外では級数は発散するため、漸近展開など別の手法が必要
次のステップとして、以下の記事も参考にしてください。
- 熱伝導方程式の導出と解法 — 偏微分方程式への展開
- 常微分方程式の数値解法まとめ — 数値計算手法
- フーリエ級数の導出と収束条件 — 直交関数展開の基礎