是否有函数调用可以替换此代码中的for循环？

这段代码片段是我的一个项目的瓶颈。有函数调用可以替代for循环并加速吗？是否有函数调用可以替换此代码中的for循环？

D = np.zeros((nOcc,nOcc,nVir,nVir)) 
for i in range(nOcc): 
    for j in range(i+1): 
     tmp = Ew[i] + Ew[j] 
     for a in range(nVir): 
     tmp2 = tmp - Ew[a+nOcc] 
     for b in range(a+1): 
      tmp3 = 1.0/(tmp2 - Ew[b+nOcc]) 
      D[i,j,a,b] = Iiajb[i,a,j,b]*tmp3 
      D[i,j,b,a] = Iiajb[i,b,j,a]*tmp3 
      D[j,i,a,b] = D[i,j,b,a] 
      D[j,i,b,a] = D[i,j,a,b]

来源

2013-10-30 John Montgomery

应该做什么功能？ “nOcc”和“nVir”有多大？ – GWW

什么是nOcc，nVir，Ew，Iiajb ... – voithos

可能更多的codereview问题 – samrap

要开始让产生一些任意的数据，这就是服从少数需要原则：

nOcc = 30 
nVir = 120 
Ew = np.random.rand(nOcc+nVir) 
Ew[:nOcc]*=-1 
Ia = np.random.rand(nOcc) 
Ib = np.random.rand(nVir) 
I = np.einsum('a,b,c,d->abcd',Ia,Ib,Ia,Ib)

让包装你的基本代码为例：

def oldcalc_D(Iiajb,nOcc,nVir,Ew): 
    D = np.zeros((nOcc,nOcc,nVir,nVir)) 
    for i in range(nOcc): 
     for j in range(i+1): 
      tmp = Ew[i] + Ew[j] 
      for a in range(nVir): 
      tmp2 = tmp - Ew[a+nOcc] 
      for b in range(a+1): 
       tmp3 = 1.0/(tmp2 - Ew[b+nOcc]) 
       D[i,j,a,b] = Iiajb[i,a,j,b]*tmp3 
       D[i,j,b,a] = Iiajb[i,b,j,a]*tmp3 
       D[j,i,a,b] = D[i,j,b,a] 
       D[j,i,b,a] = D[i,j,a,b] 
    return D

以整体对称的优点是一般一个好策略;然而，仅在numpy的是不值得的成本，所以我们会忽略这一点，简单地向量化代码：

def newcalc_D(I,nOcc,nVir,Ew): 
    O = Ew[:nOcc] 
    V = Ew[nOcc:] 
    D = O[:,None,None,None] - V[:,None,None] + O[:,None] - V 
    return (I/D).swapaxes(1,2)

一些计时：

np.allclose(oldcalc_D(I,nOcc,nVir,Ew),newcalc_D(I,nOcc,nVir,Ew)) 
True 

%timeit newcalc_D(I,nOcc,nVir,Ew) 
1 loops, best of 3: 142 ms per loop 

%timeit oldcalc_D(I,nOcc,nVir,Ew) 
1 loops, best of 3: 15 s per loop

所以只有大约〜100倍的速度更快，因为我说这是一个相当简单的过关，让你知道该怎么做。这可以做得好得多，但应该是计算的一个微不足道的部分，因为整数变换是（O）N^5而在（O）N^4处是这个。对于这些操作，我使用numba的autojit功能：

from numba import autojit 

numba_D = autojit(oldcalc_D) 

%timeit numba_D(I,nOcc,nVir,Ew) 
10 loops, best of 3: 55.1 ms per loop

来源

2013-10-31 04:29:13 Daniel

正是我想要的，非常感谢Ophion！ –

，除非这是Python3，您可能希望通过与xrange更换range启动：前创建整个列表，而后者仅仅是一个迭代器，这是你在这种情况下需要。
对于大型N s，速度差异应该是显而易见的。

此外，看作为你的使用numpy，可能有一个向量化的方法来实现算法。如果是这样的话，矢量化的实现应该快几个数量级。但除非你解释变量和算法，否则我们无法帮助你朝那个方向发展。

来源

2013-10-30 21:40:17 LeartS

是否有函数调用可以替换此代码中的for循环？

回答

相关问题