如何通过MATLAB向量化简化双重循环代码以提升运行效率？

原代码逻辑拆解 #

原代码通过双重循环遍历所有上三角索引对(i,j)（i<j），构造基向量外积矩阵并按权重w累加得到矩阵B。核心逻辑可拆解为：

• 对每个i<j，计算权重w = D1(i,j)/D2(i,j)（若D2(i,j)<=0则w=0）
• 向B中累加w * (e_i - e_j)(e_i - e_j)^T，其中e_k是第k个标准基向量

向量化实现步骤 #

利用矩阵运算直接替代循环，避免逐元素迭代的开销：

1. 构造权重矩阵W
   先计算全局权重，再仅保留上三角有效部分：

N = 500;
D_1 = ones(N,N);
D_2 = ones(N,N);

% 计算初始权重，D2<=0时置0
W = D_1 ./ D_2;
W(D_2 <= 0) = 0;
% 仅保留上三角（i<j）的权重，下三角和对角线置0
tril(W, 0) = 0;

2. 直接计算矩阵B
   根据外积累加的数学特性，B可通过对角矩阵与对称矩阵的组合快速得到：

% 计算对角线元素：每个位置k的权重和 = 所有j>k的W(k,j) + 所有i<k的W(i,k)
diag_vals = sum(W, 2) + sum(W, 1)';
% 构造对角矩阵部分
B_diag = diag(diag_vals);
% 构造非对角线部分：负的权重对称矩阵
B_offdiag = -(W + W');
% 合并得到最终B
B = B_diag + B_offdiag;

等价性验证 #

以N=2为例：

• 原循环得到B = [[1,-1],[-1,1]]
• 向量化代码中W = [[0,1],[0,0]]，diag_vals = [1;1]，B_offdiag = [[0,-1],[-1,0]]，合并后结果与原代码完全一致。

性能优势 #

原双重循环的时间复杂度为O(N²)且存在大量MATLAB循环开销，向量化后利用底层BLAS优化的矩阵运算，运行速度可提升数倍至数十倍，尤其当N较大时（如N=500）效果显著。

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