5

我想加速用 Python 和 NumPy 编写的代码。我使用 Gray-Skott 算法(http://mrob.com/pub/comp/xmorphia/index.html)作为反应扩散模型,但使用 Numba 和 Cython 会更慢!有没有可能加快速度?提前致谢!

Python+NumPy

def GrayScott(counts, Du, Dv, F, k):
    n = 300
    U = np.zeros((n+2,n+2), dtype=np.float_)
    V = np.zeros((n+2,n+2), dtype=np.float_)
    u, v = U[1:-1,1:-1], V[1:-1,1:-1]

    r = 20
    u[:] = 1.0
    U[n/2-r:n/2+r,n/2-r:n/2+r] = 0.50
    V[n/2-r:n/2+r,n/2-r:n/2+r] = 0.25
    u += 0.15*np.random.random((n,n))
    v += 0.15*np.random.random((n,n))

    for i in range(counts):
        Lu = (                 U[0:-2,1:-1] +
              U[1:-1,0:-2] - 4*U[1:-1,1:-1] + U[1:-1,2:] +
                               U[2:  ,1:-1] )
        Lv = (                 V[0:-2,1:-1] +
              V[1:-1,0:-2] - 4*V[1:-1,1:-1] + V[1:-1,2:] +
                               V[2:  ,1:-1] )
        uvv = u*v*v
        u += Du*Lu - uvv + F*(1 - u)
        v += Dv*Lv + uvv - (F + k)*v

    return V

努巴

from numba import jit, autojit

@autojit
def numbaGrayScott(counts, Du, Dv, F, k):
    n = 300
    U = np.zeros((n+2,n+2), dtype=np.float_)
    V = np.zeros((n+2,n+2), dtype=np.float_)
    u, v = U[1:-1,1:-1], V[1:-1,1:-1]

    r = 20
    u[:] = 1.0
    U[n/2-r:n/2+r,n/2-r:n/2+r] = 0.50
    V[n/2-r:n/2+r,n/2-r:n/2+r] = 0.25
    u += 0.15*np.random.random((n,n))
    v += 0.15*np.random.random((n,n))

    Lu = np.zeros_like(u)
    Lv = np.zeros_like(v)

    for i in range(counts):
        for row in range(n):
            for col in range(n):
                Lu[row,col] = U[row+1,col+2] + U[row+1,col] + U[row+2,col+1] + U[row,col+1] - 4*U[row+1,col+1]
                Lv[row,col] = V[row+1,col+2] + V[row+1,col] + V[row+2,col+1] + V[row,col+1] - 4*V[row+1,col+1]

        uvv = u*v*v
        u += Du*Lu - uvv + F*(1 - u)
        v += Dv*Lv + uvv - (F + k)*v

    return V

赛通

%%cython
cimport cython
import numpy as np
cimport numpy as np

cpdef cythonGrayScott(int counts, double Du, double Dv, double F, double k):
    cdef int n = 300
    cdef np.ndarray U = np.zeros((n+2,n+2), dtype=np.float_)
    cdef np.ndarray V = np.zeros((n+2,n+2), dtype=np.float_)
    cdef np.ndarray u = U[1:-1,1:-1]
    cdef np.ndarray v = V[1:-1,1:-1]

    cdef int r = 20
    u[:] = 1.0
    U[n/2-r:n/2+r,n/2-r:n/2+r] = 0.50
    V[n/2-r:n/2+r,n/2-r:n/2+r] = 0.25
    u += 0.15*np.random.random((n,n))
    v += 0.15*np.random.random((n,n))

    cdef np.ndarray Lu = np.zeros_like(u)
    cdef np.ndarray Lv = np.zeros_like(v)
    cdef int i, row, col
    cdef np.ndarray uvv

    for i in range(counts):
        for row in range(n):
            for col in range(n):
                Lu[row,col] = U[row+1,col+2] + U[row+1,col] + U[row+2,col+1] + U[row,col+1] - 4*U[row+1,col+1]
                Lv[row,col] = V[row+1,col+2] + V[row+1,col] + V[row+2,col+1] + V[row,col+1] - 4*V[row+1,col+1]

        uvv = u*v*v
        u += Du*Lu - uvv + F*(1 - u)
        v += Dv*Lv + uvv - (F + k)*v

    return V

使用示例:

GrayScott(4000, 0.16, 0.08, 0.04, 0.06)
4

2 回答 2

7

以下是加速 cython 版本的步骤:

  • cdef np.ndarray不要让元素访问更快,你需要在 cython: 中使用 memoryview: cdef double[:, ::1] bU = U

  • 关闭boundscheckwraparound

  • 在 for 循环中进行所有计算。

这是修改后的cython代码:

%%cython
#cython: boundscheck=False
#cython: wraparound=False
cimport cython
import numpy as np
cimport numpy as np

cpdef cythonGrayScott(int counts, double Du, double Dv, double F, double k):
    cdef int n = 300
    cdef np.ndarray U = np.zeros((n+2,n+2), dtype=np.float_)
    cdef np.ndarray V = np.zeros((n+2,n+2), dtype=np.float_)
    cdef np.ndarray u = U[1:-1,1:-1]
    cdef np.ndarray v = V[1:-1,1:-1]

    cdef int r = 20
    u[:] = 1.0
    U[n/2-r:n/2+r,n/2-r:n/2+r] = 0.50
    V[n/2-r:n/2+r,n/2-r:n/2+r] = 0.25
    u += 0.15*np.random.random((n,n))
    v += 0.15*np.random.random((n,n))

    cdef np.ndarray Lu = np.zeros_like(u)
    cdef np.ndarray Lv = np.zeros_like(v)
    cdef int i, c, r1, c1, r2, c2
    cdef double uvv

    cdef double[:, ::1] bU = U
    cdef double[:, ::1] bV = V
    cdef double[:, ::1] bLu = Lu
    cdef double[:, ::1] bLv = Lv

    for i in range(counts):
        for r in range(n):
            r1 = r + 1
            r2 = r + 2
            for c in range(n):
                c1 = c + 1
                c2 = c + 2
                bLu[r,c] = bU[r1,c2] + bU[r1,c] + bU[r2,c1] + bU[r,c1] - 4*bU[r1,c1]
                bLv[r,c] = bV[r1,c2] + bV[r1,c] + bV[r2,c1] + bV[r,c1] - 4*bV[r1,c1]

        for r in range(n):
            r1 = r + 1
            for c in range(n):
                c1 = c + 1
                uvv = bU[r1,c1]*bV[r1,c1]*bV[r1,c1]
                bU[r1,c1] += Du*bLu[r,c] - uvv + F*(1 - bU[r1,c1])
                bV[r1,c1] += Dv*bLv[r,c] + uvv - (F + k)*bV[r1,c1]

    return V

它比 numpy 版本快大约 11 倍。

于 2014-11-09T01:41:26.463 回答
5

除了循环和所涉及的大量操作之外,在您的情况下最有可能扼杀性能的是数组分配。我不知道为什么您的 Numba 和 Cython 版本没有达到您的期望,但是您可以通过就地执行所有操作来使您的 numpy 代码快 2 倍(以牺牲一些可读性为代价),即替换您当前的循环和:

Lu, Lv, uvv = np.empty_like(u), np.empty_like(v), np.empty_like(u)

for i in range(counts):
    Lu[:] = u
    Lu *= -4
    Lu += U[:-2,1:-1]
    Lu += U[1:-1,:-2]
    Lu += U[1:-1,2:]
    Lu += U[2:,1:-1]
    Lu *= Du

    Lv[:] = v
    Lv *= -4
    Lv += V[:-2,1:-1]
    Lv += V[1:-1,:-2]
    Lv += V[1:-1,2:]
    Lv += V[2:,1:-1]
    Lv *= Dv

    uvv[:] = u
    uvv *= v
    uvv *= v
    Lu -= uvv
    Lv += uvv

    u *= 1 - F
    u += F
    u += Lu

    v *= 1 - F - k
    v += Lv
于 2014-11-09T01:08:21.967 回答