Cython化の第1歩は、julia_py.py を julia_cy0.pyxに変更し、コンパイルしてみて、どれだけ高速化するかを調べる。
中身は何も変わらないので、プログラムは省略。
In [1]: import julia_cy0
In [2]: julia_cy0.julia()
time=14.09321284294128425.32秒から14.09秒に改善された。何も努力しなくて1.8倍速になた。
では、次々に型宣言をすることで、高速化していこう。
まずは、一番内側の、とても小さい関数（複素数zの絶対値の2乗を返す）に手をつける。
この関数は、Cythonから呼び出されるだけなので、defではなくcdefで関数の宣言をする。
型は、戻り値をdouble、引数をdouble complexとした。
この関数norm2は非常に短く、これを呼び出すのも一箇所だけである。そして、この関数は頻繁に呼び出されるので、呼び出しのオーバーヘッドが高いと思われる。そのため、inlineを加えて、呼び出し元のnorm2()を呼び出す部分を、コンパイル時にnorm2()の中身に置き換えるように指示した。
cdef inline double norm2(double complex z):
return z.real * z.real + z.imag * z.imagこの、たった1行だけの変更した後に実行してみた。
In [1]: import julia_cy1
In [2]: julia_cy1.julia()
time=5.78563666343689さらに高速化され、5.78秒に縮まり、25.32/5.78=4.38倍の速度になった。
次は、数列の発散までの項数を求めるescape()関数について型宣言する。
cdef int escape( double complex z, double complex c,
double z_max, int n_max):
cdef:
int i = 0
double z_max2 = z_max * z_max
while norm2(z) < z_max2 and i < n_max:
z = z * z + c
i += 1
return iIn [1]: import julia_cy2
In [2]: julia_cy2.julia()
time=0.5478911399841309ちょっと型宣言を加えただけだが、高速化の効果は大きく、実行時間がこれだけで1/10になった。
全体では、25.32/0.547=46.2倍とかなり高速になった。
1秒を切ったので、もう待たされる感覚は無いかもしれないが、まだ100倍速には達していないので、まだ高速化を進めていこう。
def calc_julia( int resolution, double complex c,
double bound=1.2, double z_max=2.0, int n_max=1000 ):
cdef:
double step = 2.0 * bound / resolution
double complex z
double real, imag
int[:,:] counts
counts = np.zeros((resolution+1,resolution+1), dtype=np.int32)
for i in range(resolution+1):
imag = bound - i * step
for j in range(resolution+1):
real = -bound + j * step
z = real + imag *1j
counts[i,j] = escape(z,c,z_max,n_max)
return np.asarray(counts)今回は、calc_juliaの中で型宣言を行った。
countsは、numpyの配列のままになっていたが、counts配列を２次元配列として宣言することによって、効率UPを図った。
n [1]: import julia_cy3
In [2]: julia_cy3.julia()
time=0.347775936126709実行結果は、さらに向上し、0.347秒になった。
全体では、25.32/0.347=72.9倍とかなり高速になった。順調である。
countsという配列を使っているのだが、こういう場合は、境界チェックとラップアラウンドチェックが行われて遅くなる。配列の要素指定のインデックス値が範囲をハズレないし、負の値を使うこともないので、境界チェックとラップアラウンドチェックを無効にできると、オライリーのCythonの本に出ていた。
それで、つぎの2行をプログラムの先頭に加えた。
# cython: boundscheck=False
# cython: wraparound=Falseあるいは、次の組み合わせでも同じである。
from cython cimport boundscheck, wraparound
@boundscheck(False)
@wraparound(False)これで、このCythonのファイルでは、境界チェックとラップアラウンドチェックが無効になり速くなる。
ということで、コンパイルし、実行してみた。
In [1]: import julia_cy４
In [2]: julia_cy4.julia()
time=0.32135868072509766若干速くなって、0.321秒であった。
全体では、25.32/0.321=78.8倍とかなり高速になった。順調である。
配列の添字範囲チェックが無効になるのは、escape関数の結果代入時だけなので、頻度はかなり少ないので、効果は限られるようだ。
型を指定するだけでは、高速化は、もう限界のような気がする。
残るは、最後の切り札、並列処理、マルチスレッドであろうが、説明は次回にする。
＃人気の記事
＃タグ
