【发布时间】:2017-03-26 23:07:03
【问题描述】:
在下面的代码中,我们计算所有给定点对之间的向量大小。为了在 NumPy 中加速这个操作,我们可以使用广播
import numpy as np
points = np.random.rand(10,3)
pair_vectors = points[:,np.newaxis,:] - points[np.newaxis,:,:]
pair_dists = np.linalg.norm(pair_vectors,axis=2).shape
或外部产品迭代
it = np.nditer([points,points,None], flags=['external_loop'], op_axes=[[0,-1,1],[-1,0,1],None])
for a,b,c in it:
c[...] = b - a
pair_vectors = it.operands[2]
pair_dists = np.linalg.norm(pair_vectors,axis=2)
我的问题是如何使用广播或外部产品迭代来创建一个格式为 10x10x6 的数组,其中最后一个轴包含一对中两个点的坐标(扩展)。并且以一种相关的方式,是否可以直接使用广播或外积迭代来计算对距离,即生成一个 10x10 形式的矩阵,而无需首先计算差异向量(约简)。
为了澄清,以下代码使用慢循环创建所需的矩阵。
pair_coords = np.zeros(10,10,6)
pair_dists = np.zeros(10,10)
for i in range(10):
for j in range(10):
pair_coords[i,j,0:3] = points[i,:]
pair_coords[i,j,3:6] = points[j,:]
pair_dists[i,j] = np.linalg.norm(points[i,:]-points[j,:])
这是使用外积迭代计算距离(或应用任何其他函数,该函数采用一对中两个点的 6 个坐标并产生一个标量)的失败尝试。
res = np.zeros((10,10))
it = np.nditer([points,points,res], flags=['reduce_ok','external_loop'], op_axes=[[0,-1,1],[-1,0,1],None])
for a,b,c in it: c[...] = np.linalg.norm(b-a)
pair_dists = it.operands[2]
【问题讨论】:
-
pair_vectors是 10x10x3。我不明白 10x10x6 会包含什么。 -
这个
norm是:np.sqrt(np.sum(pair_vectors**2,axis=2))。它消除了所有差异;对大小 3 维度求和并取 sqrt。 -
nditer有点像列表zip。它协调多个数组的迭代。它在处理broadcasting方面也超越了zip。但是您已经在c代码 (cython) 中使用它来获得真正的速度提升。 -
我最终在 fortran 和 f2py 中使用了循环,它真的很快。我认为 NumPy 是一个适合后处理的环境,但结果证明对简单的循环操作进行矢量化非常复杂。我想 NumPy 在有一种直接的方法来向量化你的函数时很有用。
标签: python numpy array-broadcasting