<p>我正在使用Python库MDAnalysis编写代码。
我有一个原子位置的数组(502,3)
我想得到一系列键(原子(I+1)-原子(I)的位置向量)
然后我想得到一个平均张量qab=它本质上是anp.外部(ua、ub)
所有原子的平均值。在</p>
<p>我可以用fortran子例程重写这段代码,但我认为看到一个美味的numpy切片解决方案会更令人愉快:)这是我写的代码,我想让它更快更漂亮。提前谢谢你的帮助。在</p>
<pre><code>bb = u.selectAtoms("all")
coor = bb.positions
print coor.shape # return is (502,3)
# the coordinates of the atoms that we split on 3 dimensions
bbx = coor[:,0]
bby = coor[:,1]
bbz = coor[:,2]
#the bond vectors we obtain by subtructing atoms positions
bbx_ave = bbx[:-1] - bbx[1:]
bby_ave = bby[:-1] - bby[1:]
bbz_ave = bbz[:-1] - bbz[1:]
#now we concentrate everything back so we have one big array
bb_res = np.vstack((bbx_ave, bby_ave, bbz_ave))
# print bb_res.shape # the return is (3,501)
# now we have an array of bonds
# qab - the result tensor which we want to get
nq = len(bb_res)
qab = np.zeros((3,3))
count = 0.
for i in xrange(0, nq):
for j in xrange(0, i+1):
qab = qab + np.outer(bb_res[:,i], bb_res[:,j])
count += 1.
qab = qab/count
print qab
[[ 0.21333394 0.5333341 0. ]
[ 0.5333341 4. 0. ]
[ 0. 0. 0. ]]
</code></pre>