2

我目前正在解决一个问题,我必须对 scipy/numpy 数组的特定条目进行求和,并且我正在寻找一种方法来完全摆脱所有 Python for 循环。我在 Mac OS X 上使用 Python 3.3。下面是我正在做的总和之一的示例,为了示例的缘故,我对填充有随机整数的数组的条目进行求和。

from scipy import ones, conjugate, sum, random

n = 5
M = random.randint(5,size=(4*n**2,4*n**2))
H = sum((M[i+1,:2*n**2]*M[i,:2*n**2].conjugate()).sum() * (-M[i,:2*n**2]*M[i+1,:2*n**2].conjugate()).sum() for i in range(0,2*n**2,2))

我首先计算两个矩阵条目的乘积,然后对一半的列求和。我这样做了两次,然后以两步对一半以上的行求和。

这可能看起来很奇怪,但我正在为格子上的系统使用哈密顿量,其中每一行对应于某个格点,偶数行和奇数行代表向上或向下旋转。n 最终会很大,我需要加快这个总结速度。

现在,我无法弄清楚如何摆脱 for 循环。我尝试在行索引中使用 range() 参数来做到这一点,但这并没有得到相同的结果。

谢谢!

4

2 回答 2

1

这个怎么样:

# We need a complex test case to make sure the conjugate works properly
M = (random.randint(5,size=(4*n**2,4*n**2)) +
     random.randint(5,size=(4*n**2,4*n**2))*1j)

H_bis = np.sum(M[1:2*n**2:2, :2*n**2] * M[:2*n**2:2, :2*n**2].conjugate(),
               axis=-1)
H_bis *= np.sum(-M[:2*n**2:2, :2*n**2] * M[1:2*n**2:2, :2*n**2].conjugate(),
                axis=-1)
H_bis = np.sum(H_bis)
>>> np.allclose(H, H_bis)
True
于 2013-09-20T20:06:48.167 回答
1

在我看来你可以做到

H_bis = np.sum(M[1:2*n**2:2, :2*n**2] * M[:2*n**2:2, :2*n**2].conjugate(), axis=1)
H_bis = H_bis * H_bis.conjugate()
H_bis = -np.sum(H_bis)
于 2013-09-20T20:28:05.103 回答