Cythonで高速化#
Pythonで配列の要素に対して大量のループ計算を行う場合、Pythonの実行効率はC言語に比べて数百倍遅くなることがあります。まず、配列演算関数の例を通じて、Cythonがどのように配列の高速演算を実現するかを示します。
ベクトル集合の距離行列の計算#
一組のベクトルの中で、各ペアのベクトル間の距離を計算する関数を実装します。以下の配列X
の形状は(200, 3)
で、200個の3次元空間中の点と見なすことができます:
import numpy as np
np.random.seed(42)
X = np.random.rand(200, 3)
各点間の距離を計算するために、NumPyのブロードキャスト機能を使用するか、SciPyのscipy.spatial.distance.pdist()
を使用することができます。ただし、ここではPythonとCythonの性能を比較するために、三重ループを使用して各要素を個別に計算します:
def pairwise_dist_python(X):
m, n = X.shape
D = np.empty((m, m), dtype=np.float64)
for i in range(m):
for j in range(i, m):
d = 0.0
for k in range(n):
tmp = X[i, k] - X[j, k]
d += tmp * tmp
D[i, j] = D[j, i] = d**0.5
return D
プログラムの結果が正しいことを確認するために、SciPyのpdist()
を使用して同じ計算を行い、両者の計算速度を比較します。両者の計算速度は約300倍異なることがわかります。
from scipy.spatial.distance import pdist, squareform
%timeit squareform(pdist(X))
%timeit pairwise_dist_python(X)
np.allclose(squareform(pdist(X)), pairwise_dist_python(X))
204 μs ± 8.66 μs per loop (mean ± std. dev. of 7 runs, 10,000 loops each)
51 ms ± 834 μs per loop (mean ± std. dev. of 7 runs, 10 loops each)
True
%%cythonマジックコマンド#
Cythonプログラムはコンパイルが必要で、通常はsetup.py
プログラムを書く必要があります。これについては後で詳しく説明します。まず、Notebookの%%cython
マジックコマンドを使用してCythonプログラムを迅速にコンパイルして実行します。
%%cython
マジックコマンドはセルコマンドで、セル内のプログラム全体がCythonによって拡張モジュールにコンパイルされ、コンパイル後のモジュールからすべてのオブジェクトが自動的にロードされます。Cythonプログラムは独立したモジュールであるため、このモジュール内でnumpy
ライブラリを再ロードする必要があります。以下のプログラムはpairwise_dist_python()
と完全に同じで、その計算速度もpairwise_dist_python()
とほぼ同じで、わずかな改善しか見られません:
%load_ext helper.cython
%%cython
import numpy as np
def pairwise_dist_cython(X):
m, n = X.shape
D = np.empty((m, m), dtype=np.float64)
for i in xrange(m):
for j in xrange(i, m):
d = 0.0
for k in xrange(n):
tmp = X[i, k] - X[j, k]
d += tmp * tmp
D[i, j] = D[j, i] = d ** 0.5
return D
以下でその計算速度をテストします:
%timeit pairwise_dist_cython(X)
np.allclose(pairwise_dist_cython(X), pairwise_dist_python(X))
40.3 ms ± 698 μs per loop (mean ± std. dev. of 7 runs, 10 loops each)
True
以下のプログラムでは、Cythonのすべての最適化手段を使用しています。cdef
は変数の型を宣言するキーワードで、@cython
はCythonのコンパイルを指示するコマンドです。
%%cython
import numpy as np
import cython
from libc.math cimport sqrt
@cython.boundscheck(False)
@cython.wraparound(False)
def pairwise_dist_cython2(double[:, ::1] X):
cdef int m, n, i, j, k
cdef double tmp, d
m, n = X.shape[0], X.shape[1]
cdef double[:, ::1] D = np.empty((m, m), dtype=np.float64)
for i in range(m):
for j in range(i, m):
d = 0.0
for k in range(n):
tmp = X[i, k] - X[j, k]
d += tmp * tmp
D[i, j] = D[j, i] = sqrt(d)
return np.asarray(D)
Cythonでコンパイルされた関数は、SciPyのpdist()
よりも少し速いです:
%timeit pairwise_dist_cython2(X)
np.allclose(pairwise_dist_cython2(X), pairwise_dist_python(X))
108 μs ± 1.61 μs per loop (mean ± std. dev. of 7 runs, 10,000 loops each)
True
Cythonで拡張モジュール作成#
NotebookでCythonプログラムの速度と正確性をテストした後、他のPythonプログラムから呼び出せるように拡張モジュールとしてコンパイルしたいと考えます。実際、%%cython
マジックコマンドは自動的にCythonコードを拡張モジュールにコンパイルし、そのモジュールからすべてのグローバルオブジェクトをロードします。以下では、sys.modules
を使用してpairwise_dist_cython2()
関数が含まれるモジュールのファイルパスを見つけます。
import sys
sys.modules[pairwise_dist_cython2.__module__]
<module '_cython_magic_e995555757c69007b2c0463c54f193e7e1c94190' from 'C:\\Users\\ruoyu\\.ipython\\cython\\_cython_magic_e995555757c69007b2c0463c54f193e7e1c94190.cp313-win_amd64.pyd'>
Cythonプログラムを拡張モジュールとしてコンパイルするには、以下のようなsetup_fast_pdist.py
インストールスクリプトを書く必要があります。その中で、"fast_pdist"
という名前の拡張モジュールを表すExtension
オブジェクトを定義します。これにはCythonソースプログラムfast_pdist.pyx
が含まれており、このファイルの内容は前に定義したpairwise_dist_cython2()
と同じです。NumPyの機能を使用しているため、numpy.get_include()
を使用してNumPyのヘッダーファイルのパスを指定する必要があります。
%%writefile setup_fast_pdist.py
from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
import numpy as np
ext_modules = [
Extension("fast_pdist", ["fast_pdist.pyx"],
include_dirs = [np.get_include()]),
]
setup(
name = 'a faster version of pdist',
cmdclass = {'build_ext': build_ext},
ext_modules = ext_modules
)
Overwriting setup_fast_pdist.py
%%writefile fast_pdist.pyx
# cython: language_level=3
import numpy as np
import cython
from libc.math cimport sqrt
@cython.boundscheck(False)
@cython.wraparound(False)
def pairwise_dist_cython2(double[:, ::1] X):
cdef int m, n, i, j, k
cdef double tmp, d
m, n = X.shape[0], X.shape[1]
cdef double[:, ::1] D = np.empty((m, m), dtype=np.float64)
for i in range(m):
for j in range(i, m):
d = 0.0
for k in range(n):
tmp = X[i, k] - X[j, k]
d += tmp * tmp
D[i, j] = D[j, i] = sqrt(d)
return np.asarray(D)
Overwriting fast_pdist.pyx
コマンドラインで以下のコマンドを実行すると、Cythonファイルがコンパイルされ、拡張モジュールfast_pdist.pyd
が生成されます。コマンドライン引数--inplace
は、生成された拡張モジュールがソースプログラムと同じパスに配置されることを意味します。
!python setup_fast_pdist.py build_ext --inplace
Compiling fast_pdist.pyx because it changed.
[1/1] Cythonizing fast_pdist.pyx
次に、この拡張モジュールをロードし、その中の関数を呼び出すことができます:
import fast_pdist
np.allclose(fast_pdist.pairwise_dist_cython2(X), pairwise_dist_python(X))
True
setup.py
を使用する以外に、cythonize
コマンドを使用することもできます。引数-i
の意味は、前述の--inplace
と同じです。
Warning
Notebook環境で既にロードされている拡張モジュールを再コンパイルすると、ファイルが上書きできないエラーが発生することがあります。この場合、Notebookのカーネルプロセスを再起動して再コンパイルする必要があります。
!cythonize -i fast_pdist.pyx
error: could not delete 'fast_pdist.cp313-win_amd64.pyd': 拒\u7edd\u8bbf\u95ee。
Cythonの基本文法#
このセクションでは、Cythonの出力コードを確認しながら、基本的な文法について説明します。
C言語でのPythonオブジェクト型#
Cython言語で書かれたプログラムは、CythonコンパイラによってC言語プログラムにコンパイルされ、その後C言語コンパイラによってPythonから呼び出せる拡張モジュールにコンパイルされます。Cython言語では、Pythonのオブジェクト型とC言語の型を処理することができます。読者がCython言語をより深く理解するために、このセクションでは、C言語でPythonオブジェクトがどのように表現されているかを紹介します。
PythonはC言語で書かれており、すべてのオブジェクトはC言語の構造体で表現されます。どのようなオブジェクトの構造体でも、最初の2つのフィールドの意味は固定されています:
ob_refcnt
:このオブジェクトの参照カウントで、参照カウントが0になると、この構造体が占めるメモリが解放されます。ob_type
:型オブジェクトへのポインタ。
PythonのC言語プログラムでは、PyObject
構造体が定義されており、これには上記の2つのフィールドしかありません。他のオブジェクト型の構造体は、その後新しいフィールドを追加します。例えば、float
型に対応する構造体は以下のようになります:
typedef struct {
Py_ssize_t ob_refcnt;
struct _typeobject *ob_type;
double ob_fval;
} PyFloatObject;
64ビットシステムでは、Py_ssize_t
とポインタは8バイトで表現され、double
型の長さも8バイトです。したがって、Pythonのfloat
オブジェクトは24バイトを占めます。以下では、sys.getsizeof()
を使用して浮動小数点数オブジェクト1.0のバイト数を取得します:
import sys
sys.getsizeof(1.0)
24
すべてのPythonオブジェクトの最初の2つのフィールドの型は同じであるため、C言語ではPyObject *
型のポインタを使用してPythonオブジェクトを表現します。Cythonでは、型が宣言されていないすべての変数はPythonオブジェクトを表すため、C言語にコンパイルされると、これらの変数はすべてPyObject *
型のポインタになります。Pythonは、C言語で拡張モジュールを書くためのPython/C APIを提供しています。
以下では、%%cython
マジックコマンドに-a
引数を追加して、Cythonプログラムとそのコンパイル後のC言語プログラムの関係を表示します。Notebookで以下のプログラムを実行すると、次のような折りたたみ可能なコードが出力されます。各行のCythonコードがコンパイルされたC言語コードの量を色で表示します。色が濃いほど、その行のコードの実行速度が遅くなる可能性があります。これらのC言語コードはすべてCythonによって自動生成されているため、読むのが難しいかもしれませんが、読者はこれらのコードを完全に理解する必要はなく、大まかな動作原理を理解するだけで十分です。
%%cython -a
a = 1.0
b = 2.0
c = a + b
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: a = 1.0
if (PyDict_SetItem(__pyx_d, __pyx_n_s_a, __pyx_float_1_0) < 0) __PYX_ERR(0, 1, __pyx_L1_error) /* … */ __pyx_t_4 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_4) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
+2: b = 2.0
if (PyDict_SetItem(__pyx_d, __pyx_n_s_b, __pyx_float_2_0) < 0) __PYX_ERR(0, 2, __pyx_L1_error)
+3: c = a + b
__Pyx_GetModuleGlobalName(__pyx_t_2, __pyx_n_s_a); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __Pyx_GetModuleGlobalName(__pyx_t_3, __pyx_n_s_b); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_3); __pyx_t_4 = PyNumber_Add(__pyx_t_2, __pyx_t_3); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; if (PyDict_SetItem(__pyx_d, __pyx_n_s_c, __pyx_t_4) < 0) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
上記のプログラムでは、変数a
、b
、c
はすべてPythonオブジェクトであるため、c = a + b
は次のようなC言語コードにコンパイルされます。__Pyx_GOTREF
と__Pyx_DECREF
は、ガベージコレクションが正常に動作することを保証するために、Pythonオブジェクトの参照カウンタを増減します。
__Pyx_GetModuleGlobalName(__pyx_t_2, __pyx_n_s_a); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 3, __pyx_L1_error)
__Pyx_GOTREF(__pyx_t_2);
__Pyx_GetModuleGlobalName(__pyx_t_3, __pyx_n_s_b); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 3, __pyx_L1_error)
__Pyx_GOTREF(__pyx_t_3);
__pyx_t_4 = PyNumber_Add(__pyx_t_2, __pyx_t_3); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 3, __pyx_L1_error)
__Pyx_GOTREF(__pyx_t_4);
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
if (PyDict_SetItem(__pyx_d, __pyx_n_s_c, __pyx_t_4) < 0) __PYX_ERR(0, 3, __pyx_L1_error)
__Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
__pyx_n_s_a
は文字列オブジェクト"a"
です。__Pyx_GetModuleGlobalName(__pyx_t_1, __pyx_n_s_a)
は、グローバルディクショナリから"a"
で表されるPythonオブジェクトを取得し、変数__pyx_t_1
に代入します。PyNumber_Add()
はPython/C API関数で、2つのPythonオブジェクトを数値として加算し、新しいPythonオブジェクト__pyx_t_3
を返します。最後に、PyDict_SetItem()
を呼び出して、__pyx_t_3
をグローバルディクショナリ__pyx_d
に追加し、対応するキーとして__pyx_n_s_c
で表される文字列オブジェクト"c"
を使用します。
したがって、単純な加算演算でも、グローバルディクショナリを2回検索して変数a
とb
で表されるオブジェクトを取得し、API関数を1回呼び出して加算演算を行い、結果をグローバルディクショナリに書き込む必要があります。
cdef
を用いた変数型の宣言#
プログラムの実行速度を大幅に向上させるためには、Pythonプログラムで頻繁に使用される変数に対してcdef
キーワードを使用して変数の型を宣言する必要があります。cdef
キーワードを使用して変数の型を宣言すると、以下の2つの効果があります:
変数のアクセス速度の向上:Pythonでは、グローバル変数やオブジェクトの属性とそれに対応する値の関係はディクショナリに保存されており、これらの変数や属性を読み書きするたびにディクショナリへのアクセスが必要です。
cdef
キーワードを使用して変数を宣言すると、その変数と値の対応関係はコンパイル時に決定されるため、ディクショナリの検索に必要な時間を節約できます。数値処理速度の向上:PythonではすべてのオブジェクトがC言語の構造体であり、それらを演算するためにはPythonが提供するAPI関数を呼び出す必要があります。値の型がわかっている場合、Cythonは可能な限りC言語が提供する演算機能を使用するため、計算速度が大幅に向上します。
以下の例では、cdef
を使用してa
、b
、c
の3つの変数をdouble
型として定義しています。生成されるC言語コードは以下のようになります:
static double a, b, c;
...
PyMODINIT_FUNC init(){
...
a = 1.0;
b = 2.0;
c = a + b;
...
}
このコードは、3つのdouble
型のグローバル変数を作成し、モジュールの初期化関数でこれらのグローバル変数に値を代入します。加算演算と変数の代入はC言語の操作であり、前節でAPI関数を呼び出すコードに比べて非常に簡潔です。
Warning
ここでcdef
を使用して定義した3つのグローバル変数はC言語のグローバル変数であり、コンパイル後の拡張モジュールからPythonでこれらの値を取得することはできません。
%%cython -a
cdef double a = 1.0
cdef double b = 2.0
cdef double c = a + b
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: cdef double a = 1.0
__pyx_v_54_cython_magic_a9eb7207565ddb950a621eff032ce32dee5975d5_a = 1.0; /* … */ __pyx_t_2 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
+2: cdef double b = 2.0
__pyx_v_54_cython_magic_a9eb7207565ddb950a621eff032ce32dee5975d5_b = 2.0;
+3: cdef double c = a + b
__pyx_v_54_cython_magic_a9eb7207565ddb950a621eff032ce32dee5975d5_c = (__pyx_v_54_cython_magic_a9eb7207565ddb950a621eff032ce32dee5975d5_a + __pyx_v_54_cython_magic_a9eb7207565ddb950a621eff032ce32dee5975d5_b);
C言語の変数とPythonオブジェクトを演算する場合、C言語の変数をPythonオブジェクトに変換してからAPI関数を呼び出して演算を行います。以下の例では、s
はC言語の変数で、a
はPythonのfloat
オブジェクトです。
❶s
とa
を直接加算すると、以下のようなC言語コードが実行されます。_v_s
はdouble
型の変数で、_t_1
から_t_3
はPythonオブジェクトのポインタ、_s_a
は文字列オブジェクト"a"
です:
_t_1 = PyFloat_FromDouble(_v_s); //sをPythonオブジェクトに変換
__Pyx_GetModuleGlobalName(_t_2, _s_a); //グローバルディクショナリから"a"に対応するオブジェクトを取得
_t_3 = PyNumber_Add(_t_1, _t_2); //加算演算
_v_s = __pyx_PyFloat_AsDouble(_t_3); //加算演算の結果をdouble型に変換
❷<double>
を使用してPythonオブジェクトをC言語のdouble
型に強制変換するため、ここでの加算演算はC言語の加算演算子を使用します。以下のようなC言語コードが実行されます。_t_1
はPythonオブジェクトのポインタ、_t_2
はdouble
型です:
__Pyx_GetModuleGlobalName(_t_1, _s_a);
_t_2 = __pyx_PyFloat_AsDouble(_t_1);
_v_s = _v_s + _t_2;
%%cython -a
cdef double s = 0
a = 3.0
s = s + a #❶
s = s + <double>a #❷
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: cdef double s = 0
__pyx_v_54_cython_magic_8c23bbe70ae22a6bfd9ef5ac40a39af6522c6c77_s = 0.0; /* … */ __pyx_t_4 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_4) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
+2: a = 3.0
if (PyDict_SetItem(__pyx_d, __pyx_n_s_a, __pyx_float_3_0) < 0) __PYX_ERR(0, 2, __pyx_L1_error)
+3: s = s + a #❶
__pyx_t_2 = PyFloat_FromDouble(__pyx_v_54_cython_magic_8c23bbe70ae22a6bfd9ef5ac40a39af6522c6c77_s); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __Pyx_GetModuleGlobalName(__pyx_t_3, __pyx_n_s_a); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_3); __pyx_t_4 = PyNumber_Add(__pyx_t_2, __pyx_t_3); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; __pyx_t_5 = __pyx_PyFloat_AsDouble(__pyx_t_4); if (unlikely((__pyx_t_5 == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0; __pyx_v_54_cython_magic_8c23bbe70ae22a6bfd9ef5ac40a39af6522c6c77_s = __pyx_t_5;
+4: s = s + <double>a #❷
__Pyx_GetModuleGlobalName(__pyx_t_4, __pyx_n_s_a); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 4, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); __pyx_t_5 = __pyx_PyFloat_AsDouble(__pyx_t_4); if (unlikely((__pyx_t_5 == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 4, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0; __pyx_v_54_cython_magic_8c23bbe70ae22a6bfd9ef5ac40a39af6522c6c77_s = (__pyx_v_54_cython_magic_8c23bbe70ae22a6bfd9ef5ac40a39af6522c6c77_s + ((double)__pyx_t_5));
C言語のさまざまなデータ型に加えて、CythonはPythonの多くの組み込みデータ型、例えばリスト、辞書、タプルなどを認識できます。以下の例では、list
を使用してclist
変数をリストオブジェクトとして宣言し、int
型のcindex
変数をインデックスとしてリストの要素を取得します。比較のために、2つの動的変数pylist
とpyindex
を使用して同じ操作を行います。
%%cython -a
cdef list clist = [1000, 2, 3]
cdef int cindex = 0
clist[cindex] #❶
pylist = [1000, 2, 3]
pyindex = 0
pylist[pyindex] #❷
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: cdef list clist = [1000, 2, 3]
__pyx_t_2 = PyList_New(3); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __Pyx_INCREF(__pyx_int_1000); __Pyx_GIVEREF(__pyx_int_1000); if (__Pyx_PyList_SET_ITEM(__pyx_t_2, 0, __pyx_int_1000)) __PYX_ERR(0, 1, __pyx_L1_error); __Pyx_INCREF(__pyx_int_2); __Pyx_GIVEREF(__pyx_int_2); if (__Pyx_PyList_SET_ITEM(__pyx_t_2, 1, __pyx_int_2)) __PYX_ERR(0, 1, __pyx_L1_error); __Pyx_INCREF(__pyx_int_3); __Pyx_GIVEREF(__pyx_int_3); if (__Pyx_PyList_SET_ITEM(__pyx_t_2, 2, __pyx_int_3)) __PYX_ERR(0, 1, __pyx_L1_error); __Pyx_XGOTREF(__pyx_v_54_cython_magic_30b42740991fe54adc6243f328268e0f7d234184_clist); __Pyx_DECREF_SET(__pyx_v_54_cython_magic_30b42740991fe54adc6243f328268e0f7d234184_clist, ((PyObject*)__pyx_t_2)); __Pyx_GIVEREF(__pyx_t_2); __pyx_t_2 = 0; /* … */ __pyx_t_4 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_4) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
+2: cdef int cindex = 0
__pyx_v_54_cython_magic_30b42740991fe54adc6243f328268e0f7d234184_cindex = 0;
+3: clist[cindex] #❶
if (unlikely(__pyx_v_54_cython_magic_30b42740991fe54adc6243f328268e0f7d234184_clist == Py_None)) { PyErr_SetString(PyExc_TypeError, "'NoneType' object is not subscriptable"); __PYX_ERR(0, 3, __pyx_L1_error) } __pyx_t_2 = __Pyx_GetItemInt_List(__pyx_v_54_cython_magic_30b42740991fe54adc6243f328268e0f7d234184_clist, __pyx_v_54_cython_magic_30b42740991fe54adc6243f328268e0f7d234184_cindex, int, 1, __Pyx_PyInt_From_int, 1, 1, 1); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 3, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
4:
+5: pylist = [1000, 2, 3]
__pyx_t_2 = PyList_New(3); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 5, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __Pyx_INCREF(__pyx_int_1000); __Pyx_GIVEREF(__pyx_int_1000); if (__Pyx_PyList_SET_ITEM(__pyx_t_2, 0, __pyx_int_1000)) __PYX_ERR(0, 5, __pyx_L1_error); __Pyx_INCREF(__pyx_int_2); __Pyx_GIVEREF(__pyx_int_2); if (__Pyx_PyList_SET_ITEM(__pyx_t_2, 1, __pyx_int_2)) __PYX_ERR(0, 5, __pyx_L1_error); __Pyx_INCREF(__pyx_int_3); __Pyx_GIVEREF(__pyx_int_3); if (__Pyx_PyList_SET_ITEM(__pyx_t_2, 2, __pyx_int_3)) __PYX_ERR(0, 5, __pyx_L1_error); if (PyDict_SetItem(__pyx_d, __pyx_n_s_pylist, __pyx_t_2) < 0) __PYX_ERR(0, 5, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
+6: pyindex = 0
if (PyDict_SetItem(__pyx_d, __pyx_n_s_pyindex, __pyx_int_0) < 0) __PYX_ERR(0, 6, __pyx_L1_error)
+7: pylist[pyindex] #❷
__Pyx_GetModuleGlobalName(__pyx_t_2, __pyx_n_s_pylist); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 7, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __Pyx_GetModuleGlobalName(__pyx_t_3, __pyx_n_s_pyindex); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 7, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_3); __pyx_t_4 = __Pyx_PyObject_GetItem(__pyx_t_2, __pyx_t_3); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 7, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
❶clist
とcindex
の型が宣言されているため、clist[cindex]
は以下のコードにコンパイルされ、リストオブジェクト内の特定のインデックスに対応する要素を取得するためにヘルパー関数__Pyx_GetItemInt_List()
を直接呼び出します:
__Pyx_GetItemInt_List(clist, cindex, int, 1, __Pyx_PyInt_From_int, 1, 1, 1)
❷pylist
とpyindex
変数には型宣言がないため、これらをPythonのオブジェクトとして扱います。pylist[pyindex]
は以下のコードにコンパイルされます:
__Pyx_GetModuleGlobalName(__pyx_t_1, __pyx_n_s_pylist);
__Pyx_GetModuleGlobalName(__pyx_t_2, __pyx_n_s_pyindex);
__pyx_t_3 = __Pyx_PyObject_GetItem(__pyx_t_1, __pyx_t_2);
まず、__Pyx_GetModuleGlobalName()
を呼び出してグローバル変数辞書からpylist
とpyindex
に対応するオブジェクトを取得し、次に__Pyx_PyObject_GetItem()
を呼び出してインデックスに対応する要素を取得します。
[]
を使用したインデックスアクセス操作に加えて、Cythonは多くの内部型のメソッドと属性を認識できます。次のテーブルには、現在Cythonが認識できるPythonデータ型と関連する演算操作がリストされています。表内の既知の属性やメソッドに遭遇した場合、Cythonは効率的なコードを生成します。
データ型 |
認識可能な演算 |
---|---|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
前のプログラムにclist.append(pyindex)
を追加して、append()
がコンパイルされた後のプログラムを確認できます。未知の属性やメソッドに遭遇した場合、Pythonオブジェクトが提供する汎用APIインターフェースが呼び出されます。例えば、clist.count(0)
などです。
def
を用いた関数定義#
通常、プログラムの実行速度を向上させるために、Cythonでdef
キーワードを使用して関数を定義し、Pythonからそれらを呼び出します。Pythonから呼び出す必要があるため、def
で定義された関数の引数と戻り値はPythonオブジェクトです。しかし、Cythonでは関数の引数に型定義を追加することができ、Cythonはこれらの引数を型チェックし、対応するC言語の型に自動的に変換します。
以下のpy_square_add()
関数の2つの引数はdouble
型です。Pythonから呼び出す際には、2つのPythonのfloat
オブジェクトが渡され、py_square_add()
内部ではこれら2つのfloat
オブジェクトがC言語のdouble
型の変数に変換され、double
型の結果を計算した後、それをPythonのfloat
オブジェクトに変換して返します。
対応するC言語のコードは以下のようになります。ここで、x
、y
、_r
はPythonのオブジェクトポインタ型であり、_v_x
と_v_y
はC言語のdouble
型変数です。API関数PyFloat_AsDouble()
とPyFloat_FromDouble()
を使用して、float
オブジェクトとdouble
型変数の間で変換を行います。
double _v_x = PyFloat_AsDouble(x);
double _v_y = PyFloat_AsDouble(y);
PyObject * _r;
_r = PyFloat_FromDouble(((_v_x * _v_x) + (_v_y * _v_y)));
%%cython -a
def py_square_add(double x, double y):
return x*x + y*y
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: def py_square_add(double x, double y):
/* Python wrapper */ static PyObject *__pyx_pw_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_1py_square_add(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ); /*proto*/ static PyMethodDef __pyx_mdef_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_1py_square_add = {"py_square_add", (PyCFunction)(void*)(__Pyx_PyCFunction_FastCallWithKeywords)__pyx_pw_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_1py_square_add, __Pyx_METH_FASTCALL|METH_KEYWORDS, 0}; static PyObject *__pyx_pw_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_1py_square_add(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ) { double __pyx_v_x; double __pyx_v_y; #if !CYTHON_METH_FASTCALL CYTHON_UNUSED Py_ssize_t __pyx_nargs; #endif CYTHON_UNUSED PyObject *const *__pyx_kwvalues; PyObject *__pyx_r = 0; __Pyx_RefNannyDeclarations __Pyx_RefNannySetupContext("py_square_add (wrapper)", 0); #if !CYTHON_METH_FASTCALL #if CYTHON_ASSUME_SAFE_MACROS __pyx_nargs = PyTuple_GET_SIZE(__pyx_args); #else __pyx_nargs = PyTuple_Size(__pyx_args); if (unlikely(__pyx_nargs < 0)) return NULL; #endif #endif __pyx_kwvalues = __Pyx_KwValues_FASTCALL(__pyx_args, __pyx_nargs); { PyObject **__pyx_pyargnames[] = {&__pyx_n_s_x,&__pyx_n_s_y,0}; PyObject* values[2] = {0,0}; if (__pyx_kwds) { Py_ssize_t kw_args; switch (__pyx_nargs) { case 2: values[1] = __Pyx_Arg_FASTCALL(__pyx_args, 1); CYTHON_FALLTHROUGH; case 1: values[0] = __Pyx_Arg_FASTCALL(__pyx_args, 0); CYTHON_FALLTHROUGH; case 0: break; default: goto __pyx_L5_argtuple_error; } kw_args = __Pyx_NumKwargs_FASTCALL(__pyx_kwds); switch (__pyx_nargs) { case 0: if (likely((values[0] = __Pyx_GetKwValue_FASTCALL(__pyx_kwds, __pyx_kwvalues, __pyx_n_s_x)) != 0)) { (void)__Pyx_Arg_NewRef_FASTCALL(values[0]); kw_args--; } else if (unlikely(PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) else goto __pyx_L5_argtuple_error; CYTHON_FALLTHROUGH; case 1: if (likely((values[1] = __Pyx_GetKwValue_FASTCALL(__pyx_kwds, __pyx_kwvalues, __pyx_n_s_y)) != 0)) { (void)__Pyx_Arg_NewRef_FASTCALL(values[1]); kw_args--; } else if (unlikely(PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) else { __Pyx_RaiseArgtupleInvalid("py_square_add", 1, 2, 2, 1); __PYX_ERR(0, 1, __pyx_L3_error) } } if (unlikely(kw_args > 0)) { const Py_ssize_t kwd_pos_args = __pyx_nargs; if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_kwvalues, __pyx_pyargnames, 0, values + 0, kwd_pos_args, "py_square_add") < 0)) __PYX_ERR(0, 1, __pyx_L3_error) } } else if (unlikely(__pyx_nargs != 2)) { goto __pyx_L5_argtuple_error; } else { values[0] = __Pyx_Arg_FASTCALL(__pyx_args, 0); values[1] = __Pyx_Arg_FASTCALL(__pyx_args, 1); } __pyx_v_x = __pyx_PyFloat_AsDouble(values[0]); if (unlikely((__pyx_v_x == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) __pyx_v_y = __pyx_PyFloat_AsDouble(values[1]); if (unlikely((__pyx_v_y == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) } goto __pyx_L6_skip; __pyx_L5_argtuple_error:; __Pyx_RaiseArgtupleInvalid("py_square_add", 1, 2, 2, __pyx_nargs); __PYX_ERR(0, 1, __pyx_L3_error) __pyx_L6_skip:; goto __pyx_L4_argument_unpacking_done; __pyx_L3_error:; { Py_ssize_t __pyx_temp; for (__pyx_temp=0; __pyx_temp < (Py_ssize_t)(sizeof(values)/sizeof(values[0])); ++__pyx_temp) { __Pyx_Arg_XDECREF_FASTCALL(values[__pyx_temp]); } } __Pyx_AddTraceback("_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1.py_square_add", __pyx_clineno, __pyx_lineno, __pyx_filename); __Pyx_RefNannyFinishContext(); return NULL; __pyx_L4_argument_unpacking_done:; __pyx_r = __pyx_pf_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_py_square_add(__pyx_self, __pyx_v_x, __pyx_v_y); int __pyx_lineno = 0; const char *__pyx_filename = NULL; int __pyx_clineno = 0; /* function exit code */ { Py_ssize_t __pyx_temp; for (__pyx_temp=0; __pyx_temp < (Py_ssize_t)(sizeof(values)/sizeof(values[0])); ++__pyx_temp) { __Pyx_Arg_XDECREF_FASTCALL(values[__pyx_temp]); } } __Pyx_RefNannyFinishContext(); return __pyx_r; } static PyObject *__pyx_pf_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_py_square_add(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_x, double __pyx_v_y) { PyObject *__pyx_r = NULL; /* … */ /* function exit code */ __pyx_L1_error:; __Pyx_XDECREF(__pyx_t_1); __Pyx_AddTraceback("_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1.py_square_add", __pyx_clineno, __pyx_lineno, __pyx_filename); __pyx_r = NULL; __pyx_L0:; __Pyx_XGIVEREF(__pyx_r); __Pyx_RefNannyFinishContext(); return __pyx_r; } /* … */ __pyx_tuple_ = PyTuple_Pack(2, __pyx_n_s_x, __pyx_n_s_y); if (unlikely(!__pyx_tuple_)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_tuple_); __Pyx_GIVEREF(__pyx_tuple_); /* … */ __pyx_t_2 = __Pyx_CyFunction_New(&__pyx_mdef_54_cython_magic_6708ce801a2c4515da3cb93bf03e1b748ea510a1_1py_square_add, 0, __pyx_n_s_py_square_add, NULL, __pyx_n_s_cython_magic_6708ce801a2c4515da, __pyx_d, ((PyObject *)__pyx_codeobj__2)); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_py_square_add, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; __pyx_t_2 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
+2: return x*x + y*y
__Pyx_XDECREF(__pyx_r); __pyx_t_1 = PyFloat_FromDouble(((__pyx_v_x * __pyx_v_x) + (__pyx_v_y * __pyx_v_y))); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_1); __pyx_r = __pyx_t_1; __pyx_t_1 = 0; goto __pyx_L0;
Pythonの型(例えばlist
)を使用して引数を宣言する場合、Cythonは型チェックを行い、認識可能な演算をより効率的なコードに変換します。以下の例では、alist
引数をlist
型として宣言しているため、len(alist)
はPyList_GET_SIZE(__pyx_v_alist)
にコンパイルされ、alist[i]
は__Pyx_GetItemInt_List()
関数の呼び出しにコンパイルされます:
%%cython -a
def sum_list(list alist):
cdef double s = 0
cdef int i = 0
for i in range(len(alist)):
s += <double>alist[i]
return s
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: def sum_list(list alist):
/* Python wrapper */ static PyObject *__pyx_pw_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_1sum_list(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ); /*proto*/ static PyMethodDef __pyx_mdef_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_1sum_list = {"sum_list", (PyCFunction)(void*)(__Pyx_PyCFunction_FastCallWithKeywords)__pyx_pw_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_1sum_list, __Pyx_METH_FASTCALL|METH_KEYWORDS, 0}; static PyObject *__pyx_pw_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_1sum_list(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ) { PyObject *__pyx_v_alist = 0; #if !CYTHON_METH_FASTCALL CYTHON_UNUSED Py_ssize_t __pyx_nargs; #endif CYTHON_UNUSED PyObject *const *__pyx_kwvalues; PyObject *__pyx_r = 0; __Pyx_RefNannyDeclarations __Pyx_RefNannySetupContext("sum_list (wrapper)", 0); #if !CYTHON_METH_FASTCALL #if CYTHON_ASSUME_SAFE_MACROS __pyx_nargs = PyTuple_GET_SIZE(__pyx_args); #else __pyx_nargs = PyTuple_Size(__pyx_args); if (unlikely(__pyx_nargs < 0)) return NULL; #endif #endif __pyx_kwvalues = __Pyx_KwValues_FASTCALL(__pyx_args, __pyx_nargs); { PyObject **__pyx_pyargnames[] = {&__pyx_n_s_alist,0}; PyObject* values[1] = {0}; if (__pyx_kwds) { Py_ssize_t kw_args; switch (__pyx_nargs) { case 1: values[0] = __Pyx_Arg_FASTCALL(__pyx_args, 0); CYTHON_FALLTHROUGH; case 0: break; default: goto __pyx_L5_argtuple_error; } kw_args = __Pyx_NumKwargs_FASTCALL(__pyx_kwds); switch (__pyx_nargs) { case 0: if (likely((values[0] = __Pyx_GetKwValue_FASTCALL(__pyx_kwds, __pyx_kwvalues, __pyx_n_s_alist)) != 0)) { (void)__Pyx_Arg_NewRef_FASTCALL(values[0]); kw_args--; } else if (unlikely(PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) else goto __pyx_L5_argtuple_error; } if (unlikely(kw_args > 0)) { const Py_ssize_t kwd_pos_args = __pyx_nargs; if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_kwvalues, __pyx_pyargnames, 0, values + 0, kwd_pos_args, "sum_list") < 0)) __PYX_ERR(0, 1, __pyx_L3_error) } } else if (unlikely(__pyx_nargs != 1)) { goto __pyx_L5_argtuple_error; } else { values[0] = __Pyx_Arg_FASTCALL(__pyx_args, 0); } __pyx_v_alist = ((PyObject*)values[0]); } goto __pyx_L6_skip; __pyx_L5_argtuple_error:; __Pyx_RaiseArgtupleInvalid("sum_list", 1, 1, 1, __pyx_nargs); __PYX_ERR(0, 1, __pyx_L3_error) __pyx_L6_skip:; goto __pyx_L4_argument_unpacking_done; __pyx_L3_error:; { Py_ssize_t __pyx_temp; for (__pyx_temp=0; __pyx_temp < (Py_ssize_t)(sizeof(values)/sizeof(values[0])); ++__pyx_temp) { __Pyx_Arg_XDECREF_FASTCALL(values[__pyx_temp]); } } __Pyx_AddTraceback("_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39.sum_list", __pyx_clineno, __pyx_lineno, __pyx_filename); __Pyx_RefNannyFinishContext(); return NULL; __pyx_L4_argument_unpacking_done:; if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_alist), (&PyList_Type), 1, "alist", 1))) __PYX_ERR(0, 1, __pyx_L1_error) __pyx_r = __pyx_pf_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_sum_list(__pyx_self, __pyx_v_alist); int __pyx_lineno = 0; const char *__pyx_filename = NULL; int __pyx_clineno = 0; /* function exit code */ goto __pyx_L0; __pyx_L1_error:; __pyx_r = NULL; __pyx_L0:; { Py_ssize_t __pyx_temp; for (__pyx_temp=0; __pyx_temp < (Py_ssize_t)(sizeof(values)/sizeof(values[0])); ++__pyx_temp) { __Pyx_Arg_XDECREF_FASTCALL(values[__pyx_temp]); } } __Pyx_RefNannyFinishContext(); return __pyx_r; } static PyObject *__pyx_pf_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_sum_list(CYTHON_UNUSED PyObject *__pyx_self, PyObject *__pyx_v_alist) { double __pyx_v_s; int __pyx_v_i; PyObject *__pyx_r = NULL; /* … */ /* function exit code */ __pyx_L1_error:; __Pyx_XDECREF(__pyx_t_4); __Pyx_AddTraceback("_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39.sum_list", __pyx_clineno, __pyx_lineno, __pyx_filename); __pyx_r = NULL; __pyx_L0:; __Pyx_XGIVEREF(__pyx_r); __Pyx_RefNannyFinishContext(); return __pyx_r; } /* … */ __pyx_tuple_ = PyTuple_Pack(3, __pyx_n_s_alist, __pyx_n_s_s, __pyx_n_s_i); if (unlikely(!__pyx_tuple_)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_tuple_); __Pyx_GIVEREF(__pyx_tuple_); /* … */ __pyx_t_2 = __Pyx_CyFunction_New(&__pyx_mdef_54_cython_magic_d4c804714e9354a89c70268af2fc7cd550bdae39_1sum_list, 0, __pyx_n_s_sum_list, NULL, __pyx_n_s_cython_magic_d4c804714e9354a89c, __pyx_d, ((PyObject *)__pyx_codeobj__2)); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_sum_list, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; __pyx_t_2 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
+2: cdef double s = 0
__pyx_v_s = 0.0;
+3: cdef int i = 0
__pyx_v_i = 0;
+4: for i in range(len(alist)):
if (unlikely(__pyx_v_alist == Py_None)) { PyErr_SetString(PyExc_TypeError, "object of type 'NoneType' has no len()"); __PYX_ERR(0, 4, __pyx_L1_error) } __pyx_t_1 = __Pyx_PyList_GET_SIZE(__pyx_v_alist); if (unlikely(__pyx_t_1 == ((Py_ssize_t)-1))) __PYX_ERR(0, 4, __pyx_L1_error) __pyx_t_2 = __pyx_t_1; for (__pyx_t_3 = 0; __pyx_t_3 < __pyx_t_2; __pyx_t_3+=1) { __pyx_v_i = __pyx_t_3;
+5: s += <double>alist[i]
if (unlikely(__pyx_v_alist == Py_None)) { PyErr_SetString(PyExc_TypeError, "'NoneType' object is not subscriptable"); __PYX_ERR(0, 5, __pyx_L1_error) } __pyx_t_4 = __Pyx_GetItemInt_List(__pyx_v_alist, __pyx_v_i, int, 1, __Pyx_PyInt_From_int, 1, 1, 1); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 5, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); __pyx_t_5 = __pyx_PyFloat_AsDouble(__pyx_t_4); if (unlikely((__pyx_t_5 == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 5, __pyx_L1_error) __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0; __pyx_v_s = (__pyx_v_s + ((double)__pyx_t_5)); }
+6: return s
__Pyx_XDECREF(__pyx_r); __pyx_t_4 = PyFloat_FromDouble(__pyx_v_s); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 6, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_4); __pyx_r = __pyx_t_4; __pyx_t_4 = 0; goto __pyx_L0;
cdef
を用いたC関数定義#
def
関数はPythonの呼び出しインターフェースを使用するため、Cythonプログラム内部でこれらの関数を呼び出す場合でも、引数と戻り値の型変換が必要です。大量のループでこのような関数を呼び出すと、大きなオーバーヘッドが発生します。cdef
キーワードを使用して、Cythonプログラム内部でのみ呼び出せる関数を定義することができ、呼び出しのオーバーヘッドはC言語の関数と同じです。
以下のc_square_add()
はcdef
を使用して定義され、引数と戻り値はdouble
型です。a = c_square_add(1.0, 2.0)
は以下のコードにコンパイルされます:
__pyx_v_a = __pyx_f_c_square_add(1.0, 2.0);
Tip
cdef
関数の戻り値の型を宣言しない場合、その型はPythonオブジェクトになります。
%%cython -a
cdef double c_square_add(double x, double y):
return x*x + y*y
cdef double a = c_square_add(1.0, 2.0)
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: cdef double c_square_add(double x, double y):
static double __pyx_f_54_cython_magic_3719867d67f48edc1093c858a4ef3d73694a22f6_c_square_add(double __pyx_v_x, double __pyx_v_y) { double __pyx_r; /* … */ /* function exit code */ __pyx_L0:; return __pyx_r; } /* … */ __pyx_t_3 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_3); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_3) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
+2: return x*x + y*y
__pyx_r = ((__pyx_v_x * __pyx_v_x) + (__pyx_v_y * __pyx_v_y)); goto __pyx_L0;
3:
+4: cdef double a = c_square_add(1.0, 2.0)
__pyx_t_2 = __pyx_f_54_cython_magic_3719867d67f48edc1093c858a4ef3d73694a22f6_c_square_add(1.0, 2.0); if (unlikely(__pyx_t_2 == ((double)-1) && PyErr_Occurred())) __PYX_ERR(0, 4, __pyx_L1_error) __pyx_v_54_cython_magic_3719867d67f48edc1093c858a4ef3d73694a22f6_a = __pyx_t_2;
同じ関数をCython内で高速に呼び出し、かつPythonからも呼び出せるようにしたい場合、cpdef
キーワードを使用できます。これにより、C言語関数とPythonから呼び出すためのラッパー関数の両方が生成されます。Cython内でcpdef
関数を呼び出す場合、cdef
関数よりも若干時間がかかりますが、def
関数を呼び出すよりもはるかに高速です。
%%cython -a
cpdef double cp_square_add(double x, double y):
return x*x + y*y
cp_square_add(1.0, 2.0)
Generated by Cython 3.0.12
Yellow lines hint at Python interaction.
Click on a line that starts with a "+
" to see the C code that Cython generated for it.
+1: cpdef double cp_square_add(double x, double y):
static PyObject *__pyx_pw_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_1cp_square_add(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ); /*proto*/ static double __pyx_f_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_cp_square_add(double __pyx_v_x, double __pyx_v_y, CYTHON_UNUSED int __pyx_skip_dispatch) { double __pyx_r; /* … */ /* function exit code */ __pyx_L0:; return __pyx_r; } /* Python wrapper */ static PyObject *__pyx_pw_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_1cp_square_add(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ); /*proto*/ static PyMethodDef __pyx_mdef_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_1cp_square_add = {"cp_square_add", (PyCFunction)(void*)(__Pyx_PyCFunction_FastCallWithKeywords)__pyx_pw_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_1cp_square_add, __Pyx_METH_FASTCALL|METH_KEYWORDS, 0}; static PyObject *__pyx_pw_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_1cp_square_add(PyObject *__pyx_self, #if CYTHON_METH_FASTCALL PyObject *const *__pyx_args, Py_ssize_t __pyx_nargs, PyObject *__pyx_kwds #else PyObject *__pyx_args, PyObject *__pyx_kwds #endif ) { double __pyx_v_x; double __pyx_v_y; #if !CYTHON_METH_FASTCALL CYTHON_UNUSED Py_ssize_t __pyx_nargs; #endif CYTHON_UNUSED PyObject *const *__pyx_kwvalues; PyObject *__pyx_r = 0; __Pyx_RefNannyDeclarations __Pyx_RefNannySetupContext("cp_square_add (wrapper)", 0); #if !CYTHON_METH_FASTCALL #if CYTHON_ASSUME_SAFE_MACROS __pyx_nargs = PyTuple_GET_SIZE(__pyx_args); #else __pyx_nargs = PyTuple_Size(__pyx_args); if (unlikely(__pyx_nargs < 0)) return NULL; #endif #endif __pyx_kwvalues = __Pyx_KwValues_FASTCALL(__pyx_args, __pyx_nargs); { PyObject **__pyx_pyargnames[] = {&__pyx_n_s_x,&__pyx_n_s_y,0}; PyObject* values[2] = {0,0}; if (__pyx_kwds) { Py_ssize_t kw_args; switch (__pyx_nargs) { case 2: values[1] = __Pyx_Arg_FASTCALL(__pyx_args, 1); CYTHON_FALLTHROUGH; case 1: values[0] = __Pyx_Arg_FASTCALL(__pyx_args, 0); CYTHON_FALLTHROUGH; case 0: break; default: goto __pyx_L5_argtuple_error; } kw_args = __Pyx_NumKwargs_FASTCALL(__pyx_kwds); switch (__pyx_nargs) { case 0: if (likely((values[0] = __Pyx_GetKwValue_FASTCALL(__pyx_kwds, __pyx_kwvalues, __pyx_n_s_x)) != 0)) { (void)__Pyx_Arg_NewRef_FASTCALL(values[0]); kw_args--; } else if (unlikely(PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) else goto __pyx_L5_argtuple_error; CYTHON_FALLTHROUGH; case 1: if (likely((values[1] = __Pyx_GetKwValue_FASTCALL(__pyx_kwds, __pyx_kwvalues, __pyx_n_s_y)) != 0)) { (void)__Pyx_Arg_NewRef_FASTCALL(values[1]); kw_args--; } else if (unlikely(PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) else { __Pyx_RaiseArgtupleInvalid("cp_square_add", 1, 2, 2, 1); __PYX_ERR(0, 1, __pyx_L3_error) } } if (unlikely(kw_args > 0)) { const Py_ssize_t kwd_pos_args = __pyx_nargs; if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_kwvalues, __pyx_pyargnames, 0, values + 0, kwd_pos_args, "cp_square_add") < 0)) __PYX_ERR(0, 1, __pyx_L3_error) } } else if (unlikely(__pyx_nargs != 2)) { goto __pyx_L5_argtuple_error; } else { values[0] = __Pyx_Arg_FASTCALL(__pyx_args, 0); values[1] = __Pyx_Arg_FASTCALL(__pyx_args, 1); } __pyx_v_x = __pyx_PyFloat_AsDouble(values[0]); if (unlikely((__pyx_v_x == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) __pyx_v_y = __pyx_PyFloat_AsDouble(values[1]); if (unlikely((__pyx_v_y == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L3_error) } goto __pyx_L6_skip; __pyx_L5_argtuple_error:; __Pyx_RaiseArgtupleInvalid("cp_square_add", 1, 2, 2, __pyx_nargs); __PYX_ERR(0, 1, __pyx_L3_error) __pyx_L6_skip:; goto __pyx_L4_argument_unpacking_done; __pyx_L3_error:; { Py_ssize_t __pyx_temp; for (__pyx_temp=0; __pyx_temp < (Py_ssize_t)(sizeof(values)/sizeof(values[0])); ++__pyx_temp) { __Pyx_Arg_XDECREF_FASTCALL(values[__pyx_temp]); } } __Pyx_AddTraceback("_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c.cp_square_add", __pyx_clineno, __pyx_lineno, __pyx_filename); __Pyx_RefNannyFinishContext(); return NULL; __pyx_L4_argument_unpacking_done:; __pyx_r = __pyx_pf_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_cp_square_add(__pyx_self, __pyx_v_x, __pyx_v_y); int __pyx_lineno = 0; const char *__pyx_filename = NULL; int __pyx_clineno = 0; /* function exit code */ { Py_ssize_t __pyx_temp; for (__pyx_temp=0; __pyx_temp < (Py_ssize_t)(sizeof(values)/sizeof(values[0])); ++__pyx_temp) { __Pyx_Arg_XDECREF_FASTCALL(values[__pyx_temp]); } } __Pyx_RefNannyFinishContext(); return __pyx_r; } static PyObject *__pyx_pf_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_cp_square_add(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_x, double __pyx_v_y) { PyObject *__pyx_r = NULL; __Pyx_XDECREF(__pyx_r); __pyx_t_1 = __pyx_f_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_cp_square_add(__pyx_v_x, __pyx_v_y, 0); if (unlikely(__pyx_t_1 == ((double)-1) && PyErr_Occurred())) __PYX_ERR(0, 1, __pyx_L1_error) __pyx_t_2 = PyFloat_FromDouble(__pyx_t_1); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); __pyx_r = __pyx_t_2; __pyx_t_2 = 0; goto __pyx_L0; /* function exit code */ __pyx_L1_error:; __Pyx_XDECREF(__pyx_t_2); __Pyx_AddTraceback("_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c.cp_square_add", __pyx_clineno, __pyx_lineno, __pyx_filename); __pyx_r = NULL; __pyx_L0:; __Pyx_XGIVEREF(__pyx_r); __Pyx_RefNannyFinishContext(); return __pyx_r; } /* … */ __pyx_tuple_ = PyTuple_Pack(2, __pyx_n_s_x, __pyx_n_s_y); if (unlikely(!__pyx_tuple_)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_tuple_); __Pyx_GIVEREF(__pyx_tuple_); /* … */ __pyx_t_2 = __Pyx_CyFunction_New(&__pyx_mdef_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_1cp_square_add, 0, __pyx_n_s_cp_square_add, NULL, __pyx_n_s_cython_magic_e6d158981d6f99fd80, __pyx_d, ((PyObject *)__pyx_codeobj__2)); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_cp_square_add, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; /* … */ __pyx_t_2 = __Pyx_PyDict_NewPresized(0); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_GOTREF(__pyx_t_2); if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_2) < 0) __PYX_ERR(0, 1, __pyx_L1_error) __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
+2: return x*x + y*y
__pyx_r = ((__pyx_v_x * __pyx_v_x) + (__pyx_v_y * __pyx_v_y)); goto __pyx_L0;
3:
+4: cp_square_add(1.0, 2.0)
__pyx_t_3 = __pyx_f_54_cython_magic_e6d158981d6f99fd8029338235c32e9bcbcf298c_cp_square_add(1.0, 2.0, 0); if (unlikely(__pyx_t_3 == ((double)-1) && PyErr_Occurred())) __PYX_ERR(0, 4, __pyx_L1_error)
配列の効率的な処理#
科学計算プログラムでは、多くの配列操作が行われます。前節では、Cythonのメモリビュー(MemoryView)を使用して配列の要素に高速にアクセスする方法を簡単に紹介しました。本節では、その詳細な使い方を説明します。
Cythonのメモリビュー#
Cythonのメモリビューは以下の構文で宣言されます:
cdef 要素型[次元宣言] 変数名
次元宣言の形式は以下の通りで、:
は軸を表し、::1
は対応する軸上の要素が連続して格納されていることを示します:
[:]
: 1次元配列[::1]
: 1次元連続配列[:, :]
:2次元配列[:, ::1]
:2次元配列、第1軸の要素が連続して格納
メモリビューはNumPy配列と同様に、shape
、base
、stride
などの情報を保持していますが、データストレージ自体は持たず、他のPythonオブジェクトやC言語の配列からデータストレージを取得します。
Cythonでは、メモリビューには2つの形式があり、使用法に応じてCythonが自動的に切り替えます。C言語レベルでは構造体ですが、Pythonオブジェクトとして使用する必要がある場合、Cythonは自動的にMemoryView
オブジェクトに変換します。以下に、メモリビューの使用例を示します。
以下のプログラムでは、❶buf
はC言語のグローバル2次元配列で、❷view
はメモリビューであり、そのデータストレージをグローバル配列buf
に初期化しています。
❸numpy.asarray()
を使用して、MemoryView
オブジェクトをNumPy配列に変換し、配列のメソッドや関数を呼び出すことができます。❹また、MemoryView
オブジェクトは直接NumPyの関数に渡すこともできます。これらの関数内部では、asarray()
のような関数が呼び出され、引数が配列に変換されます。
❺メモリビューをPython環境に返す場合、CythonはそれをMemoryView
オブジェクトに変換します。したがって、get_view()
を使用して変換後のMemoryView
オブジェクトを取得し、その内容を確認することができます。
❻Cythonでメモリビューのshape
属性を取得することは、メモリビュー構造体のshape
フィールドのデータを取得することに相当します。また、メモリビューの添字操作は、データストレージの対応するアドレスへの直接アクセスにコンパイルされるため、非常に効率的です。
%%cython
import numpy as np
cdef double buf[3][4] #❶
cdef double[:, ::1] view = buf #❷
def fill_value(double value):
np.asarray(view).fill(value) #❸
def sum_view():
return np.sum(view) #❹
def get_view():
return view #❺
def square_view(): #❻
cdef int i, j
for i in range(view.shape[0]):
for j in range(view.shape[1]):
view[i, j] *= view[i, j]
以下に、get_view()
を使用してMemoryView
オブジェクトを取得し、その各種属性を確認します。これらはNumPy配列と全く同じです:
import helper.magics
view = get_view()
%C view.shape; view.strides; view.itemsize; view.nbytes
view.shape view.strides view.itemsize view.nbytes
---------- ------------ ------------- -----------
(3, 4) (32, 8) 8 96
以下に、fill_view()
、square_view()
、sum_view()
などの関数を使用して、C言語のグローバル配列を操作します。MemoryView
オブジェクトは添字演算もサポートしています:
fill_value(3.0)
print(sum_view())
square_view()
print(sum_view())
view[1, 2] = 10
print(sum_view())
36.0
108.0
109.0
以下に、MemoryView
オブジェクトをNumPy配列に変換します。これらはデータストレージを共有していることがわかります:
arr = np.asarray(view)
arr[1, 0] = 11
print(sum_view())
print(arr)
111.0
[[ 9. 9. 9. 9.]
[11. 9. 10. 9.]
[ 9. 9. 9. 9.]]
メモリビューのbase
属性を使用して、実際のデータを保持するオブジェクトを取得できます。データはC言語のグローバル配列に保存されているため、Cythonはそのグローバル配列を表すarray
オブジェクトを作成し、その唯一のmemview
属性を使用して新しいMemoryView
オブジェクトを取得できます。
print(view.base.__class__)
print(view.base.memview)
<class '_cython_magic_d99ab4ab18f2536ba0f0ee5152cd6ccca941a457.array'>
<MemoryView of 'array' object>
NumPy配列とメモリビュー#
以下に、メモリビューを使用してNumPy配列を操作する方法を示します:
❶double[:]
を使用してx
を1次元メモリビューとして宣言します。要素が連続して格納されていることを指定していないため、x[i]
はx.data + i * x.strides[0]
にコンパイルされます。ここで、x.data
はx
のデータストレージの先頭アドレスで、単一バイトのポインタ型です。x.strides[0]
は第0軸上の要素間のバイト間隔数です。
❷res
の要素が連続して格納されていることを宣言しているため、res[i]
は((double *) res.data) + i
にコンパイルされます。添字変数i
と要素間のバイト間隔数の乗算が省略されるため、演算速度が向上します。
❸res.base
を使用してNumPy配列を返します。base
属性はメモリビューのC言語構造体のフィールドではないため、CythonはまずそれをMemoryView
オブジェクトに変換し、Pythonのgetattr()
関数を使用してbase
属性を取得します。
関数内で配列の要素をループ処理する場合、boundscheck
とwraparound
の2つのコンパイルオプションを無効にすることができます。これにより、生成されるC言語コードは配列の添字の範囲チェックを行わず、負の添字もサポートしないため、配列要素へのアクセス速度が向上します。
%%cython
import numpy as np
import cython
@cython.boundscheck(False)
@cython.wraparound(False)
def square(double[:] x): #❶
cdef int i
cdef double[::1] res = np.empty_like(x)
for i in range(x.shape[0]):
res[i] = x[i] * x[i] #❷
return res.base #❸
さらに、メモリビューはNumPy配列と同じスライス添字アクセス機能をサポートしています。メモリビューをスライスアクセスする場合、元のメモリビューとメモリを共有する新しい構造体が作成されるため、多少の演算オーバーヘッドがありますが、配列のスライス操作よりもはるかに高速です。
以下の例では、cpdef
を使用して、Cython内で高速に呼び出せ、かつPythonからも呼び出せるnorm()
関数を定義しています。これは1次元ベクトルをその場で正規化します。norm_axis()
は2次元配列の指定された軸を正規化します。inplace
パラメータがTrue
の場合、その場で正規化し、そうでない場合は新しい正規化された配列を返します。
❶メモリビューはスライス代入演算をサポートしており、この場合、一時的なメモリビュー構造体が作成され、2つのメモリビュー構造体間でデータがコピーされます。❷メモリビューをスライス添字で読み取る場合、新しいビューが作成され、このビューがnorm()
関数に渡されます。
Tip
この例では、#cython
コメントを使用してboundscheck
とwraparound
コンパイルオプションを設定しています。
%%cython
#cython: boundscheck=False, wraparound=False
import numpy as np
import cython
cpdef norm(double[:] x):
cdef double s
cdef int i
s = 0
for i in range(x.shape[0]):
s += x[i]*x[i]
s = 1 / s**0.5
for i in range(x.shape[0]):
x[i] *= s
def norm_axis(double[:, :] x, int axis=0, bint inplace=True):
cdef int i
cdef double[:, :] data
if not inplace:
data = np.empty_like(x)
data[:] = x #❶
else:
data = x
if axis == 1:
for i in range(data.shape[0]):
norm(data[i, :]) #❷
elif axis == 0:
for i in range(data.shape[1]):
norm(data[:, i]) #❷
return data.base
また、アドレス演算子を使用してメモリビューのデータのアドレスを取得し、それをC言語の関数にポインタとして渡すこともできます。以下の例では、❶cimport
を使用してC言語の標準ヘッダーファイルstring.h
からmemcpy()
関数をインポートしています。この関数のシグネチャは以下の通りです:src
が指すアドレスを開始アドレスとする連続したn
バイトを、dst
が指すアドレスを開始アドレスとする領域にコピーします。
void *memcpy(void *dst, const void *src, size_t n);
❷アドレス演算子&
を使用して、メモリビューの最初の要素を指すポインタを取得します。&dst[0]
はdouble *
型のポインタを返します。memcpy()
はバイト長を受け取るため、ここではsizeof(double)
を使用して倍精度浮動小数点数のバイト数を計算し、メモリビューの第0軸の長さを掛けます。
ここで使用しているメモリビューは連続しているため、memcpy()
を直接使用してメモリをコピーできます。連続していない場合、strides
属性もC言語関数に渡す必要があり、これによりメモリビューの要素に正しくアクセスできます。
%%cython
from libc.string cimport memcpy #❶
def copy_memview(double[::1] src, double[::1] dst):
memcpy(&dst[0], &src[0], sizeof(double)*dst.shape[0]) #❷
a = np.random.rand(10)
b = np.zeros_like(a)
copy_memview(a, b)
assert np.all(a == b)
また、型変換操作を使用してC言語のポインタをメモリビューに変換することもできます。例えば、addr
がポインタの場合、<double[:10]>addr
を使用して、長さ10の倍精度浮動小数点数のメモリビューに変換できます。長さは変数で指定することもできます。ただし、このような変換を行う際は、メモリの割り当てと解放に十分注意する必要があります。そうしないと、ダングリングポインタが発生し、プログラム全体がクラッシュする可能性があります。
PythonのオブジェクトとAPI#
Cythonは、Pythonの多くの組み込み型の一般的な操作を認識し、計算速度を向上させることができます。さらに、CythonでPython/C API関数を呼び出すことで、C言語でしか実現できない操作を実現することもできます。
list
オブジェクトの操作#
Pythonでは、まず空のリストを作成し、append()
メソッドを使用して要素を追加するか、[None]*n
を使用してn個の要素を持つリストを作成し、インデックスを使用してリストの内容を設定することができます。Cythonでは、リストを操作するAPI関数を呼び出すことで、リストの作成速度を向上させることができます。
以下のコードは、API関数とappend()
メソッドを使用してリストを作成する速度を比較しています。my_range()
では、以下の3つのAPI関数を呼び出します:
object PyList_New(Py_ssize_t len)
: 指定されたサイズのリストを作成し、リストの要素をNULLに設定します。ここでのNULLはC言語の空アドレスであり、PythonのNoneオブジェクトではありません。つまり、このリストは作成されていますが、その内容は初期化されておらず、Pythonでは使用できません。void PyList_SET_ITEM(object list, Py_ssize_t i, object o)
: リストlist
の指定されたインデックスi
の内容をo
に設定します。これは実際にはC言語のマクロであり、要素がNULLのリストの要素を迅速に設定するために使用できます。void Py_INCREF(object o)
: オブジェクトo
の参照カウンタを1増やします。PyList_SET_ITEM()
を呼び出してオブジェクトo
をリストに追加する場合、リストはオブジェクトo
の参照カウンタを増やす必要がありますが、PyList_SET_ITEM()
は自動的に追加されたオブジェクトの参照カウンタを増やしません。そのため、Py_INCREF()
を呼び出す必要があります。このような参照カウンタを増やさない関数は、Python/C APIのドキュメントで「参照を盗む」と注記されています。
%%cython --compile-args=-w
#cython: boundscheck=False, wraparound=False
from cpython.list cimport PyList_New, PyList_SET_ITEM #❶
from cpython.ref cimport Py_INCREF
def my_range(int n):
cdef int i
cdef object obj #❷
cdef list result
result = PyList_New(n)
for i in range(n):
obj = i
PyList_SET_ITEM(result, i, obj)
Py_INCREF(obj)
return result
def my_range2(int n):
cdef int i
cdef list result
result = []
for i in range(n):
result.append(i)
return result
❶cimport
を使用して、CythonのヘッダーファイルからPyList_New
とPyList_SET_ITEM
の2つの関数宣言を読み込みます。Cythonのヘッダーファイルは.pxd
拡張子を持ち、C言語のヘッダーファイルと同様に、さまざまな関数と型の定義を含むことができます。
CythonのインストールディレクトリのIncludes
ディレクトリに、すべてのCythonヘッダーファイルがあります。numpy
ディレクトリにはNumPy関連の型とAPI関数の宣言が含まれ、libc
ディレクトリにはC言語標準ライブラリの関数の宣言が含まれ、cpython
ディレクトリにはPython/C APIの関数宣言が含まれています。
❷obj
変数は、C言語の整数変数i
をPythonの整数オブジェクトに変換するために一時的に保存されます。
my_range()
は直接目標サイズのリストを作成するため、段階的な拡張によるオーバーヘッドを省き、append()
を使用する場合の1.5倍以上の速度を実現します。
%timeit list(range(100))
%timeit my_range(100)
%timeit my_range2(100)
857 ns ± 28.8 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)
715 ns ± 4.23 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)
1.46 μs ± 32 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)
tuple
オブジェクトの作成#
Pythonではtuple
オブジェクトは不変ですが、CythonでAPI関数を呼び出してtuple
オブジェクトを作成する際にその内容を設定することで、高速にtuple
オブジェクトを作成することができます。以下のto_tuple_list()
は、API関数を呼び出して2次元配列をタプルのリストに変換します。
リストと同様に、タプルにも対応する初期化関数があります:PyTuple_New
とPyTuple_SET_ITEM
で、その使用方法はリストと同じですので、ここでは繰り返しません。
%%cython --compile-args=-w
#cython: boundscheck=False, wraparound=False
from cpython.list cimport PyList_New, PyList_SET_ITEM
from cpython.tuple cimport PyTuple_New, PyTuple_SET_ITEM
from cpython.ref cimport Py_INCREF
def to_tuple_list(double[:, :] arr):
cdef int m, n
cdef int i, j
cdef list result
cdef tuple t
cdef object obj
m, n = arr.shape[0], arr.shape[1]
result = PyList_New(m)
for i in range(m):
t = PyTuple_New(n)
for j in range(n):
obj = arr[i, j]
PyTuple_SET_ITEM(t, j, obj)
Py_INCREF(obj)
PyList_SET_ITEM(result, i, t)
Py_INCREF(t)
return result
以下は、NumPy配列のtolist()
メソッドとto_tuple_list()
の速度を比較します:
import numpy as np
arr = np.random.randint(0, 10, (5, 2)).astype(np.double)
print(to_tuple_list(arr))
arr = np.random.rand(100, 5)
%timeit to_tuple_list(arr)
%timeit arr.tolist()
[(8.0, 9.0), (8.0, 8.0), (5.0, 7.0), (0.0, 9.0), (3.0, 0.0)]
14.4 μs ± 399 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)
17 μs ± 116 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)
array.array
を動的配列として使用#
NumPyの章で紹介したように、Pythonの標準モジュールarray
のarray
オブジェクトを1次元の動的配列として使用することができます。Cythonコードでも、array.array
オブジェクトを動的配列として使用することができます。ただし、Cythonコードでarray.append()
を呼び出して要素を追加する場合、Pythonの関数オブジェクトを呼び出す必要があるため、速度向上の効果はありません。Cythonのcpython/array.pxd
ヘッダーファイルには、array
オブジェクトを迅速に操作するためのcdef
関数が提供されています。
以下のin_circle()
は、2次元座標配列points
の中で、cx
、cy
、r
で表される円の内部にあるすべての座標を収集します。ここで、(cx, cy)
は円の中心座標、r
は円の半径です。円の内部にある点の数が事前にわからないため、プログラムではarray.array
動的配列を使用して条件を満たす点を逐次追加します。
%%cython -c-Ofast --compile-args=-w
#cython: boundscheck=False, wraparound=False
import numpy as np
from cpython cimport array
def in_circle(double[:, :] points, double cx, double cy, double r):
cdef array.array[double] res = array.array("d") #❶
cdef double r2 = r * r
cdef double p[2] #❷
cdef int i
for i in range(points.shape[0]):
p[0] = points[i, 0]
p[1] = points[i, 1]
if (p[0] - cx)**2 + (p[1] - cy)**2 < r2:
array.extend_buffer(res, <char*>p, 2) #❸
return np.frombuffer(res, np.double).copy().reshape(-1, 2) #❹
❶cimport
キーワードを使用してarray
のヘッダーファイルを読み込んだ後、その中のarray.array
型を使用してCython変数res
を定義し、新しいarray.array
オブジェクトを作成してそれに割り当てます。"d"
は要素の型が倍精度浮動小数点数であることを示します。❷p
は2つの要素を持つC言語の配列で、現在処理中の点の座標を一時的に保存するために使用します。❸array.extend_buffer()
を呼び出してp
をres
に追加します。extend_buffer()
の最初の引数はarray.array
オブジェクトで、2番目の引数は追加するデータの先頭アドレスを指すchar*
ポインタ、3番目の引数は追加する要素の数(バイト数ではない)です。❹最後に、numpy.frombuffer()
を使用してres
とメモリを共有するNumPy配列を作成し、Pythonがres
オブジェクトをガベージコレクトできるように配列をコピーします。
以下は、in_circle()
とNumPy関連メソッドの計算速度を比較します。ほとんどの点が円の外側にある場合、in_circle()
の計算速度はより高速になります。
points = np.random.rand(10000, 2)
cx, cy, r = 0.3, 0.5, 0.05
%timeit points[(points[:, 0] - cx)**2 + (points[:, 1] - cy)**2 < r**2, :]
%timeit in_circle(points, cx, cy, r)
82.9 μs ± 1.4 μs per loop (mean ± std. dev. of 7 runs, 10,000 loops each)
25.5 μs ± 690 ns per loop (mean ± std. dev. of 7 runs, 10,000 loops each)