如何无循环实现Numpy矩阵的特定乘法计算以提升效率？

针对你需要对矩阵u的每一列计算二次型u_i^T M u_i的需求，这里提供两种无需显式循环的numpy向量化实现方法，大幅提升大尺寸矩阵场景下的计算效率：

原始循环实现（对照参考） #

import numpy as np

M = np.array([[1,2,3],[3,4,5],[6,7,8]])
u = np.array([[1,2,3],[4,5,7],[2,4,9]])
res = np.zeros((3,))

for i in range(3):
      res[i] = np.matmul(np.matmul(u[:,i].T, M), u[:,i])
# 输出结果: res = array([ 231.,  594., 1957.])

高效向量化方案 #

方案1：使用np.einsum（可读性最强）

einsum通过索引直接描述运算逻辑，精准匹配二次型的计算需求：

res = np.einsum('ji,jk,ki->i', u, M, u)

• 索引说明：ji代表取u的转置（对应列向量u_i^T），jk对应矩阵M，ki对应u的列向量；收缩重复的j、k索引后，保留列维度的i，最终得到每一列的二次型结果。

方案2：矩阵乘法+逐元素求和

利用矩阵乘法和广播机制完成批量计算：

res = (u.T @ M * u.T).sum(axis=1)

• 原理：u.T @ M的第i行是u_i^T @ M，与u.T的第i行（即u_i的转置）逐元素相乘后求和，等价于u_i^T @ M @ u_i的计算结果。

性能优势 #

两种方案均依赖numpy底层优化的C代码执行运算，完全规避了Python显式循环的性能开销，在大矩阵（如维度为1000+）场景下，速度提升可达数十倍甚至上百倍。

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