<p>你可以做的一件很明显的事情就是更换线路</p>
<pre><code>r_test_fast = reshape_vector(r_test)
</code></pre>
<p>与</p>
^{pr2}$
<p>也许在性能上不会有什么大的不同,但无论如何,使用纽普的内置而不是重新设计轮子是有意义的。在</p>
<p>一般来说,正如您现在可能已经注意到的,优化numpy的诀窍是借助numpy整个数组操作来表达算法,或者至少用切片来表示,而不是迭代python代码中的每个元素。阻止这种“矢量化”的是所谓的循环携带依赖关系,即每次迭代都依赖于前一次迭代的结果的循环。简单地看一下你的代码,你没有这样的东西,应该可以很好地将你的代码矢量化。在</p>
<p><strong>编辑:一个解决方案</strong></p>
<p>我还没有证实这是正确的,但应该告诉你如何处理它。在</p>
<p>首先,取<a href="https://stackoverflow.com/questions/1208118/using-numpy-to-build-an-array-of-all-combinations-of-two-arrays">cartesian() function, which we'll use</a>。那么</p>
<pre>
<code>
def calculate_dipole_vect(mus, r_i, mom_i):
# Treat each mu sequentially
Bs = []
omega = []
for mu in mus:
rel = mu - r_i
r_norm = np.sqrt((rel * rel).sum(1))
r_unit = rel / r_norm[:, np.newaxis]
A = 1e-7
num = A*(3*np.sum(mom_i * r_unit, 0)*r_unit - mom_i)
den = r_norm ** 3
B = np.sum(num / den[:, np.newaxis], 0)
Bs.append(B)
omega.append(gamma_mu * np.sqrt(np.dot(B, B)))
return Bs, omega
# Transpose to get more "natural" ordering with row-major numpy
r_i = r_i.T
mom_i = mom_i.T
t_start = time.clock()
r_frac = cartesian((np.arange(n[0]) / float(n[0]),
np.arange(n[1]) / float(n[1]),
np.arange(n[2]) / float(n[2])))
r_test = np.dot(r_frac, a)
B, omega = calculate_dipole_vect(r_test, r_i, mom_i)
print 'Total time for vectorized: %f s' % (time.clock() - t_start)
</code>
</pre>
<p>好吧,在我的测试中,这实际上比我开始的基于循环的方法稍微慢一些。问题是,在这个问题的原始版本中,它已经通过对shape(20000,3)数组的整个数组操作进行矢量化,所以任何进一步的矢量化都不会带来更多的好处。事实上,如上所述,它可能会使性能恶化,可能是因为大型临时阵列。在</p>