ベッセル関数の級数展開とベッセルの微分方程式の導出【特殊関数】

ベッセル関数と呼ばれる特殊関数について解説します。

ガンマ関数やベータ関数と比べると式は複雑で、出会う場面も稀な特殊関数です。

とは言え、ベッセル関数は有用な関数であることも事実です。

例えば、波動方程式の解の表示や人工衛星の軌道計画にベッセル関数は用いられており、現代文明を支える重要な関数の一つです。

今回は現代文明を陰ながら支えるベッセル関数について解説します。

ケプラー方程式とベッセル関数

惑星の運動を時間の関数として表すとケプラー方程式と呼ばれる超越方程式が現れます。

$M$は平均近点離角と呼ばれ、惑星の平均運動から簡単に求められます。一方、$E$は離心近点離角であり、実際に惑星が位置する場所を表します。

実用上は$E$の方が重要であるため、$M$が与えられたとき$E$を求めたい場合がほとんどです。

しかしながら、ケプラー方程式を見れば分かるように、$E$が三角関数の中に含まれているため、簡単に離心近点離角を求めることはできません。

ケプラー方程式を解くとは、$M$が与えられたときに$E$を求めることです。

多くの人々がケプラー方程式を解こうと挑戦しましたが、ケプラー方程式は超越方程式の一種であり、解析的には解けないことが後に判明しました。

とは言え、解く方法を何とか探さなければなりません。

この要請を受け、ケプラー方程式を解く方法として生み出されたのがベッセル関数です。

さて、$E$は（第一種）ベッセル関数により次のように求められることが分かっています。

歴史的には、ケプラー方程式からベッセル関数が導かれましたが、現在では天体力学以外の多くの分野でも利用される汎用的な関数です。

ベッセル関数の定義

実際にケプラー方程式からベッセル関数を見出してみましょう。

まず、ケプラー方程式を次のように変形すると、

\begin{eqnarray}
e\sin E = E\,-M \tag{1}
\end{eqnarray}

とできて、左辺は三角関数であるため、右辺も周期関数であることが分かります。

周期関数であるため、右辺をフーリエ級数展開できます。（→フーリエ級数展開の公式）

この周期関数の周期は$2\pi$であり、$E$を$M$により表すことを考慮して、フーリエ級数展開を行うと、次のように展開できます。（$a_0, a_k, b_k$は定数）

さて、式(1)について改めて考えると、左辺が正弦関数であることから、$E\,-M$も正弦関数となるはずです。したがって、上式の余弦関数は無視でき、以下のように簡単に出来ます。

$a_0, b_k$は定数のため、両辺を$M$で微分して、

となります。（$a_0$の項を消去できました）

ここから$b_k$を取り出すため、テクニックを使います。

そのテクニックとは、三角関数の直交性を利用することです。（→三角関数の直交性とは？）

すなわち、$\cos lM$を両辺に掛けて$0$から$\pi$まで積分を実行すると、

となって、$b_k$の項に関しては、$k=l$を満たす項だけが残ります。

従って$k=l$の場合にのみ注目して式を整理すると、

とできて、$b_k$を求める目途が付きました。

しかし、このままでは計算ができないため、積分変数を$M$から$E$に変更します。

$M$と$E$の対応関係は表のようになり、

\begin{array}{c|ccc}
\,M & 0 & \cdots & \pi\,\,\, \\\hline
\,E & 0 & \cdots & \pi\,\,\,
\end{array}

また変数変換後も変数の範囲は変わらないため、次のように上式は変換できて、

\begin{split}
\ff{\pi }{2}kb_k &= \int_0^{\pi} \cos kM \diff E
\end{split}

となります。

この式に関して$M = E \,- e\sin E$であることを思い出すと、

\begin{eqnarray}
\ff{\pi }{2}kb_k = \int_0^{\pi} \cos k( E \,- e\sin E) \diff E
\end{eqnarray}

とできます。以上より$b_k$が次のように求められます。

\begin{eqnarray}
b_k = \ff{2}{\pi k} \int_0^{\pi} \cos ( kE \,- ek\sin E) \diff E
\end{eqnarray}

ここで、$x = ek$として次のようにベッセル関数$J_k(x)$なるものを定義します。

ベッセル関数を用いると$b_k$は、

\begin{eqnarray}
b_k = \ff{2}{k} J_k(x)
\end{eqnarray}

と表せます。

さて、$a_0$については$M=0$のとき、$E=0$であることから$0$となります。

以上より、ケプラー方程式はベッセル関数を用いて、

と級数展開できます。

ベッセル関数の展開

ベッセル関数を導出できましたが、依然として$E$が中に含まれているため、具体的には計算できません。

このままでは不完全なため、ベッセル関数について考察を深めていきましょう。

まず、先程のベッセル関数を加法定理に従って展開します。

右辺第一項と第二項に分けて、計算を進めていきましょう。なお、$k$は整数とします。

右辺第一項の計算

右辺第一項について、積分区間を次のように分割します。

ここで右辺第二項について$E’ = \pi \,- E$として積分変数を変換します。

積分区間は表のようになり、$\diff E = -\diff E’$となります。

\begin{array}{c|ccc}
\,E & \DL{\ff{\pi}{2}} & \cdots & \pi\,\,\, \\\hline
\,E’ & \DL{\ff{\pi}{2}} & \cdots & 0 \,\,\,
\end{array}

表より、次のように変形が行えます。

$\cos k\pi = (-1)^k$、$\sin k\pi = 0$なので、右辺第一項は最終的に次のように整理できます。（$E’$を$E$に書き換えてます）

上式より$k$が奇数のとき右辺第一項は$0$となることが分かります。

右辺第二項の計算

右辺第二項に関しても積分区間を$\DL{\ff{\pi}{2}}$で分割して、同様に$E’ = \pi \,- E$として式変形を行います。

すると、

これより、$k$が偶数のとき右辺第二項は$0$になることが分かります。

ベッセル関数の級数表示

先ほどの計算結果から、ベッセル関数は$k$が奇数と偶数の場合で次のように分類できることが分かります。

ベッセル関数と三角関数の級数展開

まず、$k$が奇数の場合でのベッセル関数の級数展開を検討します。

$k$が奇数のとき、ベッセル関数は右辺第二項のみになります。

\begin{eqnarray}
J_k(x) &=& \ff{1}{\pi}\int_0^{\pi}\sin kE \cdot\sin(x\sin E) \diff E \tag{a} \EE
\end{eqnarray}

この積分を実行したいのですが、$\sin(x\sin E)$が邪魔で積分を行えません。

そのため、$\sin(x\sin E)$の級数展開を行い、積分が実行できるように工夫します。

$\sin(x\sin E)$を級数展開すると次式になります。（→三角関数の級数展開公式）

この式に取り組む前に、$m$を奇数として$\sin^m E$に関して考察しましょう。（$m = 2j-1$）

さて、オイラーの公式より、$\sin \theta$は$\DL{\ff{e^{i\theta}\,- e^{-i\theta} }{2i}}$と表せるので、$\sin^m E$は

となります。（→オイラーの公式を利用した三角関数の表示）

さて、$m$は奇数のため、次のように計算できて、

\begin{eqnarray}
\left( \ff{1}{2i} \right)^m &=& \ff{1}{2^m}\ff{1}{i\cdot i^{m-1}} \EE
&=& \ff{1}{2^m}\ff{i^{m-1}}{i\cdot i^{2(m-1)}} \EE
&=& \ff{1}{2^m}\ff{(-1)^{\ff{m-1}{2}}}{i} \EE
\end{eqnarray}

これより、

となります。

二項定理の適用

式(2)の累乗に関して二項定理を用いると、

とできて、$s$を$0$以上の整数として、組み合わせ（combination）の性質より$\DL{{}_mC_{s} = {}_mC_{m-s}}$であり、また、$(-1)^{m-s} = -(-1)^s $ であることを利用して式(2)に代入すると、

$\sin(x\sin E)$の級数展開

さて、二項定理による計算結果を式(1)に適用すると、

となります。（ただし、$m=2j-1$）

さて、式(a)の$\DL{\ff{1}{\pi}\int_0^{\pi}\sin kE \cdot\sin(x\sin E) \diff E}$の積分計算を行うために、$\sin(x\sin E)$を級数展開しましたが、ここで三角関数の直交性を利用すると、$m \,- 2s = k$を満たす項のみが残ることが事前に判断できるので、$\sin(x\sin E)$の項の内、

となる多項式をベッセル関数の計算に用いれば良いことが分かります。

ベッセル関数の級数展開

以上より、ベッセル関数の中身は次のように書き下せます。

式変形において、${}_{k+2s}C_s = \DL{\ff{(k+2s)!}{s!(k+s)!}}$であることを利用しています。

以上より、$k$が奇数のときのベッセル関数の級数表示は、次のようになります。

ベッセル関数の級数表示

次に$k$が偶数の場合について考えましょう。

$k$が偶数のとき、ベッセル関数は、

\begin{eqnarray}
J_k(x) = \ff{1}{\pi}\int_0^{\pi} \cos kE\cdot\cos(x\sin E) \diff E \tag{3}
\end{eqnarray}

となります。

先程と同様に、$\cos(x\sin E)$についても級数展開を実行すると、

となって、以降も先程と同様に二項定理を利用して愚直に計算を行うと、(3)式は、

と級数展開できます。

$k$が偶数の場合でも奇数の場合でも、級数展開の結果が一致することが分かりました。

以上より、ベッセル関数の級数展開は次のようになります。

※$k$が負の整数の場合でもベッセル関数は計算でき、$J_{-k}(x) = (-1)^kJ_k(x)$の関係が成立します。

ベッセル関数のグラフ

ベッセル関数を具体的な級数で表すことに成功したので、グラフを描画してみましよう。

縦軸に$J_k(x)$の値、横軸に$x$を取り$k=0,1,2,3,4$までのベッセル関数を描画すると次のようになります。

ベッセル関数に三角関数が含まれていることから、ベッセル関数も周期関数であると見当が付きますが、実際に描画することでベッセル関数が周期関数であることが視覚的に理解できます。

さて、ベッセル関数のグラフから、$x$の増加に伴い振幅が小さくなることが分かります。この様子は減衰振動の様子を連想させます。（→減衰振動とその様子）

さらに、$k$が増加するにつれ、振幅の最大値も減少することが分かります。

このグラフの形状を観察すると、水滴が落下したときに、水面に発生させる同心円状の波の断面に似ていることに気が付きます。

実際、この直感は正しく、同心円状に広がる波の様子や、太鼓の膜の振動の様子はベッセル関数を用いることで記述することができます。

この例以外にも、振り子の周期について一般解を求めると、ベッセル関数が現れることが知られています。

ベッセル関数と微分方程式

いささか唐突ですが、次の微分方程式を考えます。

ベッセル関数がこの微分方程式の解になることを示します。

ベッセル関数を次のように一般的な形に改めます。

\begin{eqnarray}
J_k(x) &=& \ff{1}{\pi} \int_0^{\pi} \cos (k\theta \,-\, x \sin \theta) \diff \theta \EE
\end{eqnarray}

$y=J_k(x)$として、ベッセル関数を$x$の関数と見ると、一階微分は次のようになります。

\begin{eqnarray}
\ff{\diff J_k(x)}{\diff x} &= \ff{1}{\pi}\int_0^{\pi} \sin \theta\cdot\sin ( k\theta \,- x\sin \theta) \diff \theta \EE \tag{4}
\end{eqnarray}

ここで部分積分を利用すると、

となります。

次に二階微分について考えます。式(4)より、

となり、以上より、

となります。

最後に、次の微分方程式を考えます。

ここで$k\theta \,- x\sin\theta = \phi$とすると次の表のような対応関係が成立し、

\begin{array}{c|ccc}
\,\theta & 0 & \cdots & \pi\,\,\, \\\hline
\,\phi & 0 & \cdots & k\pi \,\,\,
\end{array}

$\DL{\diff \theta = \ff{\diff \phi}{k \,- x\cos\theta} }$となるので先ほどの積分は次のように計算できます。

\begin{split}
&\quad\ff{\diff^2 J_k(x)}{\diff x^2} + \ff{1}{x}\ff{\diff J_k(x)}{\diff x} + \left( 1\,-\ff{n^2}{x^2} \right)J_k(x) \EE
&= -\ff{k}{\pi x^2}\int_0^{\pi} (k-x\cos \theta)\cdot\cos( k\theta \,- x\sin \theta)\diff \theta \EE
&= -\ff{k}{\pi x^2}\int_0^{k\pi} (k-x\cos \theta)\cdot\cos\phi\cdot \ff{\diff \phi}{k \,- x\cos\theta} \EE
&= -\ff{k}{\pi x^2}\int_0^{k\pi}\cos\phi \,\diff \phi\EE
&= 0
\end{split}

となり、ベッセル関数は最初の微分方程式の解となることが示せました。

ベッセル関数がこの微分方程式のなるため、この微分方程式をベッセルの微分方程式と呼びます。