从二维等距网格生成距离矩阵：求更Pythonic的实现方案

我有一个包含N个点的n×n二维等距网格，希望生成一个N×N的距离矩阵D，其中元素D_ij表示点i与点j之间的欧氏距离。例如，网格[[a b],[c d]]对应的距离矩阵应为[[0 1 1 1.414],[1 0 1.414 1],[1 1.414 0 1],[1.414 1 1 0]]。

我目前用暴力循环的方式实现了需求，代码如下：

import numpy as np

N=2**2
n=int(np.sqrt(N))

row = np.arange(0,n)
grid = np.tile(row,n)

#East matrix. For each point, the matrix yields the distance in the x direction(east) to the other points. The first row covers each point in the first row in the grid, and so on.

D_E = np.zeros((N,N))
for i in range(N):
    for j in range(N):
        
        D_E[i,j] = grid[j]-grid[i]
#Same logic as above, but this time in the y_direction
D_S = np.zeros((N,N))

row2 = np.arange(0,n)
grid2 = np.repeat(row2,n)
for i in range(N):
    for j in range(N):
        D_S[i,j] = grid2[j]-grid2[i]

D = np.sqrt(D_E**2 + D_S**2)

请问是否有更Pythonic的实现方法？

补充说明：我的最终目标是利用该距离矩阵生成网格的N×N协方差矩阵，因为协方差矩阵的每个元素由输入网格点间欧氏距离的协方差函数决定。

当然有，利用numpy的向量化操作和广播机制可以完全避免嵌套循环，让代码更简洁、高效，这也是Pythonic的核心体现。

方法一：基于广播的坐标差计算 #

import numpy as np

n = 2
N = n ** 2

# 生成每个点的x、y坐标
x = np.tile(np.arange(n), n)  # 横向重复的x坐标：[0,1,0,1]（n=2时）
y = np.repeat(np.arange(n), n)  # 纵向重复的y坐标：[0,0,1,1]（n=2时）

# 利用广播计算所有点对的x、y坐标差
dx = x[np.newaxis, :] - x[:, np.newaxis]
dy = y[np.newaxis, :] - y[:, np.newaxis]

# 计算欧氏距离矩阵
D = np.sqrt(dx ** 2 + dy ** 2)

方法二：用numpy内置函数简化操作 #

可以用np.subtract.outer直接生成差矩阵，再用np.hypot快速计算欧氏距离（比手动开平方更高效）：

import numpy as np

n = 2
N = n ** 2

x = np.tile(np.arange(n), n)
y = np.repeat(np.arange(n), n)

# 生成所有点对的坐标差矩阵
dx = np.subtract.outer(x, x)
dy = np.subtract.outer(y, y)

# 计算欧氏距离
D = np.hypot(dx, dy)

为什么这更Pythonic？ #

• 无冗余循环：完全抛弃嵌套for循环，利用numpy的向量化操作，代码行数大幅减少，执行效率提升明显（n越大，优势越显著）
• 语义清晰：代码逻辑直观，从生成坐标到计算距离的流程一目了然
• 利用内置函数：np.subtract.outer和np.hypot都是numpy针对这类场景优化的函数，既简洁又高效

内容的提问来源于stack exchange，提问作者Joe