倭算数理研究所

科学・数学・学習関連の記事を、「倭マン日記」とは別に書いていくのだ！

トップ > Scala > spire で多項式 Polynomial を使ってみる (5) ： 多項式で表される特殊関数

spire で多項式 Polynomial を使ってみる (5) ： 多項式で表される特殊関数

Scala 特殊関数

spire を使ってみるシリーズ（目次）。　今回は多項式で表される特殊関数に関連する spire.math.poly.SpecialPolynomials クラスのメソッドを見ていきます。

SpecialPolynomials クラスでは

  • ルジャンドルの多項式： legendre()
  • ラゲールの多項式： laguerre()
  • チェビシェフの多項式
    • 第1種： chebyshevsFirstKind()
    • 第2種： chebyshevsSecondKind()
  • エルミートの多項式
    • 物理学で使うバージョン： physHermite()
    • 確率論で使うバージョン： probHermite()

の6つが使えます (spire 0.10)。　使い方はどれも同じで、上記のメソッドに Int 値を与えて、長さがその Int 値の Stream を返します。　この Stream の n 番目の要素は各特殊関数のパラメータの値が n となっています。　なぜ Stream を返すのに n 個で必ず切らないといけなくしているのかよく分かりませんが*1、まぁ使い方は別に難しくはないと思います。

各特殊関数の箇所に漸化式を書いているものは spire の実装で使われていることのメモでもあります (spire 0.10)。

以下のサンプルコードでは下記の import 文が書かれているものとします：

import spire.math._
import spire.math.poly.SpecialPolynomials

【参考】

ルジャンドルの多項式

ルジャンドル (Legendre) の多項式 { P_n(x) } は

　　{ \displaystyle\begin{align*} P_n(x) &= \frac{1}{2^n n!}\frac{d^n}{dx^n}(x^2 - 1)^n \end{align*}}

で与えられ、以下を満たします：

　　{ \displaystyle\begin{align*} P_0(x) &= 1 \\ P_1(x) &= x \\ P_{n+1}(x) &= \frac{1}{n+1}\Big\{(2n+1) xP_n(x) - nP_{n-1}(x)\Big\} \end{align*}}

spire では SpecialPolynomials.legendre() メソッドで取得できます：

  // Legendre Polynomials
  val ps: Stream[Polynomial[Rational]] = SpecialPolynomials.legendres(5)

  assert( ps(0) == Polynomial("1") )
  assert( ps(1) == Polynomial("x") )
  assert( ps(2) == Polynomial("3/2x^2 - 1/2") )
  assert( ps(3) == Polynomial("5/2x^3 - 3/2x") )
  assert( ps(4) == Polynomial("35/8x^4 - 15/4x^2 + 3/8") )

ラゲールの多項式

ラゲール (Laguerre) の多項式 { L_n(x) } は

　　{ \displaystyle\begin{align*} L_n(x) &= \sum_{r=0}^n(-1)^r\binom{n}{r}\frac{x^r}{r!} \\[2mm] &= \frac{e^x}{n!}\frac{d^n}{dx^n}\left(e^{-x}x^n\right) = \frac{1}{n!}\left(\frac{d}{dx} - 1\right)^nx^n \end{align*}}

で定義され、以下の漸化式を満たします：

　　{ \displaystyle\begin{align*} L_{n+1}(x) = \frac{1}{n+1}\Big\{(2n+1-x)L_n(x) - nL_{n-1}(x)\Big\} \end{align*}}

spire では SpecialPolynomial.laguerre() メソッドで取得できます：

  // Laguerre Polynomials
  val ls: Stream[Polynomial[Rational]] = SpecialPolynomials.laguerres(5)

  assert( ls(0) == Polynomial("1") )
  assert( ls(1) == Polynomial("-x + 1") )
  assert( ls(2) == Polynomial("1/2x^2 - 2x + 1") )
  assert( ls(3) == Polynomial("-1/6x^3 + 3/2x^2 - 3x + 1") )
  assert( ls(4) == Polynomial("1/24x^4 - 2/3x^3 + 3x^2 - 4x + 1") )

チェビシェフの多項式

第1種チェビシェフの多項式 { T_n(x) } は余弦の { n } 倍角の公式が

　　{ \displaystyle\begin{align*} \cos n\theta = T_n(\cos\theta) \end{align*}}

と表されるように定義されます（余弦の { n } 倍角の公式は『三角関数の n 倍角の公式を導く ～余弦編～』参照）。　 { T_n(x) } は以下を満たします：

　　{ \displaystyle\begin{align*} T_0(x) &= 1 \\ T_1(x) &= x \\ T_{n+1}(x) &= 2x T_n(x) - T_{n-1}(x) \end{align*}}

spire では SpecialPolynomial.chebyshevsFirstKind() メソッドで取得できます：

  // Chebyshev Polynomials of the First Kind
  val ts: Stream[Polynomial[Rational]] = SpecialPolynomials.chebyshevsFirstKind(5)

  assert( ts(0) == Polynomial("1") )
  assert( ts(1) == Polynomial("x") )
  assert( ts(2) == Polynomial("2x^2 - 1") )
  assert( ts(3) == Polynomial("4x^3 - 3x") )
  assert( ts(4) == Polynomial("8x^4 - 8x^2 + 1") )

第2種チェビシェフの多項式 { U_n(x) } は正弦の { n+1 } 倍角の公式が

　　{ \displaystyle\begin{align*} \frac{\sin (n+1)\theta}{\sin\theta} = U_n(\cos \theta) \end{align*}}

と表されるように定義されます。　 { U_n(x) } は以下を満たします：

　　{ \displaystyle\begin{align*} U_0(x) &= 1 \\ U_1(x) &= 2x \\ U_{n+1}(x) &= 2x U_n(x) - U_{n-1}(x) \end{align*}}

これらは { n = 1 } の場合だけが { T_n(x) = x } と異なるだけで、特に漸化式は同じです。

spire では SpecialPolynomial.chebyshevsSecondKind() メソッドで取得できます：

  // Chebyshev Polynomials of the Second Kind
  val us: Stream[Polynomial[Rational]] = SpecialPolynomials.chebyshevsSecondKind(5)

  assert( us(0) == Polynomial("1") )
  assert( us(1) == Polynomial("2x") )
  assert( us(2) == Polynomial("4x^2 - 1") )
  assert( us(3) == Polynomial("8x^3 - 4x") )
  assert( us(4) == Polynomial("16x^4 - 12x^2 + 1") )

エルミートの多項式

エルミート (Hermite) の多項式には分野によって2種類の定義があります。　物理学で使われるものは

　　{ \displaystyle\begin{align*} H_n(x) = (-1)^ne^{x^2}\frac{d^n}{dx^n}e^{-x^2} = \left(2x - \frac{d}{dx}\right)^n1 \end{align*}}

で定義されます。　漸化式の形で書けば以下のように書けます：

　　{ \displaystyle\begin{align*} H_0(x) &= 1 \\ H_1(x) &= 2x \\ H_{n+1}(x) &= \left(2x - \frac{d}{dx}\right)H_n(x) \end{align*}}

エルミートの多項式がエルミートの微分方程式の解として求まることは『エルミートの微分方程式とその級数解』で示しました。

spire では SpecialPolynomials.physHermites() メソッドで取得できます：

  // Hermite Polynomials (physics)
  val hs: Stream[Polynomial[Rational]] = SpecialPolynomials.physHermites(5)

  assert( hs(0) == Polynomial("1") )
  assert( hs(1) == Polynomial("2x") )
  assert( hs(2) == Polynomial("4x^2 - 2") )
  assert( hs(3) == Polynomial("8x^3 - 12x") )
  assert( hs(4) == Polynomial("16x^4 - 48x^2 + 12") )

確率論で使われるものは

　　{ \displaystyle\begin{align*} He_n(x) = (-1)^ne^\frac{x^2}{2}\frac{d^n}{dx^n}e^\frac{-x^2}{2} = \left(x - \frac{d}{dx}\right)^n1 \end{align*}}

で定義され、漸化式の形で以下のようにも書けます：

　　{ \displaystyle\begin{align*} He_0(x) &= 1 \\ He_1(x) &= x \\ He_{n+1}(x) &= \left(x - \frac{d}{dx}\right)He_n(x) \end{align*}}

spire では SpecialPolynomials.probHermites() メソッドで取得できます：

  // Hermite Polynomials (prob)
  val hes: Stream[Polynomial[Rational]] = SpecialPolynomials.probHermites(5)

  assert( hes(0) == Polynomial("1") )
  assert( hes(1) == Polynomial("x") )
  assert( hes(2) == Polynomial("x^2 - 1") )
  assert( hes(3) == Polynomial("x^3 - 3x") )
  assert( hes(4) == Polynomial("x^4 - 6x^2 + 3") )

これで多項式関連の API は大体見てきました。　spire.math.poly パッケージには他にもいくつかクラスなどが定義されていますが、直接使う必要はなさそうなもの（根を求めるなど）なので、このシリーズはこの辺で終了。

Scalaスケーラブルプログラミング第3版

Scalaスケーラブルプログラミング第3版

*1:ソースコードを見ると、漸化式で無限 Stream を作った後に単に take() を読んでるだけですが・・・