N 体積分器の二重ループを最適化しようとしていますが、コードの問題は、格納された変数をメモリ ビューの場所に書き込むときに大量のオーバーヘッドが発生することであることがわかりました。
私は元々、このコードを numpy でベクトル化していましたが、パーティクルの位置を更新するために別の for ループ内で呼び出され、オーバーヘッドが非常に大きくなりました。位置 (X) の np.ndarray Nx2 ベクトルがあり、運動量 (XOut) の Nx2 ベクトルを返したいです。以下にリストされている現在のコードはメモリ ビューを返しますが、最終的にはこれを埋め込みたいので問題ありません。このボトルネックをデバッグしたら、他の Cython 関数に関数を追加します。
cython -a "name.pyx" コマンドを試してみたところ、多かれ少なかれすべてが C タイプであることがわかりました。ただし、ループの下部に向かって、 XOut[ii,0] -= valuex のメモリビューへの書き込みが実行時間のほとんどを占めていることがわかりました。これを定数に変更して XOut[ii,0] -= 5 にすると、コードは最大 40 倍速くなります。これは、その行である種のコピー操作を行っているため、速度が低下していることを意味すると思います。私の Cython/C++ のバックグラウンドは初歩的なものですが、ポインターからメモリビューに書き込むように構文を変更する必要があると思います。どんなアドバイスでも大歓迎です。ありがとう！
import numpy as np
cimport numpy as np
from cython.view cimport array as cvarray
cimport cython
from libc.math cimport sinh, cosh, sin, cos, acos, exp, sqrt, fabs, M_PI
DTYPE = np.float64
ctypedef np.float64_t DTYPE_t
cdef DTYPE_t pi = 3.141592653589793
@cython.cdivision(True)
@cython.boundscheck(False) # turn off bounds-checking for entire function
@cython.wraparound(False) # turn off negative index wrapping for entire function
def intTerms(const DTYPE_t[:,:] X, DTYPE_t epsilon, DTYPE_t[:,:] XOut):
cdef Py_ssize_t ii,jj,N
N = X.shape[0]
cdef DTYPE_t valuex,valuey,r2,xvec,yvec
for ii in range(0,N):
for jj in range(ii+1,N):
xvec = X[ii,0]-X[jj,0]
yvec = X[ii,1]-X[jj,1]
r2 = max(xvec**2+yvec**2,epsilon)
valuex = xvec/r2**2
valuey = yvec/r2**2
XOut[ii,0] -= valuex
XOut[ii,1] -= 5 #valuey
XOut[jj,0] += 5 #valuex
XOut[jj,1] += 5 #valuey
XOut[ii,0] /= 2*pi
XOut[ii,1] /= 2*pi
return XOut