ヌメロフ法

\( y''(z) = P(z) = g(z) y(z) \)というようなタイプの微分方程式の数値解を得ようとする時には、ヌメロフ法が便利だと思う。

ヌメロフ法の簡単な説明としては、まずテイラー展開で \begin{align} y_{n+1} + y_{n-1} - 2 y_{n} = \frac{2}{2!} y''(z) h^{2} + \frac{2}{4!} y''''(z) h^{4} + \mathcal{O}(h^{6}) \; , \end{align} とする。 \( y_{n+1} = y(z+h), \; y_{n-1} = y(z-h) \)
そして\( y''''(z) \)については、 \( y''(z) = P(z) \)であることを使って、 \begin{align} P_{n+1} + P_{n-1} - 2 P_{n} = \frac{2}{2!} P''(z) h^{2} + \mathcal{O}(h^{4} ) \; , \\ h^{2} y''''(z) = h^{2} P''(z) = P_{n+1} + P_{n-1} - 2P_{n} - \mathcal{O}(h^{4 } ) \; , \end{align} となることがわかる。それを利用して、 \begin{align} y_{n+1} + y_{n-1} - 2 y_{n} & = P_{n} h^{2} + \frac{ 2 h^{2} }{4! } \left( P_{n+1} + P_{n-1} - 2P_{n} - \mathcal{O}(h^{4 } ) \right) + \mathcal{O}(h^{6 } ) \; , \\ & = \frac{h^{2} }{12 } P_{n+1} + \frac{5h^{2}}{6} P_{n} + \frac{h^{2}}{12 } P_{n-1} + \mathcal{O}(h^{6} ) \; \end{align} となることがわかる。よって、 \begin{align} \left( 1 - \frac{h^{2} }{12} g_{n+1} \right) y_{n+1} + \left( 1 - \frac{h^{2} }{12} g_{n-1} \right) y_{n-1} = \left( 2 - \frac{5 h^{2 } }{6 } g_{n} \right) y_{n} + \mathcal{O}(h^{6} ) \; \end{align} を得る。これがヌメロフ法である。

ヌメロフ法を使ってこの手の微分方程式を解いてみる。 \begin{align} \frac{d^{2} \chi_{n} }{d z^{2}} + \left( \frac{ - \frac{\alpha^{2} }{4} + \frac{1}{4} }{z^{2} } + \frac{ n + \frac{\alpha}{2} + \frac{1 }{2} }{z} - \frac{1}{4} \right) \chi_{n} = 0 \; \end{align} といった方程式はどうだろう。この解析解は \begin{align} \chi_{n } (z) = z^{ \frac{\alpha + 1}{2} } e^{ \frac{-z}{2} } L^{(\alpha ) }_{n } (z) \end{align} である。
例として、 \( \alpha=1.0 \; , \; N=2 \)と \( \alpha=2.5 \; , N=3 \), \( \alpha=4.5, \; N=5 \), \( \alpha=5, \; N=7 \) の場合では、

\( \alpha=1.0, \; N=10 \)　と　\( \alpha=1.0, \; N=20 \)の場合は

のようになり、Nが大きくなるほど収束は遅くなることがわかる。
全体的に解析解と数値解はかなり近い値となっている。計算速度もかなり速かったと思う。