使用 Eigen + Intel MKL + Pardiso
Use Eigen + Intel MKL + Pardiso
我正在尝试使用 Eigen 对 MKL 和 Pardiso 的支持(参见下面的示例)。我已经使用 Intel link line advisor 来提供编译器选项,但我尝试的一切都不成功。特别是,编译:
$ icpc -I${HOME}/src/eigen -DEIGEN_USE_MKL_ALL -DMKL_ILP64 -I${MKLROOT}/include main.cpp -L${MKLROOT}/lib/intel64 -lmkl_intel_ilp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl
导致以下类型错误:
~/src/eigen/Eigen/src/PardisoSupport/PardisoSupport.h(50): error: argument of type "int *" is incompatible with parameter of type "const long long *"
::pardiso(pt, &maxfct, &mnum, &type, &phase, &n, a, ia, ja, perm, &nrhs, iparm, &msglvl, b, x, &error);
在没有 64 位支持的情况下编译时,使用:
$ icpc -I${HOME}/src/eigen -DEIGEN_USE_MKL_ALL -I${MKLROOT}/include main.cpp -L${MKLROOT}/lib/intel64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl
结果:
main.cpp:(.text+0x759): undefined reference to `pardiso'
我也试过用 OpenMP 编译,同样失败。
我应该如何编译?
(我使用的是 Intel 2016.3.210 和最新的 Eigen,3.3.3)
示例 (main.cpp)
本例求解泊松方程,具有两个狄利克雷边界条件 (u_{0} = u_{N+1} = 0)
#include <iostream>
#include <iomanip>
#include <math.h>
#include <Eigen/Dense>
#include <Eigen/PardisoSupport>
int main()
{
typedef Eigen::SparseMatrix<double> SpMat;
typedef Eigen::Triplet <double> Trip;
typedef Eigen::PardisoLU <SpMat > Solver;
size_t N = 11;
SpMat A(N,N);
Eigen::VectorXd f(N);
f *= 0.0;
f((N-1)/2) = 1.0;
std::vector<Trip> Atr;
Atr.push_back(Trip(0,0,+2.0));
Atr.push_back(Trip(0,1,-1.0));
for ( size_t i=1; i<N-1; ++i ) {
Atr.push_back(Trip(i,i-1,-1.0));
Atr.push_back(Trip(i,i ,+2.0));
Atr.push_back(Trip(i,i+1,-1.0));
}
Atr.push_back(Trip(N-1,N-2,-1.0));
Atr.push_back(Trip(N-1,N-1,+2.0));
A.setFromTriplets(Atr.begin(),Atr.end());
Solver solver;
solver.analyzePattern(A);
solver.factorize(A);
Eigen::VectorXd u = solver.solve(f);
return 0;
}
原来eigen3 documentation很清楚这一点:
Using Intel MKL through Eigen is easy:
...
- on a 64bits system, you must use the LP64 interface (not the ILP64 one)
我现在已经成功编译
icpc -I${HOME}/src/eigen -DEIGEN_USE_MKL_ALL -DMKL_LP64 -I${MKLROOT}/include main.cpp -L${MKLROOT}/lib/intel64 -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl
我正在尝试使用 Eigen 对 MKL 和 Pardiso 的支持(参见下面的示例)。我已经使用 Intel link line advisor 来提供编译器选项,但我尝试的一切都不成功。特别是,编译:
$ icpc -I${HOME}/src/eigen -DEIGEN_USE_MKL_ALL -DMKL_ILP64 -I${MKLROOT}/include main.cpp -L${MKLROOT}/lib/intel64 -lmkl_intel_ilp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl
导致以下类型错误:
~/src/eigen/Eigen/src/PardisoSupport/PardisoSupport.h(50): error: argument of type "int *" is incompatible with parameter of type "const long long *" ::pardiso(pt, &maxfct, &mnum, &type, &phase, &n, a, ia, ja, perm, &nrhs, iparm, &msglvl, b, x, &error);
在没有 64 位支持的情况下编译时,使用:
$ icpc -I${HOME}/src/eigen -DEIGEN_USE_MKL_ALL -I${MKLROOT}/include main.cpp -L${MKLROOT}/lib/intel64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl
结果:
main.cpp:(.text+0x759): undefined reference to `pardiso'
我也试过用 OpenMP 编译,同样失败。
我应该如何编译?
(我使用的是 Intel 2016.3.210 和最新的 Eigen,3.3.3)
示例 (main.cpp)
本例求解泊松方程,具有两个狄利克雷边界条件 (u_{0} = u_{N+1} = 0)
#include <iostream>
#include <iomanip>
#include <math.h>
#include <Eigen/Dense>
#include <Eigen/PardisoSupport>
int main()
{
typedef Eigen::SparseMatrix<double> SpMat;
typedef Eigen::Triplet <double> Trip;
typedef Eigen::PardisoLU <SpMat > Solver;
size_t N = 11;
SpMat A(N,N);
Eigen::VectorXd f(N);
f *= 0.0;
f((N-1)/2) = 1.0;
std::vector<Trip> Atr;
Atr.push_back(Trip(0,0,+2.0));
Atr.push_back(Trip(0,1,-1.0));
for ( size_t i=1; i<N-1; ++i ) {
Atr.push_back(Trip(i,i-1,-1.0));
Atr.push_back(Trip(i,i ,+2.0));
Atr.push_back(Trip(i,i+1,-1.0));
}
Atr.push_back(Trip(N-1,N-2,-1.0));
Atr.push_back(Trip(N-1,N-1,+2.0));
A.setFromTriplets(Atr.begin(),Atr.end());
Solver solver;
solver.analyzePattern(A);
solver.factorize(A);
Eigen::VectorXd u = solver.solve(f);
return 0;
}
原来eigen3 documentation很清楚这一点:
Using Intel MKL through Eigen is easy:
...
- on a 64bits system, you must use the LP64 interface (not the ILP64 one)
我现在已经成功编译
icpc -I${HOME}/src/eigen -DEIGEN_USE_MKL_ALL -DMKL_LP64 -I${MKLROOT}/include main.cpp -L${MKLROOT}/lib/intel64 -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl