ガンマ関数の計算

ガンマ関数Γ(z)は Γ(z)=∫[t=0,∞]tz-1exp(-t)dt で定義される. Re(z)が正でない整数を除く複素数 z で定義されている. ガンマ関数は, 次のような性質を持っている.

  • Γ(z+1)=zΓ(z). 特に, 負でない整数 n について Γ(n+1)=n!
  • Γ(z)Γ(1-z)=π csc(πz)

特に, 前者の性質については階乗の計算を計算機上で効率よく行う方法として知られている. 計算機で十分に大きな n については ln(Γ(n)) を直接計算する方法が多数知られているため, 投げると表よりも裏の方が2倍出やすいコインを 1000 枚同時に投げてちょうど 300 枚が表になる確率を計算機で求めるときに, 定義どおり 1000C300(1/3)300(1-1/3)1000-300を計算するよりも, ln(Γ(1000+1))-ln(Γ(300+1))-ln(Γ(700+1))+300ln(1/3)+700ln(2/3) を計算して, それを指数とし, 底をネイピア数とする指数関数を計算したほうが,
中間値に大きな数が出てこないので, オーバーフローを気にしなくてもよい.

なお, 十分に大きな z に対して, ln(Γ(z))≒(ln(2π)-ln(z))/2+z(ln(z+1/(12z-1/(10z)))-1) なる近似がよく知られている (この近似による誤差は大よそ O((1/z)3/2).

一方, 絶対値が小さな z に先述の近似式を直接適用できない. その場合は二重指数変換などを利用した数値積分を利用することになる. Γ(z)の定義式を二重指数変換する (t=exp(u-exp(-u)) ) と, ∫[u=-∞,∞]exp(u-exp(u))z-1exp(u-exp(-u))(1+exp(-u))exp(-exp(u-exp(-u)))du となる. 被積分関数は u が十分に大きいと急速に 0 に正の方向から近づくため, この積分を台形則で計算すると高速に真値へ収束する.

コメントを残す