数组与数组列之间的最小欧几里德距离及相应的指数

Question

更新

以防万一这对任何人都有用。如下文所述，欧几里德距离计算和np.argmin几乎占了整个运行时。通过用numba重写距离计算，与已经很快的np.einsum相比，在大多数情况下，我至少可以减少20%。在

@jit(nopython=True)
def calculateDistances_numba(currentLocation, traces)

    deltaX = traces[:, 0, :] - currentLocation[0]
    deltaY = traces[:, 1, :] - currentLocation[1]
    deltaZ = traces[:, 2, :] - currentLocation[2]

    distances = (deltaX**2 + deltaY**2 + deltaZ**2)*0.5

    return distances

~~~~

问题

我有一个大数组，vertices.shape = (N, 3); N ~ 5e6，用来描述非结构化网格的三维顶点。我有n更小的坐标和数据数组，我想将它们线性插值到vertices上。这些存储在另一个数组traces.shape = (L, 3, n); L ~ 2e4; n ~ 2e3的第三个轴上。对于每个顶点（在vertices中的行），我希望快速地找到来自不同小数组的两个最近的点（traces中的页，即它们沿着axis=2的索引不同）。我指的是欧几里德距离d = (deltaX**2 + deltaY**2 + deltaZ**2)。此函数的目的是在两个已知值之间线性插值到顶点中的点。在

我当前的函数运行得相当好，但是对于上面给出的预期数组大小（8小时以上）来说，速度变得非常慢。我已经检查了我的整个代码，可以肯定地说这是一个昂贵的计算。在

当前功能

^{pr2}$

%timeit输出

%timeit interpolater(currentLocation, streamlineBlock, nStreamlines)
10 loops, best of 3: 42.8 ms per loop

由于数据的结构，我只能选择traces的一块来搜索（L~2e3），这明显减少了运行时间。要搜索的括号是当前位置的函数。在

%timeit interpolaterNew(...)
100 loops, best of 3: 6.27 ms per loop

cProfile输出

cProfile告诉我的np.einsum公司以及np.argmin公司是最慢的-事实上他们是计算的绝大多数。请注意，这是针对我的一小部分数据的代码，因此可能不能准确地反映上述函数。在

4251460 function calls (4151427 primitive calls) in 17.907 seconds

Ordered by: internal time

ncalls  tottime  percall  cumtime  percall filename:lineno(function)
50012    7.236    0.000   17.063    0.000 vc3d_traces.py:554(interpolateToVertices)
50012    3.940    0.000    3.940    0.000 {built-in method numpy.core.multiarray.c_einsum}
50012    3.505    0.000    3.505    0.000 {method 'argmin' of 'numpy.ndarray' objects}
100025    0.291    0.000    0.291    0.000 {method 'reduce' of 'numpy.ufunc' objects}
50012    0.289    0.000   17.352    0.000 frame.py:4238(f)
50012    0.223    0.000    0.223    0.000 {method 'searchsorted' of 'numpy.ndarray' objects}
100024    0.191    0.000    0.346    0.000 indexing.py:1815(_convert_key)
1    0.190    0.190   17.905   17.905 {pandas._libs.lib.reduce}
100024    0.159    0.000    0.504    0.000 fromnumeric.py:1730(sum)
100024    0.155    0.000    0.155    0.000 {method 'get_value' of 'pandas._libs.index.IndexEngine' objects}

问题

我现在对如何提高业绩有点不知所措。考虑到距离计算和argmin排序是最昂贵的，是否可以将这些步骤“矢量化”，将计算应用到整个“顶点”数组中？我确实尝试了这个方法，在axis=4中广播，但没有成功——计算机死机了。cProfile报告是否指向其他方面，或者在我的代码中有什么明显不好的地方？有谁能给我一个更好的方法？最后，使用TQM，每秒的迭代次数会有一个大的、快速的减少（在最初的几分钟内减少到175次），这是意料之中的吗？在