在MKL/CBLAS中如何高效实现逐元素向量与矩阵相乘运算？

首先明确你的运算本质：目标结果等价于原3×3矩阵的每一列分别和行向量对应位置的标量做逐元素乘法，不需要做广播扩维操作，MKL/CBLAS中有两种远优于你当前方案的实现方式：

方案1：逐列调用cblas_Xscal（小矩阵最优方案） #

你的运算没有额外加法操作，完全可以通过BLAS level 1的scal（向量标量乘）函数实现，不需要额外分配内存存储扩展后的向量，访存效率最高。
以双精度浮点数、列优先存储矩阵为例，参考代码逻辑：

// 入参定义：
// m = 3：矩阵行数，n = 3：矩阵列数
// A：指向原矩阵首地址，运算后直接覆盖为结果
// v：指向行向量首地址，长度为n
// lda = 3：矩阵的leading dimension（列优先下为行数，行优先下为列数）
for (int col_idx = 0; col_idx < n; col_idx++) {
    // 对第col_idx列的所有元素乘以v[col_idx]
    cblas_dscal(m, v[col_idx], A + col_idx * lda, 1);
}

如果你的矩阵是行优先存储，只要改为遍历每一行调用cblas_dscal，将步长参数改为lda即可。

方案2：调用GEMM接口（大矩阵最优方案） #

如果矩阵尺寸较大（≥1000×1000），可以将运算转换为矩阵乘形式：目标矩阵 C = A * diag(v)，其中diag(v)是对角元为行向量元素的对角矩阵，直接调用BLAS level 3的cblas_dgemm接口。MKL对大矩阵GEMM做了极致的SIMD指令优化、缓存分层优化和多核心并行优化，运算效率远高于手动循环。
参考调用逻辑：

cblas_dgemm(
    CblasColMajor, CblasNoTrans, CblasNoTrans,
    m, n, n,
    1.0, A, lda,
    diag_v, lda, // diag_v为对角矩阵，非对角元置0
    0.0, C, lda
);

如果不想额外存储完整对角矩阵，可以搭配MKL稀疏矩阵乘法接口cblas_dgemm_sparse，将diag_v存储为稀疏对角格式，进一步降低内存开销。

性能对比 #

• 小尺寸矩阵（<100×100）：方案1性能更高，接口调用开销远低于GEMM的启动开销
• 大尺寸矩阵（≥1000×1000）：方案2性能更高，MKL的底层优化收益完全覆盖调用开销
• 两种方案的内存开销、运算效率都远高于现有扩维逐乘的实现，大矩阵场景下内存开销差距可达数量级。

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