🎉

Symbolics.jlでロドリゲスの公式を展開しよう

2022/12/17に公開約36,800字

JuliaLang Advent Calendar 2022 17日目です🎄

Juliaの数式処理パッケージ Symbolics.jl を使って, 量子力学の教科書によく登場するRodriguesの公式を展開していきます. 下表の直交多項式はRodriguesの公式で定義されることが多いものの, 微分が入っていて扱いづらいため, 実装には漸化式や閉形式が用いられます. Rodriguesの公式は文献によって異なることがあり, どの定義に基づく実装なのか確かめるのはなかなか大変です. そこで, Symbolics.jlに頑張ってもらいましょう.

名前 Rodriguesの公式
エルミート多項式(Hermite)
H_{n}(x) = (-1)^n \mathrm{e}^{x^2} \frac{\mathrm{d}^n}{\mathrm{d}x^n} \mathrm{e}^{-x^2}
ルジャンドル多項式(Legendre)
P_n(x) = \frac{1}{2^n n!} \frac{\mathrm{d}^n}{\mathrm{d}x ^n} \left( x^2-1 \right)^n
ルジャンドル陪多項式(associated Legendre)
P_k^m(x) = \left( 1-x^2 \right)^{m/2} \frac{\mathrm{d}^m}{\mathrm{d}x^m} P_k(x)
ラゲール多項式(Laguerre)
L_n(x) = \frac{\mathrm{e}^x}{n!} \frac{\mathrm{d}^n}{\mathrm{d}x ^n} \left( \mathrm{e}^{-x} x^n \right)
ラゲール陪多項式(associated Laguerre)
L_n^{k}(x) = \frac{\mathrm{d}^k}{\mathrm{d}x^k} L_n(x)
一般化ラゲール多項式(generalized Laguerre)
L_n^{(\alpha)}(x)= \frac{x^{-\alpha}e^x}{n!} \frac{d^n}{dx^n}\left(x^{n+\alpha}e^{-x}\right)

※ 古い量子力学の教科書では, Laguerre多項式の定義に\frac{1}{n!}の因子を含まないL_n(x) = \mathrm{e}^x \frac{\mathrm{d}^n}{\mathrm{d}x ^n} \left( \mathrm{e}^{-x} x^n \right)が採用されることが多いように思われ, 当初はこの定義を採用するつもりでした. しかし, DLMFや, 一般化ラゲール多項式の定義との整合性を取るために\frac{1}{n!}の因子を加え, 上表の定義を採用することとしました.

パッケージ・環境

Symbolics.jlに加えて, 表示を整えるためにLatexify.jlLaTeXStrings.jlを用います.

入力
# using Pkg
# Pkg.add("Symbolics")
# Pkg.add("Latexify")
# Pkg.add("LaTeXStrings")
# Pkg.add("SpecialFunctions")
using Symbolics
using Latexify
using LaTeXStrings
using SpecialFunctions
versioninfo()
出力
Julia Version 1.8.3
Commit 0434deb161 (2022-11-14 20:14 UTC)
Platform Info:
  OS: Windows (x86_64-w64-mingw32)
  CPU: 8 × 11th Gen Intel(R) Core(TM) i7-1185G7 @ 3.00GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-13.0.1 (ORCJIT, tigerlake)
  Threads: 1 on 8 virtual cores

Symbolics.jlの使い方

公式ドキュメントによれば, Differential()で微分演算子を定義して, expand_derivatives()で微分を展開できます. 例えば\exp(-x^2)の微分は下記のように計算できます.

@variables x
D = Differential(x)

a = exp(-x^2)
b = D(a)
c = expand_derivatives(b)

display(a)
display(b)
display(c)
e^{ - x^{2}}
\frac{\mathrm{d} e^{ - x^{2}}}{\mathrm{d}x}
- 2 x e^{ - x^{2}}

2階微分は次のように計算できます.

@variables x
D = Differential(x)^2

a = exp(-x^2)
b = D(a)
c = expand_derivatives(b)

display(a)
display(b)
display(c)
e^{ - x^{2}}
\frac{\mathrm{d}^{2} e^{ - x^{2}}}{\mathrm{d}x^{2}}
- 2 e^{ - x^{2}} + 4 x^{2} e^{ - x^{2}}

上記のように計算した導関数に\exp(x^2)を掛けて整理すればHermite多項式が計算できます. simplify()のデフォルトだとうまく簡略化できないので, expand=trueのオプションを追加しましょう. すると, 下記のようにn=2のHermite多項式が求められます.

@variables x
n = 2
D = Differential(x)^n               # n階の微分演算子
a = (-1)^n * exp(x^2)               # 微分演算子の左側
b = exp(-x^2)                       # 微分演算子の右側
c = a * D(b)                        # Rodriguesの公式
d = expand_derivatives(c)           # 微分演算を展開
e = simplify(d, expand=true)        # 簡略化

display(a)
display(b)
display(c)
display(d)
display(e)
e^{x^{2}}
e^{ - x^{2}}
e^{x^{2}} \frac{\mathrm{d}^{2} e^{ - x^{2}}}{\mathrm{d}x^{2}}
\left( - 2 e^{ - x^{2}} + 4 x^{2} e^{ - x^{2}} \right) e^{x^{2}}
-2 + 4 x^{2}

D = Differential(x)^0でエラーが出ますので, n=0つまり0階微分のときは恒等写像x->xになるようにするとn=0n>0で共通のコードが使えるためすっきります. 三項演算子を使うとさらにシンプルに書けます.

@variables x
n = 0
D = n==0 ? x->x : Differential(x)^n # n階の微分演算子
a = (-1)^n * exp(x^2)               # 微分演算子の左側
b = exp(-x^2)                       # 微分演算子の右側
c = a * D(b)                        # Rodriguesの公式
d = expand_derivatives(c)           # 微分演算を展開
e = simplify(d, expand=true)        # 簡略化

display(a)
display(b)
display(c)
display(d)
display(e)
e^{x^{2}}
e^{ - x^{2}}
e^{x^{2}} e^{ - x^{2}}
e^{x^{2}} e^{ - x^{2}}
1

なお, Int64だと階乗はfactorial(20)までしか使えないので, factorial(21)以降を使う予定の場合はfactorial(big(21))のようにBigIntを使う必要があります. 追加のパッケージ等は不要です.

julia> factorial(20)
2432902008176640000

julia> factorial(21)
OverflowError: 21 is too large to look up in the table; consider using `factorial(big(21))` instead
Stacktrace:...

エルミート多項式

Hermite多項式のRodriguesの公式(上段)を展開し, 閉形式(下段)と比較します. 閉形式については, いつも通りに関数を定義し, そこに@variables xで定義したSymbolics.jlの変数を渡せばSymbolics.jlのオブジェクトとして返ってきます. 展開した結果はDLMF 18.5.18等と見比べてください.

\begin{aligned} H_{n}(x) &= (-1)^n \mathrm{e}^{x^2} \frac{\mathrm{d}^n}{\mathrm{d}x^n} \mathrm{e}^{-x^2} \\ &= n! \sum_{m=0}^{\lfloor n/2 \rfloor} \frac{(-1)^m}{m! (n-2m)!}(2 x)^{n-2m} \end{aligned}
for n in 0:6
    # Rodriguesの公式の展開
    @variables x
    D = n==0 ? x->x : Differential(x)^n # n階の微分演算子 
    a = (-1)^n * exp(x^2)               # 微分演算子の左側
    b = exp(-x^2)                       # 微分演算子の右側
    c = a * D(b)                        # Rodriguesの公式
    d = expand_derivatives(c)           # 微分演算を展開
    e = simplify(d, expand=true)        # 簡略化
    # 閉形式
    H(x; n=0) = factorial(n) * sum(m -> (-1)^m // (factorial(m)  * factorial(n-2*m)) * (2*x)^(n-2*m), 0:Int(floor(n/2)))
    # LaTeX表示
    display(latexstring(
    """
    \\begin{aligned}
        H_{$n}(x)
         = $(latexify(c, env=:raw))
        &= $(latexify(e, env=:raw)) \\\\
        &= $(latexify(H(x, n=n), env=:raw))
    \\end{aligned}
    """))
end
\begin{aligned} H_{0}(x) = e^{x^{2}} e^{ - x^{2}} &= 1 \\ &= 1 \end{aligned}
\begin{aligned} H_{1}(x) = - e^{x^{2}} \frac{\mathrm{d} e^{ - x^{2}}}{\mathrm{d}x} &= 2 x \\ &= 2 x \end{aligned}
\begin{aligned} H_{2}(x) = e^{x^{2}} \frac{\mathrm{d}^{2} e^{ - x^{2}}}{\mathrm{d}x^{2}} &= -2 + 4 x^{2} \\ &= -2 + 4 x^{2} \end{aligned}
\begin{aligned} H_{3}(x) = - e^{x^{2}} \frac{\mathrm{d}^{3} e^{ - x^{2}}}{\mathrm{d}x^{3}} &= 8 x^{3} - 12 x \\ &= 8 x^{3} - 12 x \end{aligned}
\begin{aligned} H_{4}(x) = e^{x^{2}} \frac{\mathrm{d}^{4} e^{ - x^{2}}}{\mathrm{d}x^{4}} &= 12 - 48 x^{2} + 16 x^{4} \\ &= 12 - 48 x^{2} + 16 x^{4} \end{aligned}
\begin{aligned} H_{5}(x) = - e^{x^{2}} \frac{\mathrm{d}^{5} e^{ - x^{2}}}{\mathrm{d}x^{5}} &= 32 x^{5} + 120 x - 160 x^{3} \\ &= 32 x^{5} + 120 x - 160 x^{3} \end{aligned}
\begin{aligned} H_{6}(x) = e^{x^{2}} \frac{\mathrm{d}^{6} e^{ - x^{2}}}{\mathrm{d}x^{6}} &= -120 + 720 x^{2} + 64 x^{6} - 480 x^{4} \\ &= -120 + 720 x^{2} + 64 x^{6} - 480 x^{4} \end{aligned}

ちゃんと一致していますね! 多項式の順番が少し変ですが, 値は合っているので, 閉形式がRodriguesの公式による定義と一致することの確認としては問題ありません.

ルジャンドル多項式

Legendre多項式のRodriguesの公式を展開し, 閉形式と比較します. 閉形式については, Legendre陪多項式のm=0の場合を用います. 展開した結果はDLMF 18.5.16等と見比べてください.

\begin{aligned} P_n(x) &= \frac{1}{2^n n!} \frac{\mathrm{d}^n}{\mathrm{d}x ^n} \left(x^2-1 \right)^n \\ &= P_n^0(x) \end{aligned}
# using Symbolics
# using Latexify
# using LaTeXStrings

for n in 0:6
    # Rodriguesの公式の展開
    @variables x
    D = n==0 ? x->x : Differential(x)^n # n階の微分演算子
    a = 1 // (2^n * factorial(n))       # 微分演算子の左側
    b = (x^2 - 1)^n                     # 微分演算子の右側
    c = a * D(b)                        # Rodriguesの公式
    d = expand_derivatives(c)           # 微分演算を展開
    e = simplify(d, expand=true)        # 簡略化
    # 閉形式
    P(x; n=0, m=0) = (1//2)^n * (1-x^2)^(m//2) * sum(j -> (-1)^j * factorial(2*n-2*j) // (factorial(j) * factorial(n-j) * factorial(n-2*j-m)) * x^(n-2*j-m), 0:Int(floor((n-m)/2)))
    # LaTeX表示
    display(latexstring(
    """
    \\begin{aligned}
        P_{$(n)}(x)
        = $(latexify(c, env=:raw))
        &= $(latexify(e, env=:raw)) \\\\
        &= $(latexify(P(x, n=n), env=:raw))
    \\end{aligned}
    """))
end
\begin{aligned} P_{0}(x) = 1 &= 1 \\ &= 1 \end{aligned}
\begin{aligned} P_{1}(x) = \frac{1}{2} \frac{\mathrm{d}}{\mathrm{d}x} \left( -1 + x^{2} \right) &= x \\ &= x \end{aligned}
\begin{aligned} P_{2}(x) = \frac{1}{8} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \left( -1 + x^{2} \right)^{2} &= \frac{-1}{2} + \frac{3}{2} x^{2} \\ &= \frac{-1}{2} + \frac{3}{2} x^{2} \end{aligned}
\begin{aligned} P_{3}(x) = \frac{1}{48} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \left( -1 + x^{2} \right)^{3} &= - \frac{3}{2} x + \frac{5}{2} x^{3} \\ &= - \frac{3}{2} x + \frac{5}{2} x^{3} \end{aligned}
\begin{aligned} P_{4}(x) = \frac{1}{384} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \left( -1 + x^{2} \right)^{4} &= \frac{3}{8} - \frac{15}{4} x^{2} + \frac{35}{8} x^{4} \\ &= \frac{3}{8} - \frac{15}{4} x^{2} + \frac{35}{8} x^{4} \end{aligned}
\begin{aligned} P_{5}(x) = \frac{1}{3840} \frac{\mathrm{d}^{5}}{\mathrm{d}x^{5}} \left( -1 + x^{2} \right)^{5} &= \frac{15}{8} x + \frac{63}{8} x^{5} - \frac{35}{4} x^{3} \\ &= \frac{15}{8} x + \frac{63}{8} x^{5} - \frac{35}{4} x^{3} \end{aligned}
\begin{aligned} P_{6}(x) = \frac{1}{46080} \frac{\mathrm{d}^{6}}{\mathrm{d}x^{6}} \left( -1 + x^{2} \right)^{6} &= \frac{-5}{16} + \frac{105}{16} x^{2} + \frac{231}{16} x^{6} - \frac{315}{16} x^{4} \\ &= \frac{-5}{16} + \frac{105}{16} x^{2} + \frac{231}{16} x^{6} - \frac{315}{16} x^{4} \end{aligned}

ルジャンドル陪多項式

Legendre陪多項式のRodriguesの公式(上段)を展開し, 閉形式(下段)と比較します. 展開した結果はMcQuarrie&Simon(1997)等と見比べてください.

\begin{aligned} P_n^m(x) &= \left( 1-x^2 \right)^{m/2} \frac{\mathrm{d}^m}{\mathrm{d}x^m} P_n(x) \\ &= \frac{1}{2^n} (1-x^2)^{m/2} \sum_{j=0}^{\left\lfloor\frac{n-m}{2}\right\rfloor} (-1)^j \frac{(2n-2j)!}{j! (n-j)! (n-2j-m)!} x^{(n-2j-m)} \end{aligned}
# using Symbolics
# using Latexify
# using LaTeXStrings

for n in 0:4
    for m in 0:n
        # Rodriguesの公式の展開
        @variables x
        Dn = n==0 ? x->x : Differential(x)^n # n階の微分演算子
        Dm = m==0 ? x->x : Differential(x)^m # m階の微分演算子
        a = 1 // (2^n * factorial(n))         # 微分演算子の左側
        b = (x^2 - 1)^n                       # 微分演算子の右側
        c = (1 - x^2)^(m//2) * Dm(a * Dn(b)) # Rodriguesの公式
        d = expand_derivatives(c)             # 微分演算を展開
        e = simplify(d, expand=true)          # 簡略化
        # 閉形式
        P(x; n=0, m=0) = (1//2)^n * (1-x^2)^(m//2) * sum(j -> (-1)^j * factorial(2*n-2*j) // (factorial(j) * factorial(n-j) * factorial(n-2*j-m)) * x^(n-2*j-m), 0:Int(floor((n-m)/2)))
        # LaTeX表示
        display(latexstring(
        """
        \\begin{aligned}
            P_{$n}^{$m}(x)
             = $(latexify(c, env=:raw))
            &= $(latexify(e, env=:raw)) \\\\
            &= $(latexify(simplify(P(x, n=n, m=m), expand=true), env=:raw))
        \\end{aligned}
        """))
    end
end
\begin{aligned} P_{0}^{0}(x) = 1 &= 1 \\ &= 1 \end{aligned}
\begin{aligned} P_{1}^{0}(x) = \frac{1}{2} \frac{\mathrm{d}}{\mathrm{d}x} \left( -1 + x^{2} \right) &= x \\ &= x \end{aligned}
\begin{aligned} P_{1}^{1}(x) = \left( 1 - x^{2} \right)^{\frac{1}{2}} \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{2} \frac{\mathrm{d}}{\mathrm{d}x} \left( -1 + x^{2} \right) &= \left( 1 - x^{2} \right)^{\frac{1}{2}} \\ &= \left( 1 - x^{2} \right)^{\frac{1}{2}} \end{aligned}
\begin{aligned} P_{2}^{0}(x) = \frac{1}{8} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \left( -1 + x^{2} \right)^{2} &= \frac{-1}{2} + \frac{3}{2} x^{2} \\ &= \frac{-1}{2} + \frac{3}{2} x^{2} \end{aligned}
\begin{aligned} P_{2}^{1}(x) = \left( 1 - x^{2} \right)^{\frac{1}{2}} \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{8} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \left( -1 + x^{2} \right)^{2} &= 3 \left( 1 - x^{2} \right)^{\frac{1}{2}} x \\ &= 3 \left( 1 - x^{2} \right)^{\frac{1}{2}} x \end{aligned}
\begin{aligned} P_{2}^{2}(x) = \left( 1 - x^{2} \right) \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \frac{1}{8} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \left( -1 + x^{2} \right)^{2} &= 3 - 3 x^{2} \\ &= 3 - 3 x^{2} \end{aligned}
\begin{aligned} P_{3}^{0}(x) = \frac{1}{48} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \left( -1 + x^{2} \right)^{3} &= - \frac{3}{2} x + \frac{5}{2} x^{3} \\ &= - \frac{3}{2} x + \frac{5}{2} x^{3} \end{aligned}
\begin{aligned} P_{3}^{1}(x) = \left( 1 - x^{2} \right)^{\frac{1}{2}} \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{48} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \left( -1 + x^{2} \right)^{3} &= - \frac{3}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} + \frac{15}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} x^{2} \\ &= - \frac{3}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} + \frac{15}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} x^{2} \end{aligned}
\begin{aligned} P_{3}^{2}(x) = \left( 1 - x^{2} \right) \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \frac{1}{48} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \left( -1 + x^{2} \right)^{3} &= 15 x - 15 x^{3} \\ &= 15 x - 15 x^{3} \end{aligned}
\begin{aligned} P_{3}^{3}(x) = \left( 1 - x^{2} \right)^{\frac{3}{2}} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \frac{1}{48} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \left( -1 + x^{2} \right)^{3} &= 15 \left( 1 - x^{2} \right)^{\frac{3}{2}} \\ &= 15 \left( 1 - x^{2} \right)^{\frac{3}{2}} \end{aligned}
\begin{aligned} P_{4}^{0}(x) = \frac{1}{384} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \left( -1 + x^{2} \right)^{4} &= \frac{3}{8} - \frac{15}{4} x^{2} + \frac{35}{8} x^{4} \\ &= \frac{3}{8} - \frac{15}{4} x^{2} + \frac{35}{8} x^{4} \end{aligned}
\begin{aligned} P_{4}^{1}(x) = \left( 1 - x^{2} \right)^{\frac{1}{2}} \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{384} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \left( -1 + x^{2} \right)^{4} &= - \frac{15}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} x + \frac{35}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} x^{3} \\ &= - \frac{15}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} x + \frac{35}{2} \left( 1 - x^{2} \right)^{\frac{1}{2}} x^{3} \end{aligned}
\begin{aligned} P_{4}^{2}(x) = \left( 1 - x^{2} \right) \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \frac{1}{384} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \left( -1 + x^{2} \right)^{4} &= \frac{-15}{2} + 60 x^{2} - \frac{105}{2} x^{4} \\ &= \frac{-15}{2} + 60 x^{2} - \frac{105}{2} x^{4} \end{aligned}
\begin{aligned} P_{4}^{3}(x) = \left( 1 - x^{2} \right)^{\frac{3}{2}} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \frac{1}{384} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \left( -1 + x^{2} \right)^{4} &= 105 \left( 1 - x^{2} \right)^{\frac{3}{2}} x \\ &= 105 \left( 1 - x^{2} \right)^{\frac{3}{2}} x \end{aligned}
\begin{aligned} P_{4}^{4}(x) = \left( 1 - x^{2} \right)^{2} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \frac{1}{384} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \left( -1 + x^{2} \right)^{4} &= 105 \left( 1 - x^{2} \right)^{2} \\ &= 105 \left( 1 - x^{2} \right)^{2} \end{aligned}

ラゲール多項式

Laguerre多項式のRodriguesの公式を展開し, 閉形式と比較します. 閉形式についてはLaguerre陪多項式のk=0および一般化Laguerre多項式の\alpha=0の場合を用います. 展開した結果はDLMF等と見比べてください.

\begin{aligned} L_n(x) &= \frac{\mathrm{e}^x}{n!} \frac{\mathrm{d}^n}{\mathrm{d}x ^n} \left( \mathrm{e}^{-x} x^n \right) \\ &= L_n^{0}(x) \\ &= L_n^{(0)}(x) \end{aligned}
# using Symbolics
# using Latexify
# using LaTeXStrings

for n in 0:6
    # Rodriguesの公式の展開
    @variables x
    D = n==0 ? x->x : Differential(x)^n # n階の微分演算子
    a = exp(x) / factorial(n)           # 微分演算子の左側
    b = exp(-x) * x^n                   # 微分演算子の右側
    c = a * D(b)                        # Rodriguesの公式
    d = expand_derivatives(c)           # 微分演算を展開
    e = simplify(d, expand=true)        # 簡略化
    # 閉形式
    AL(x; n=0, k=0)  = sum(m -> (-1)^(m+k) * factorial(n) // (factorial(m) * factorial(m+k) * factorial(n-m-k)) * x^m, 0:n-k)
    GL(x; n=0, α=0) = sum(k -> (-1)^(k) * binomial(n+α,n-k) * x^k // factorial(k), 0:n)
    # LaTeX表示
    display(latexstring(
    """
    \\begin{aligned}
        L_{$n}(x)
         = $(latexify(c, env=:raw))
        &= $(latexify(e, env=:raw)) \\\\
        &= $(latexify(AL(x, n=n, k=0), env=:raw)) \\\\
        &= $(latexify(GL(x, n=n, α=0), env=:raw))
    \\end{aligned}
    """))
end
\begin{aligned} L_{0}(x) = e^{x} e^{ - x} &= 1 \\ &= 1 \\ &= 1 \end{aligned}
\begin{aligned} L_{1}(x) = e^{x} \frac{\mathrm{d}}{\mathrm{d}x} x e^{ - x} &= 1 - x \\ &= 1 - x \\ &= 1 - x \end{aligned}
\begin{aligned} L_{2}(x) = \frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{2} e^{ - x} &= 1 + \frac{1}{2} x^{2} - 2 x \\ &= 1 + \frac{1}{2} x^{2} - 2 x \\ &= 1 + \frac{1}{2} x^{2} - 2 x \end{aligned}
\begin{aligned} L_{3}(x) = \frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{3} e^{ - x} &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \\ &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \\ &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \end{aligned}
\begin{aligned} L_{4}(x) = \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \\ &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \\ &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \end{aligned}
\begin{aligned} L_{5}(x) = \frac{1}{120} e^{x} \frac{\mathrm{d}^{5}}{\mathrm{d}x^{5}} x^{5} e^{ - x} &= 1 + \frac{5}{24} x^{4} - \frac{1}{120} x^{5} - \frac{5}{3} x^{3} - 5 x + 5 x^{2} \\ &= 1 + \frac{5}{24} x^{4} - \frac{1}{120} x^{5} - \frac{5}{3} x^{3} - 5 x + 5 x^{2} \\ &= 1 + \frac{5}{24} x^{4} - \frac{1}{120} x^{5} - \frac{5}{3} x^{3} - 5 x + 5 x^{2} \end{aligned}
\begin{aligned} L_{6}(x) = \frac{1}{720} e^{x} \frac{\mathrm{d}^{6}}{\mathrm{d}x^{6}} x^{6} e^{ - x} &= 1 + \frac{1}{720} x^{6} + \frac{5}{8} x^{4} - \frac{1}{20} x^{5} - \frac{10}{3} x^{3} - 6 x + \frac{15}{2} x^{2} \\ &= 1 + \frac{1}{720} x^{6} + \frac{5}{8} x^{4} - \frac{1}{20} x^{5} - \frac{10}{3} x^{3} - 6 x + \frac{15}{2} x^{2} \\ &= 1 + \frac{1}{720} x^{6} + \frac{5}{8} x^{4} - \frac{1}{20} x^{5} - \frac{10}{3} x^{3} - 6 x + \frac{15}{2} x^{2} \end{aligned}

ラゲール陪多項式

Laguerre陪多項式のRodriguesの公式(上段)を展開し, 閉形式(下段)と比較します. 閉形式については, いつも通りに関数を定義し, そこに@variables xで定義したSymbolics.jlの変数を渡せばSymbolics.jlのオブジェクトとして返ってきます. 展開した結果はWikipedia等と見比べてください.

\begin{aligned} L_n^{k}(x) &= \frac{\mathrm{d}^k}{\mathrm{d}x^k} L_n(x) \\ &= \sum_{m=0}^{n-k} (-1)^{m+k} \frac{n!}{m!(m+k)!(n-m-k)!} x^m \\ &= (-1)^k L_{n-k}^{(k)}(x) \end{aligned}
# using Symbolics
# using Latexify
# using LaTeXStrings

for n in 0:4
    for k in 0:n
        # Rodriguesの公式の展開
        @variables x
        Dn = n==0 ? x->x : Differential(x)^n # n階の微分演算子
        Dk = k==0 ? x->x : Differential(x)^k # k階の微分演算子
        a = exp(x) / factorial(n)             # 微分演算子の左側
        b = exp(-x) * x^n                     # 微分演算子の右側
        c = Dk(a * Dn(b))                     # Rodriguesの公式
        d = expand_derivatives(c)             # Rodriguesの公式の微分演算を展開
        e = simplify(d, expand=true)          # 簡略化
        # 閉形式
        AL(x; n=0, k=0)  = sum(m -> (-1)^(m+k) * factorial(n) // (factorial(m) * factorial(m+k) * factorial(n-m-k)) * x^m, 0:n-k)
        GL(x; n=0, α=0) = sum(k -> (-1)^(k) * binomial(n+α,n-k) * x^k // factorial(k), 0:n)
        # LaTeX表示
        display(latexstring(
        """
        \\begin{aligned}
            L_{$n}^{$k}(x)
             = $(latexify(c, env=:raw))
            &= $(latexify(e, env=:raw)) \\\\
            &= $(latexify(AL(x, n=n, k=k), env=:raw)) \\\\
            &= $(latexify((-1)^k * GL(x, n=n-k, α=k), env=:raw))            
        \\end{aligned}
        """))
    end
end
\begin{aligned} L_{0}^{0}(x) = e^{x} e^{ - x} &= 1 \\ &= 1 \\ &= 1 \end{aligned}
\begin{aligned} L_{1}^{0}(x) = e^{x} \frac{\mathrm{d}}{\mathrm{d}x} x e^{ - x} &= 1 - x \\ &= 1 - x \\ &= 1 - x \end{aligned}
\begin{aligned} L_{1}^{1}(x) = \frac{\mathrm{d}}{\mathrm{d}x} e^{x} \frac{\mathrm{d}}{\mathrm{d}x} x e^{ - x} &= -1 \\ &= -1 \\ &= -1 \end{aligned}
\begin{aligned} L_{2}^{0}(x) = \frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{2} e^{ - x} &= 1 + \frac{1}{2} x^{2} - 2 x \\ &= 1 + \frac{1}{2} x^{2} - 2 x \\ &= 1 + \frac{1}{2} x^{2} - 2 x \end{aligned}
\begin{aligned} L_{2}^{1}(x) = \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{2} e^{ - x} &= -2 + x \\ &= -2 + x \\ &= -2 + x \end{aligned}
\begin{aligned} L_{2}^{2}(x) = \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{2} e^{ - x} &= 1 \\ &= 1 \\ &= 1 \end{aligned}
\begin{aligned} L_{3}^{0}(x) = \frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{3} e^{ - x} &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \\ &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \\ &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \end{aligned}
\begin{aligned} L_{3}^{1}(x) = \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{3} e^{ - x} &= -3 - \frac{1}{2} x^{2} + 3 x \\ &= -3 - \frac{1}{2} x^{2} + 3 x \\ &= -3 - \frac{1}{2} x^{2} + 3 x \end{aligned}
\begin{aligned} L_{3}^{2}(x) = \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{3} e^{ - x} &= 3 - x \\ &= 3 - x \\ &= 3 - x \end{aligned}
\begin{aligned} L_{3}^{3}(x) = \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{3} e^{ - x} &= -1 \\ &= -1 \\ &= -1 \end{aligned}
\begin{aligned} L_{4}^{0}(x) = \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \\ &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \\ &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \end{aligned}
\begin{aligned} L_{4}^{1}(x) = \frac{\mathrm{d}}{\mathrm{d}x} \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= -4 - 2 x^{2} + \frac{1}{6} x^{3} + 6 x \\ &= -4 - 2 x^{2} + \frac{1}{6} x^{3} + 6 x \\ &= -4 - 2 x^{2} + \frac{1}{6} x^{3} + 6 x \end{aligned}
\begin{aligned} L_{4}^{2}(x) = \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= 6 + \frac{1}{2} x^{2} - 4 x \\ &= 6 + \frac{1}{2} x^{2} - 4 x \\ &= 6 + \frac{1}{2} x^{2} - 4 x \end{aligned}
\begin{aligned} L_{4}^{3}(x) = \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= -4 + x \\ &= -4 + x \\ &= -4 + x \end{aligned}
\begin{aligned} L_{4}^{4}(x) = \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= 1 \\ &= 1 \\ &= 1 \end{aligned}

一般化ラゲール多項式

一般化Laguerre多項式のRodriguesの公式(上段)を展開し, 閉形式(下段)と比較します.

\begin{aligned} L_n^{(\alpha)}(x) &= \frac{x^{-\alpha}e^x}{n!} \frac{d^n}{dx^n}\left(x^{n+\alpha}e^{-x}\right) \\ &= \sum_{k=0}^n(-1)^k \left(\begin{array}{l} n+\alpha \\ n-k \end{array}\right) \frac{x^k}{k !} \\ &= \sum_{k=0}^n(-1)^k \frac{\Gamma(\alpha+n+1)}{\Gamma(\alpha+k+1)\Gamma(n-k+1)} \frac{x^k}{k !} \end{aligned}

二項係数\left(\begin{array}{l} a \\ b \end{array}\right)binomial(a,b)で計算できます. ただし, \alphaが整数の場合以外はbinomial(a,b)は計算できないので, \left(\begin{array}{l} a \\ b \end{array}\right) = \frac{a !}{b !(a-b) !} = \frac{\Gamma(a+1)}{\Gamma(b+1)\Gamma(a-b+1)} を利用し, 階乗表示を経由してガンマ関数から計算することにします.(そもそも, \alphaが整数なら陪関数との対応関係を使えば十分です. Morse振動子の波動関数では\alphaは非整数なので, 非整数に対応しておく必要があります.) なお, Laguerre多項式の定義によっては, 一般化Laguerre多項式はLaguerre多項式の一般化にならないことに注意してください. \alpha=0の場合はL_n^{(0)}(x) = \frac{e^x}{n!} \frac{d^n}{dx^n}\left(x^{n}e^{-x}\right)となり, Laguerre多項式の定義をL_n(x) = e^x \frac{d^n}{dx^n}\left(x^{n}e^{-x}\right)とした場合とは一致しません. このノートではLaguerre多項式をL_n(x) = \frac{e^x}{n!} \frac{d^n}{dx^n}\left(x^{n}e^{-x}\right)と定義したので, 一般化になっています.

# using Symbolics
# using Latexify
# using LaTeXStrings
# using SpecialFunctions

for n in 0:4
    for α in 0:n
        # Rodriguesの公式の展開
        @variables x
        D = n==0 ? x->x : Differential(x)^n # n階の微分演算子
        a = exp(x) * x^(-α) / factorial(n) # 微分演算子の左側
        b = exp(-x) * x^(n+α)              # 微分演算子の右側
        c = a * D(b)                        # Rodriguesの公式
        d = expand_derivatives(c)           # 微分演算を展開
        e = simplify(d, expand=true)        # 簡略化
        # 閉形式
        GL1(x; n=0, α=0) = sum(k -> (-1)^(k) * binomial(n+α,n-k) * x^k // factorial(k), 0:n)
        GL2(x; n=0, α=0) = sum(k -> (-1)^(k) * (Int(gamma(α+n+1)) // Int((gamma(α+1+k)*gamma(n-k+1)))) * x^k // factorial(k), 0:n) # 綺麗に表示されるように整数に変換している
        # GL2(x; n=0, α=0) = sum(k -> (-1)^(k) * (gamma(α+n+1) // (gamma(α+1+k)*gamma(n-k+1))) * x^k // factorial(k), 0:n)
        # LaTeX表示
        display(latexstring(
        """
        \\begin{aligned}
            L_{$n}^{($α)}(x)
             = $(latexify(c, env=:raw))
            &= $(latexify(e, env=:raw)) \\\\
            &= $(latexify(GL1(x, n=n, α=α), env=:raw)) \\\\
            &= $(latexify(GL2(x, n=n, α=α), env=:raw))
        \\end{aligned}
        """))
    end
end
\begin{aligned} L_{0}^{(0)}(x) = e^{x} e^{ - x} &= 1 \\ &= 1 \\ &= 1 \end{aligned}
\begin{aligned} L_{1}^{(0)}(x) = e^{x} \frac{\mathrm{d}}{\mathrm{d}x} x e^{ - x} &= 1 - x \\ &= 1 - x \\ &= 1 - x \end{aligned}
\begin{aligned} L_{1}^{(1)}(x) = \frac{e^{x} \frac{\mathrm{d}}{\mathrm{d}x} x^{2} e^{ - x}}{x} &= 2 - x \\ &= 2 - x \\ &= 2 - x \end{aligned}
\begin{aligned} L_{2}^{(0)}(x) = \frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{2} e^{ - x} &= 1 + \frac{1}{2} x^{2} - 2 x \\ &= 1 + \frac{1}{2} x^{2} - 2 x \\ &= 1 + \frac{1}{2} x^{2} - 2 x \end{aligned}
\begin{aligned} L_{2}^{(1)}(x) = \frac{\frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{3} e^{ - x}}{x} &= 3 + \frac{1}{2} x^{2} - 3 x \\ &= 3 + \frac{1}{2} x^{2} - 3 x \\ &= 3 + \frac{1}{2} x^{2} - 3 x \end{aligned}
\begin{aligned} L_{2}^{(2)}(x) = \frac{\frac{1}{2} e^{x} \frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}} x^{4} e^{ - x}}{x^{2}} &= 6 + \frac{1}{2} x^{2} - 4 x \\ &= 6 + \frac{1}{2} x^{2} - 4 x \\ &= 6 + \frac{1}{2} x^{2} - 4 x \end{aligned}
\begin{aligned} L_{3}^{(0)}(x) = \frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{3} e^{ - x} &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \\ &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \\ &= 1 + \frac{3}{2} x^{2} - \frac{1}{6} x^{3} - 3 x \end{aligned}
\begin{aligned} L_{3}^{(1)}(x) = \frac{\frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{4} e^{ - x}}{x} &= 4 + 2 x^{2} - \frac{1}{6} x^{3} - 6 x \\ &= 4 + 2 x^{2} - \frac{1}{6} x^{3} - 6 x \\ &= 4 + 2 x^{2} - \frac{1}{6} x^{3} - 6 x \end{aligned}
\begin{aligned} L_{3}^{(2)}(x) = \frac{\frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{5} e^{ - x}}{x^{2}} &= 10 + \frac{5}{2} x^{2} - \frac{1}{6} x^{3} - 10 x \\ &= 10 + \frac{5}{2} x^{2} - \frac{1}{6} x^{3} - 10 x \\ &= 10 + \frac{5}{2} x^{2} - \frac{1}{6} x^{3} - 10 x \end{aligned}
\begin{aligned} L_{3}^{(3)}(x) = \frac{\frac{1}{6} e^{x} \frac{\mathrm{d}^{3}}{\mathrm{d}x^{3}} x^{6} e^{ - x}}{x^{3}} &= 20 + 3 x^{2} - \frac{1}{6} x^{3} - 15 x \\ &= 20 + 3 x^{2} - \frac{1}{6} x^{3} - 15 x \\ &= 20 + 3 x^{2} - \frac{1}{6} x^{3} - 15 x \end{aligned}
\begin{aligned} L_{4}^{(0)}(x) = \frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{4} e^{ - x} &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \\ &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \\ &= 1 + 3 x^{2} + \frac{1}{24} x^{4} - \frac{2}{3} x^{3} - 4 x \end{aligned}
\begin{aligned} L_{4}^{(1)}(x) = \frac{\frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{5} e^{ - x}}{x} &= 5 + 5 x^{2} + \frac{1}{24} x^{4} - \frac{5}{6} x^{3} - 10 x \\ &= 5 + 5 x^{2} + \frac{1}{24} x^{4} - \frac{5}{6} x^{3} - 10 x \\ &= 5 + 5 x^{2} + \frac{1}{24} x^{4} - \frac{5}{6} x^{3} - 10 x \end{aligned}
\begin{aligned} L_{4}^{(2)}(x) = \frac{\frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{6} e^{ - x}}{x^{2}} &= 15 + \frac{15}{2} x^{2} + \frac{1}{24} x^{4} - x^{3} - 20 x \\ &= 15 + \frac{15}{2} x^{2} + \frac{1}{24} x^{4} - x^{3} - 20 x \\ &= 15 + \frac{15}{2} x^{2} + \frac{1}{24} x^{4} - x^{3} - 20 x \end{aligned}
\begin{aligned} L_{4}^{(3)}(x) = \frac{\frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{7} e^{ - x}}{x^{3}} &= 35 + \frac{21}{2} x^{2} + \frac{1}{24} x^{4} - \frac{7}{6} x^{3} - 35 x \\ &= 35 + \frac{21}{2} x^{2} + \frac{1}{24} x^{4} - \frac{7}{6} x^{3} - 35 x \\ &= 35 + \frac{21}{2} x^{2} + \frac{1}{24} x^{4} - \frac{7}{6} x^{3} - 35 x \end{aligned}
\begin{aligned} L_{4}^{(4)}(x) = \frac{\frac{1}{24} e^{x} \frac{\mathrm{d}^{4}}{\mathrm{d}x^{4}} x^{8} e^{ - x}}{x^{4}} &= 70 + 14 x^{2} + \frac{1}{24} x^{4} - \frac{4}{3} x^{3} - 56 x \\ &= 70 + 14 x^{2} + \frac{1}{24} x^{4} - \frac{4}{3} x^{3} - 56 x \\ &= 70 + 14 x^{2} + \frac{1}{24} x^{4} - \frac{4}{3} x^{3} - 56 x \end{aligned}

まとめ

おかげで各種の閉形式についてRodriguesの公式との整合性が簡単に確認できるようになりました. Symbolics.jlはまだ発展途中のパッケージですが, かなり活躍してくれましたね. 自分でルールを追加して機能を拡張できるようなので, これから比較にならないほど発展すると思われます. 現状のシステムに強いて注文を付けるとすれば, 多項式の各項を次数の順番でソートできると良いですね. Issueに上がっているので, 対応されるものと思われます.

参考文献

https://symbolics.juliasymbolics.org/dev/manual/derivatives/

文献
H_{n}(x)
P_n(x)
P_k^m(x)
L_n(x)
L_n^{k}(x)
L_n^{(\alpha)}(x)
cpprefjp hermite legendre assoc_legendre laguerre assoc_laguerre
The Digital Library of Mathematical Functions (DLMF) 18.3 Table1, 18.5 Table1, 18.5.13, 18.5.18 18.3 Table1, 18.5 Table1, 18.5.16 18.3 Table1, 18.5 Table1, 18.5.17 18.3 Table1, 18.5 Table1, 18.5.12
L. D. Landau, E. M. Lifshitz, Quantum Mechanics (Pergamon Press, 1965) p.595 (a.4), (a.6) p.598 (c.1) p.598 (c.4) p.603 (d.13) p.603 (d.13)
L. I. Schiff, Quantum Mechanics (McGraw-Hill Book Company, 1968) p.71 (13.12) p.79 (14.12) p.93 (16.19)
A. Messiah, Quanfum Mechanics (Dover Publications, 1999) p.491 (B.59) p.493 (B.72),p.494 Table p.493 (B.72) p.483 (B.12) p.483 (B.12)
W. Greiner, Quantum Mechanics: An Introduction Third Edition (Springer, 1994) p.152 (7.22) p.83 (4) p.83 (5) p.149 (21)
D. J. Griffiths, Introduction to Quantum Mechanics (Prentice Hall, 1995) p.41 Table 2.1, p.43 (2.70) p.126 (4.28), p.96 Table3.1 p.126 (4.27) p.139 (4.88), p.140 Table4.4 p.139 (4.87), p.140 Table4.5
D. A. McQuarrie, J. D. Simon, Physical Chemistry: A Molecular Approach (University Science Books, 1997) p.170 Table 5.2 p.195 Table6.1 p.196 (6.26), p.196 Table6.2 p.207 Table6.4
P. W. Atkins, J. De Paula, Atkins' Physical Chemistry, 8th edition (W. H. Freeman, 2008) p.293 Table 9.1
J. J. Sakurai, J. Napolitano, Modern Quantum Mechanics Third Edition (Cambridge University Press, 2021) p.524 (B.29) p.245 Problem 3.30.b
J. P. Dahl, M. Springborg, J. Chem. Phys. 88, 4535 (1988) (62), (63)
W. K. Shao, Y. He, J. Pan, J. Nonlinear Sci. Appl., 9, 5, 3388 (2016) (1.6)

https://gist.github.com/ohno/287799f1f591fb3e6541c7ecee69cd03

Discussion

ログインするとコメントできます