Gamma関数をお絵描き

本日もpythonで描いたグラフを乗っけていこうと思います。第二回はガンマ関数のグラフについてです。

Gamma関数とは。

まずガンマ関数の定義ですが、これもWikipedia（ガンマ関数 - Wikipedia）を参照されるのがベターでしょうが、一応一般的によくある定義と解析接続したものであるWeierstrassの表示を紹介したいと思います。

Gamma関数 を満たすに対して $$\Gamma (z) = \int_0^{\infty} t^{z-1} e^{-t}dt $$

これも定義域が複素平面全てではないのでグラフを描くにあたっては都合がよくありません。そこで次の解析接続後の、Weierstrassの無限積表示の定義も書いておきます。

Gamma関数（Weierstrassの無限積表示） 原点と負の整数における極を除く任意のに対して $$\Gamma (z) = \bigg[ z e^{\gamma z} \prod _{r=1}^{\infty} \bigg( 1 + \frac{z}{r} \bigg) e^{-z/r} \bigg]^{-1}$$ ただしはオイラーの定数である。

この表示の逆数、即ちは任意のにおいて正則（複素微分可能で導関数が連続）です。他にもガンマ関数に関する定理、命題は次のページによくまとまっています。

Gamma Function -- from Wolfram MathWorld

pythonを使って描く。

まずに任意の実数を代入してガンマ関数を返す時のグラフが次のコードである。

from scipy.special import gamma, zeta
from matplotlib import pyplot as plt
%matplotlib inline
import numpy as np

x = np.linspace( -5., 5., 500)
y = gamma(x)
plt.axis([-5., 5., -5., 5.])
plt.plot(x, y,"o")
plt.savefig("gamma_real.png")
plt.show()

すると得られる画像が次の通り。

さて次に絶対値付きガンマ関数の3Dのグラフのコードが次の通りです。

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from mpl_toolkits.mplot3d import Axes3D
from scipy.special import gamma

def three_D_plot(func,picname):
    x, y = np.meshgrid(np.linspace( -5, 5, 500), np.linspace( -5, 5, 500))
    z = x + y*1j
    w = func(z)
    fig = plt.figure()
    ax = Axes3D(fig)
    ax.plot_surface(x, y, np.abs(w))
    ax.set_zlim3d(0, 15.)
    ax.autoscale(False)
    plt.axis([-6., 6., -6., 6.])
    plt.savefig(picname)
    plt.show()

three_D_plot(gamma,picname="gamma3D.png")

するとグラフは次のようになります。

この図からも分かるように虚部がゼロの部分の断面図が先の図で実軸対象に折り返したものに一致することが確認できます。次にこのグラフの偏角の様子とレベルカーブを複素平面にプロットした図が次のコードが与える2つのグラフです。1つ目のカラフルな図が偏角によって色が変化していて、2つ目がレベルカーブです。

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from scipy.special import gamma

def plot_complex(complex_func,picname1,picname2):

    x, y = np.meshgrid(np.linspace(-4, 4, 400), np.linspace(-4, 4, 400))
    z = x + y*1j

    angles = (np.angle(complex_func(z)) + 2*np.pi) % (2*np.pi)
    plt.figure(figsize=(6,5))
    im_1 = plt.pcolor(x, y, angles, cmap="hsv", vmin=0, vmax=6)
    plt.colorbar(im_1)
    plt.title('argument')
    plt.savefig(picname1)

    plt.figure(figsize=(6,5))
    im_2 = plt.pcolor(x, y, np.abs(complex_func(z)), vmin=0, vmax=4)
    plt.title('absolute value')
    plt.colorbar(im_2)
    plt.savefig(picname2)

    plt.show()

plot_complex(gamma,"gamma_arg.png","gamma_abs.png")

上の図から実部が負の整数で虚部がゼロの時、黄色っぽくなっていることからここで極を取っていることが分かりますね。さて今まで書いてきたガンマ関数ですが今回はpythonのscipyというライブラリにもともとあるガンマ関数を使っています。もちろんmpmathにもガンマ関数はあって前回も書いた通り、cplotという便利なプロット用の関数も用意されてるのでそちらを使っても同じものが描けます。最後にそれを一応乗せておきましょう。

from mpmath import *
cplot(gamma, points=50000, file="gamma2D.png")

グラフとしては上の自前のものをミックスしたような形になっているのが得られます。ちなみに黒の濃淡が絶対値の大小を表しています。