python - 如何用矢量化计算距离矩阵的种类

Question

我有一个形状为 4 X 3 X 2 的 numpy 数组 A。下面的每一行都是一个节点的二维坐标。 (在我的有限元分析中，每三个节点组成一个三角形。)

array([[[0., 2.],  #node00 
        [2., 2.],  #node01
        [1., 1.]], #node02

       [[0., 2.],  #node10
        [1., 1.],  #node11
        [0., 0.]], #node12

       [[2., 2.],  #node20
        [1., 1.],  #node21
        [2., 0.]], #node22

       [[0., 0.], #node30
        [1., 1.], #node31
        [2., 0.]]]) #node32

我有另一个 numpy 数组 B 预先计算的“中心”的坐标:

array([[1.        , 1.66666667], # center0
       [0.33333333, 1.        ], # center1
       [1.66666667, 1.        ], # center2
       [1.        , 0.33333333]])# center3

我怎样才能像这样有效地计算欧几里德距离的矩阵 C

dist(center0, node00) dist(center0,node01) dist(center0, node02)
dist(center1, node10) dist(center1,node11) dist(center1, node12)
dist(center2, node20) dist(center2,node21) dist(center2, node22)
dist(center3, node30) dist(center3,node31) dist(center3, node32)

其中 dist 表示欧几里得距离公式，如 math.dist 或 numpy.linalg.norm？即结果矩阵的i,j元素为中心i到节点ij的距离。

需要矢量化代码而不是循环，因为我的实际数据来自非常大的医学成像。使用嵌套循环，可以获得如下预期输出:

In [63]: for i in range(4):
    ...:     for j in range(3):
    ...:         C[i,j]=math.dist(A[i,j], B[i]) 

In [67]: C
Out[67]:
array([[1.05409255, 1.05409255, 0.66666667],
       [1.05409255, 0.66666667, 1.05409255],
       [1.05409255, 0.66666667, 1.05409255],
       [1.05409255, 0.66666667, 1.05409255]])

[编辑] 这是与 Pairwise operations (distance) on two lists in numpy 不同的问题，因为需要在此处妥善解决诸如索引之类的问题。

Answer 1

a = np.reshape(A, [12, 2])
b = B[np.repeat(np.arange(4), 3)]
c = np.reshape(np.linalg.norm(a - b, axis=-1), (4, 3))
c
# array([[1.05409255, 1.05409255, 0.66666667],
#        [1.05409255, 0.66666667, 1.05409255],
#        [1.05409255, 0.66666667, 1.05409255],
#        [1.05409255, 0.66666667, 1.05409255]])