[VBA] 合流型超幾何関数とラゲール陪多項式

合流型超幾何関数

 次の級数で定義される関数を Kummer の合流型超幾何関数 (Confluent hypergeometric function) といいます。
 
\[F(a,c\:;x)=1+\frac{a}{c}\frac{x}{1!}+\frac{a(a+1)}{c(c+1)}\frac{x^2}{2!}+\cdots=\sum_{k=0}^{\infty}\frac{(a)_k}{(c)_k}\frac{x^k}{k!}\]
 指数関数や エルミート多項式ベッセル関数 などを包括する広いクラスの関数です(たとえば $a=c$ のときに指数関数のマクローリン級数と一致します)。定義式から明らかなように、この関数は $c$ が 0 または負整数のときは分母が 0 となって値をもちません。また、$a$ が負整数であるときには級数は有限項で打ち切られます。たとえば $a=-2$ のときは第 4 項以降はすべて 0 になるので、$F$ は 2 次関数となります。

 VBA 数値計算用最新マクロ をモジュールに貼り付ける(あるいは上書きする)と、この記事にあるマクロも含めて当サイトに載っている主要な関数を全て利用することができます。ぜひお試しください。
 

合流型超幾何関数のマクロ

 合流型超幾何級数を計算する Function Macro(ユーザー定義関数)です。

'[VBA] 合流型超幾何級数

Function CHGEO(a As Double, c As Double, x As Double) As Double

  Dim d As Double
  Dim s As Double

  d = a * x / c

  s = 1 + d

  For k = 1 To 200
    d = d * x * (a + k) / (c + k) / (k + 1)
    s = s + d
  Next k

  CHGEO = s

End Function

 この関数を呼びだすときには、

=CHGEO(a,c,x)

というように記述します。たとえば

=CHGEO(1,1,2)

と入力すると $F(1,1\:;2)=e^2$ を計算して 7.389 という値を返します。
 

複素変数の合流型超幾何関数

 一般には複素変数の合流型超幾何関数が用いられます。
 
\[F(a,c\:;z)=1+\frac{a}{c}\frac{z}{1!}+\frac{a(a+1)}{c(c+1)}\frac{z^2}{2!}+\cdots=\sum_{k=0}^{\infty}\frac{(a)_k}{(c)_k}\frac{z^k}{k!}\]
 マクロは次のようになります。

'VBA 複素変数の合流超幾何関数

Function CXCHGEO(a As Double, c As Double, z As Complex) As Complex

  Dim d As Double
  Dim k As Double
  Dim s As Complex
  Dim zk As Complex

  d = a / c

  s.re = 1 + d * z.re
  s.im = d * z.im
  zk = z

  For k = 2 To 100
    d = d * (a + k - 1) / (c + k - 1) / k
    zk = CPXCLC(zk, z, 3)
    s.re = s.re + d * zk.re
    s.im = s.im + d * zk.im
  Next k

  CXCHGEO = s

End Function

 変数 z は複素数型 (Complex型) を指定します。この関数を使って、
 
\[F(1,1,1;i)=\sum_{k=0}^{\infty}(1+i)^k\]
を数値計算してみましょう。

'VBA (1+i)^kの計算

Sub CXCHGEOtest()

 Dim fz As Complex

 fz = CXCHGEO(1, 1, CPX(1, 1))

 Debug.Print fz.re & " + " & fz.im & "i"

End Sub

 マクロを実行すると

1.46869393991589 + 2.28735528717884i

という計算結果が表示されます。
 

ラゲール陪多項式

 合流型超幾何関数には、量子力学において水素原子の波動関数を記述するときに用いられる ラゲール陪多項式 が含まれています。ラゲール陪多項式は母関数によって
 
\[g(t,x)=\frac{(-1)^k}{(1-t)^{k+1}}\exp\left(-\frac{xt}{1-t}\right)=\sum_{n=k}^{\infty}L_n^k(x)\frac{t^{n-k}}{n!}\]
と定義されています。つまり $g(t,x)$ を級数展開したときの $t^{n-k}/n!$ の係数が $L_n^k(x)$ です。合流型超幾何関数を用いると
 
\[L_n^k(x)=(-1)^k\frac{(n!)^2}{k!(n-k)!}F(k-n,k+1\:;x)\]
と表されます。特に $k=0$ のときには ラゲール多項式 とよばれ、$L_n^0(x)$ の添え字 0 を省略して単に $L_n(x)$ と書くこともあります。ラゲール陪多項式 $L_n^k(x)$ はラゲール多項式の $k$ 階微分という関係になっています。
 
\[L_n^k(x)=\frac{d^k}{dx^k}L_n(x)\]
 $L_n(x)$ を具体的に書き下してみると
 
\[\begin{align*}L_0(x)&=1\\[6pt]
L_1(x)&=-x+1\\[6pt]
L_2(x)&=x^2-4x+2\\[6pt]
L_3(x)&=-x^3+9x^2-18x+6\end{align*}\]
 それぞれ微分すると $L_n^1(x)$ が得られます。
 
\[\begin{align*}L_1^1(x)&=-1\\[6pt]
L_2^1(x)&=2x-4\\[6pt]
L_3^1(x)&=-3x^2+18x-18\end{align*}\]
 これをさらに繰り返し微分していくと、ラゲール陪多項式 $L_n^k(x)$ の系列が得られます。
 

ラゲール陪多項式のマクロ

 ラゲール陪多項式を計算するマクロは以下のようになります。

'VBA ラゲール陪多項式

Function LAG(n As Integer, k As Integer, x As Double)

  nfact = WorksheetFunction.fact(n)
  kfact = WorksheetFunction.fact(k)
  nkfact = WorksheetFunction.fact(n - k)

  LAG = (-1) ^ k * nfact ^ 2 * CHGEO(k - n, k + 1, x) / kfact / nkfact

End Function

 このユーザー定義関数は

=LAG(n,k,x)

の形で呼びだします。たとえば

=LAG(2,1,3)

と入力すると $L_2^1(3)$ を計算します。このマクロを使ってワークシートにラゲール多項式 $L_n(x)$ のグラフを描くと次のようになります。

 VBAラゲール多項式のグラフ

 $L_n(x)$ の導関数である $L_n^1(x)$ のグラフは次のようになります。

 VBAラゲールの陪多項式グラフ


 

コメント

  1. 落合史生 より:

    この関数を呼びだすときには、
    =CGHEO(a,c,x)

    =CHGEO(a,c,x)の間違いでは。
    =CGHEO(1,1,2)も。

    ラゲール多項式で検索していて、このページを見つけました。
    VBAでグラフが書けるので、多項式の特徴を知るのに助かります。

    • Blog Cat より:

       御指摘の通り、正しくは
      =CHGEO(a,c,x)
      です。記事は訂正しておきました。申し訳ありませんでした。

Blog Cat へ返信する コメントをキャンセル

メールアドレスが公開されることはありません。