Eigen+MKL集群环境下大矩阵乘法未调用MKL GEMM接口的问题排查

看起来你已经把MKL和Eigen的环境搭得有模有样了——ldd能识别MKL库，小的BLAS操作（比如DGEMV、DPOTRF）也能正常触发MKL调用，但偏偏最关键的大矩阵乘法DGEMM没动静，确实挺闹心的。我来帮你捋捋可能的问题点和解决办法：

1. 先排查Eigen表达式模板的“小脾气” #

Eigen的表达式模板会延迟计算，像转置、切片这类非连续的矩阵表达式，可能会让Eigen放弃调用MKL的DGEMM，转而用自己的实现。你瓶颈处的代码是：

S_mean.noalias()=X_alltoall*X_alltoall.transpose();

这里X_alltoall是列优先（ColMajor），它的转置X_alltoall.transpose()是行优先（RowMajor）的表达式，而非实际的连续矩阵。Eigen可能因为存储顺序的不匹配，没触发MKL的DGEMM接口。

测试方案：
把转置操作先“落地”成连续矩阵再相乘，比如：

MatrixXd Xt = X_alltoall.transpose().eval(); // eval()强制生成连续的转置矩阵
S_mean.noalias() = X_alltoall * Xt;

或者在你的小测试代码里，把C.noalias()=B.transpose()*B;改成：

MatrixType Bt = B.transpose().eval();
C.noalias() = Bt * B;

重新编译运行后，再看MKL_VERBOSE里有没有DGEMM的输出。

2. 检查Eigen的矩阵乘积阈值 #

Eigen默认对小矩阵会用自己的优化实现，只有矩阵尺寸超过某个阈值才会调用外部BLAS/MKL。虽然你的10000x10000矩阵肯定够大，但某些编译选项（比如-march=native）可能让Eigen调整了这个阈值，或者误判矩阵规模。

解决办法：
在所有Eigen头文件包含之前，手动设置阈值，强制大矩阵调用MKL。比如在代码最开头加：

#define EIGEN_DEFAULT_MATRIX_MATRIX_MULTIPLY_THRESHOLD 64 // 只要矩阵≥64x64就调用BLAS/MKL
#define EIGEN_USE_MKL_ALL
#include <Eigen/Dense>

也可以在Makefile的CXXFLAGS里加-DEIGEN_DEFAULT_MATRIX_MULTIPLY_THRESHOLD=64，确保所有源文件都能生效。

3. 确认编译宏的全局有效性 #

你在Makefile里加了-DEIGEN_USE_MKL_ALL，但要确保所有源文件在包含Eigen头文件时都能拿到这个宏。比如有没有某个cpp文件里，在包含Eigen之前自己定义了冲突的宏（比如#define EIGEN_DONT_USE_MKL）？

验证方式：
随便找一个源文件，在包含Eigen之前加一行代码：

#ifdef EIGEN_USE_MKL_ALL 
#error "MKL宏已定义" 
#endif

如果编译时报错，说明宏生效了；如果没报错，那就是这个源文件没拿到宏，得检查Makefile的编译规则是否覆盖了所有cpp文件。

4. 排除MPI编译器的潜在干扰 #

你用mpicxx作为链接器，某些MPI编译器可能会自带默认编译选项，或者修改环境变量，间接影响Eigen对MKL的调用。

测试方案：
用纯g++编译你的小测试代码，避开MPI：

g++ -std=c++20 -O3 -DEIGEN_USE_MKL_ALL -I../../LIBRARIES/eigen-3.4.0 -I${MKLROOT}/include -L${MKLROOT}/lib/intel64 -Wl,-rpath,${MKLROOT}/lib/intel64 -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl test.cpp -o test.out

运行MKL_VERBOSE=1 ./test.out，看有没有DGEMM输出。如果这个能触发，那问题就出在mpicxx的编译选项上，可以对比mpicxx -show和g++的默认选项，排查冲突点。

5. 验证MKL DGEMM符号是否真的被链接 #

有时候ldd能看到MKL库，但可执行文件里可能没实际引用DGEMM符号，Eigen自然没法调用。

验证命令：
用nm工具检查可执行文件：

nm main.out | grep -i dgemm

如果输出里有U dgemm_（U表示未定义，会在运行时从MKL库加载），或者有具体的符号地址，说明链接是对的；如果完全没输出，那说明Eigen根本没生成调用DGEMM的代码，得回到前面的宏定义和编译选项排查。

6. 暂时关闭Eigen的本地向量优化 #

你加了-march=native，Eigen会自动启用AVX/AVX512等本地向量指令的优化实现，有时候这些实现的性能和MKL相当，Eigen可能会优先用自己的。可以试试强制Eigen不用这些指令：

在Makefile的CXXFLAGS里加：

-DEIGEN_DONT_USE_AVX -DEIGEN_DONT_USE_AVX2 -DEIGEN_DONT_USE_AVX512

重新编译运行后，再看MKL_VERBOSE有没有DGEMM输出。

最后再提个小细节 #

你设置了MKL_NUM_THREADS=1来避免和MPI冲突，这个是对的，但要确保这个环境变量真的被MPI进程继承了。可以在代码里加一行输出验证：

#include <cstdlib>
// ...
std::cout << "MKL_NUM_THREADS: " << getenv("MKL_NUM_THREADS") << std::endl;

确认每个MPI进程都能读到这个值。

内容来源于stack exchange