6

我正在尝试优化我在一段代码中的循环。我认为以更 numpy 的方式编写它会使其更快,但现在更慢了!方程将长度为 n 的 numpy.array vec 作为输入:

from numpy import *

def f(vec):
    n=len(vec)
    aux=0
    for i in range(n):
        aux = aux + (1- aux)*vec[i]
    return aux

def f2(vec):
    n=len(vec)
    G=tril(array([-vec]*n),-1)+1              #numpy way!
    aux=dot(G.prod(1),vec)
    return aux


if __name__ == '__main__':
    import timeit
    print(timeit.timeit("f(ones(225)+4)", setup="from __main__ import f\nfrom numpy import ones",number=1000))
    print(timeit.timeit("f2(ones(225)+4)", setup="from __main__ import f2\nfrom numpy import ones,tril,dot",number=1000))

0.429496049881 [s]

5.66514706612 [秒]

最后我决定将整个函数插入到我的循环中,从而获得 3 倍的性能提升。我真的在寻找 100 倍的性能提升,但不知道还能做什么。这是我的最终功能:

def CALC_PROB_LOC2(int nSectors, int nZones,double[:] beta, double[:] thetaLoc,np.ndarray[double, ndim=2] h, np.ndarray[double, ndim=2] p, np.ndarray[np.float64_t, ndim=3] U_nij, np.ndarray[double, ndim=2] A_ni):
    cdef np.ndarray[double, ndim=3] Pr_nij  =np.zeros((nSectors,nZones,nZones),dtype="d")
    cdef np.ndarray[double, ndim=2] U_ni    =np.zeros((nSectors,nZones),dtype="d")
    #cdef np.ndarray[np.float64_t, ndim=1] A_ni_pos
    cdef Py_ssize_t n,i,opt
    cdef int aux_bool,options
    cdef np.ndarray[np.float64_t, ndim=1] prob,attractor,optionCosts
    cdef np.ndarray[np.float64_t, ndim=1] eq23,utilities
    cdef double disu
    cdef double eq22
    cdef double aux17
    for n in range(nSectors):
        aux_bool=1
        if n in [0,2,9,10,11,12,13,14,18,19,20]:
            for i in xrange(nZones):
                U_ni[n,i]=p[n,i]+h[n,i]
                Pr_nij[n,i,i]=1
            aux_bool=0
        if aux_bool==1:
            if beta[n]<=0:
                for i in xrange(nZones):
                    U_ni[n,i]=U_nij[n,i,i]
            else:
                A_ni_pos=A_ni[n,:]>0
                options=len(A_ni[n,:][A_ni_pos])
                attractor=A_ni[n,:][A_ni_pos]
                if options>0:
                    for i in xrange(nZones):
                        optionCosts=U_nij[n,i,A_ni_pos]
                        disu=0
                        eq22=0
                        aux17=0
                        prob=np.ones(options)/options #default value
                        if beta[n]==0:
                            Pr_nij[n,i,A_ni_pos],U_ni[n,i]= prob,0
                        if options==1:
                            Pr_nij[n,i,A_ni_pos],U_ni[n,i]= prob,optionCosts
                        else:
                            if thetaLoc[n]<=0:
                                cmin=1
                            else:
                                cmin=(optionCosts**thetaLoc[n]).min()
                                if cmin==0:
                                    cmin=100
                            utilities=optionCosts/cmin
                            eq23=-beta[n]*utilities
                            eq23=np.exp(eq23)
                            aux17=np.dot(attractor,eq23)
                            if aux17==0:
                                Pr_nij[n,i,A_ni_pos],U_ni[n,i]= 0*prob,0
                            else:
                                for opt in range(options):
                                    eq22=eq22+(1-eq22)*eq23[opt]
                                prob=attractor*eq23/aux17
                                disu=cmin*(-np.log(eq22)/beta[n])
                                Pr_nij[n,i,A_ni_pos],U_ni[n,i]= prob,disu


    return Pr_nij,U_ni
4

1 回答 1

9

这就是当线性算法被二次算法取代时会发生的情况:无论执行速度有多快,更好的算法总是会获胜(对于足够大的问题)。

很明显,f它以线性时间f2运行,并且以二次时间运行,因为这是矩阵向量点积的时间复杂度。

对数图清楚地显示了运行时间的差异(线性表示f,二次表示f2):

两种算法比较

可以解释绿线的最右边部分(即,当它不表现为直线时),因为 numpy 函数通常具有很高的开销,这对于不是很小但在运行时间占主导地位的数组来说可以忽略不计。小的。


在已经使用快速算法的 Python 中加速代码的“标准”方法是获取已编译的代码并编写扩展。Cython允许您通过使用一些类型注释来注释 Python 源代码来做到这一点,并且它可以理解 numpy 数组。

通过告诉 Cython 这vec是一个双精度数组、aux一个双精度和i一个整数,它能够生成一个 C 扩展,这对我来说快 400 倍。

def f(double[:] vec):
    n = len(vec)
    cdef double aux = 0
    cdef int i
    for i in range(n):
        aux = aux + (1- aux)*vec[i]
    return aux

如果您碰巧使用的是IPython,您可以运行%load_ext cythonmagic然后将该函数复制到以该行为前缀的单元格中%%cython进行尝试。Cython文档中解释了其他构建和编译它的方法。顺便说一句,IPython 还允许您通过%timeit在语句之前编写 timeit 代码来计时,这真的很方便。

一个完全不同的选择是使用PyPy,这是一个 Python 2.7 实现,带有 JIT 并具有一些基本的 numpy 支持。它可以通过替换import numpypyfor来运行这个小片段import numpy,但它可能无法运行你的整个程序。它比 Cython 慢一点,但它不需要编译器,也不需要注释代码。

于 2013-06-14T13:06:16.580 回答