Использование BLAS/LAPACK из Eigen

Начиная с версии Eigen 3.3 и более поздних, любые совместимые с F77 библиотеки BLAS или LAPACK могут использоваться в качестве бэкендов для операций с плотными матрицами и разложениями плотных матриц. Например, можно использовать Intel® MKL, фреймворк Apple Accelerate на OSX, OpenBLAS, Netlib LAPACK и т.д.

Не пропустите эту страницу для дальнейшего обсуждения использования Intel® MKL (также включает VML, PARDISO и т.д.).

Для использования внешней библиотеки BLAS и/или LAPACK необходимо связать собственное приложение с соответствующими библиотеками и их зависимостями. Для LAPACK необходимо также связаться со стандартной библиотекой Lapacke, которая используется в качестве удобного слоя между кодом Eigen на C++ и интерфейсом LAPACK F77. Затем необходимо активировать их использование, определив одну или несколько следующих макросов (до включения заголовков Eigen):

Примечание

Для пользователей Mac, для использования версии lapack, поставляемой с фреймворком Accelerate, также требуется библиотека lapacke. С помощью MacPorts это так же просто, как:

sudo port install lapack

и затем использовать следующие флаги линковки: -framework Accelerate /opt/local/lib/lapack/liblapacke.dylib

При этом ряд алгоритмов Eigen молча заменяется вызовами подпрограмм BLAS или LAPACK. Эти замены применяются только для динамических или достаточно больших объектов с одним из следующих четырех стандартных типов скаляров: float, double, complex<float>, и complex<double>. Операции с другими типами скаляров или смешивание вещественных и комплексных чисел будут продолжать использовать встроенные алгоритмы.

Перечень функциональности Eigen, который может быть заменён, представлен в таблице ниже.

m1*m2.transpose();
m1.selfadjointView<Lower>()*m2;
m1*m2.triangularView<Upper>();
m1.selfadjointView<Lower>().rankUpdate(m2,1.0);

?gemm
?symm/?hemm
?trmm
dsyrk/ssyrk

m1.adjoint()*b;
m1.selfadjointView<Lower>()*b;
m1.triangularView<Upper>()*b;

?gemv
?symv/?hemv
?trmv

v1 = m1.lu().solve(v2);

?getrf

v1 = m2.selfadjointView<Upper>().llt().solve(v2);

?potrf

m1.householderQr();
m1.colPivHouseholderQr();

?geqrf
?geqp3

JacobiSVD<MatrixXd> svd;
svd.compute(m1, ComputeThinV);

?gesvd

EigenSolver<MatrixXd> es(m1);
ComplexEigenSolver<MatrixXcd> ces(m1);
SelfAdjointEigenSolver<MatrixXd> saes(m1+m1.transpose());
GeneralizedSelfAdjointEigenSolver<MatrixXd>
    gsaes(m1+m1.transpose(),m2+m2.transpose());

?gees
?gees
?syev/?heev
?syev/?heev,
?potrf

RealSchur<MatrixXd> schurR(m1);
ComplexSchur<MatrixXcd> schurC(m1);

?gees

В примерах m1 и m2 — плотные матрицы, а v1 и v2 — плотные векторы.