bug #173: remove dependency to MKL for LAPACKe backend.

This commit is contained in:
Gael Guennebaud 2016-07-25 17:55:07 +02:00
parent 34b483e25d
commit 4d54e3dd33
10 changed files with 16349 additions and 42 deletions

View File

@ -40,10 +40,10 @@ namespace Eigen {
namespace internal { namespace internal {
template<typename Scalar> struct mkl_llt; template<typename Scalar> struct lapacke_llt;
#define EIGEN_MKL_LLT(EIGTYPE, MKLTYPE, MKLPREFIX) \ #define EIGEN_MKL_LLT(EIGTYPE, BLASTYPE, MKLPREFIX) \
template<> struct mkl_llt<EIGTYPE> \ template<> struct lapacke_llt<EIGTYPE> \
{ \ { \
template<typename MatrixType> \ template<typename MatrixType> \
static inline Index potrf(MatrixType& m, char uplo) \ static inline Index potrf(MatrixType& m, char uplo) \
@ -59,7 +59,7 @@ template<> struct mkl_llt<EIGTYPE> \
a = &(m.coeffRef(0,0)); \ a = &(m.coeffRef(0,0)); \
lda = convert_index<lapack_int>(m.outerStride()); \ lda = convert_index<lapack_int>(m.outerStride()); \
\ \
info = LAPACKE_##MKLPREFIX##potrf( matrix_order, uplo, size, (MKLTYPE*)a, lda ); \ info = LAPACKE_##MKLPREFIX##potrf( matrix_order, uplo, size, (BLASTYPE*)a, lda ); \
info = (info==0) ? -1 : info>0 ? info-1 : size; \ info = (info==0) ? -1 : info>0 ? info-1 : size; \
return info; \ return info; \
} \ } \
@ -69,7 +69,7 @@ template<> struct llt_inplace<EIGTYPE, Lower> \
template<typename MatrixType> \ template<typename MatrixType> \
static Index blocked(MatrixType& m) \ static Index blocked(MatrixType& m) \
{ \ { \
return mkl_llt<EIGTYPE>::potrf(m, 'L'); \ return lapacke_llt<EIGTYPE>::potrf(m, 'L'); \
} \ } \
template<typename MatrixType, typename VectorType> \ template<typename MatrixType, typename VectorType> \
static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \ static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \
@ -80,7 +80,7 @@ template<> struct llt_inplace<EIGTYPE, Upper> \
template<typename MatrixType> \ template<typename MatrixType> \
static Index blocked(MatrixType& m) \ static Index blocked(MatrixType& m) \
{ \ { \
return mkl_llt<EIGTYPE>::potrf(m, 'U'); \ return lapacke_llt<EIGTYPE>::potrf(m, 'U'); \
} \ } \
template<typename MatrixType, typename VectorType> \ template<typename MatrixType, typename VectorType> \
static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \ static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \
@ -92,8 +92,8 @@ template<> struct llt_inplace<EIGTYPE, Upper> \
EIGEN_MKL_LLT(double, double, d) EIGEN_MKL_LLT(double, double, d)
EIGEN_MKL_LLT(float, float, s) EIGEN_MKL_LLT(float, float, s)
EIGEN_MKL_LLT(dcomplex, MKL_Complex16, z) EIGEN_MKL_LLT(dcomplex, lapack_complex_double, z)
EIGEN_MKL_LLT(scomplex, MKL_Complex8, c) EIGEN_MKL_LLT(scomplex, lapack_complex_float, c)
} // end namespace internal } // end namespace internal

View File

@ -49,7 +49,7 @@
#define EIGEN_USE_LAPACKE #define EIGEN_USE_LAPACKE
#endif #endif
#if defined(EIGEN_USE_LAPACKE) || defined(EIGEN_USE_MKL_VML) #if defined(EIGEN_USE_MKL_VML)
#define EIGEN_USE_MKL #define EIGEN_USE_MKL
#endif #endif
@ -72,7 +72,6 @@
#endif #endif
#if defined EIGEN_USE_MKL #if defined EIGEN_USE_MKL
#include <mkl_lapacke.h>
#define EIGEN_MKL_VML_THRESHOLD 128 #define EIGEN_MKL_VML_THRESHOLD 128
/* MKL_DOMAIN_BLAS, etc are defined only in 10.3 update 7 */ /* MKL_DOMAIN_BLAS, etc are defined only in 10.3 update 7 */
@ -125,4 +124,8 @@ typedef int BlasIndex;
#include "../../misc/blas.h" #include "../../misc/blas.h"
#endif #endif
#ifdef EIGEN_USE_LAPACKE
#include "../../misc/lapacke.h"
#endif
#endif // EIGEN_MKL_SUPPORT_H #endif // EIGEN_MKL_SUPPORT_H

View File

@ -83,10 +83,10 @@ ComplexSchur<Matrix<EIGTYPE, Dynamic, Dynamic, EIGCOLROW> >::compute(const Eigen
\ \
} }
EIGEN_MKL_SCHUR_COMPLEX(dcomplex, MKL_Complex16, z, Z, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_SCHUR_COMPLEX(dcomplex, lapack_complex_double, z, Z, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_SCHUR_COMPLEX(scomplex, MKL_Complex8, c, C, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_SCHUR_COMPLEX(scomplex, lapack_complex_float, c, C, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_SCHUR_COMPLEX(dcomplex, MKL_Complex16, z, Z, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_SCHUR_COMPLEX(dcomplex, lapack_complex_double, z, Z, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_SCHUR_COMPLEX(scomplex, MKL_Complex8, c, C, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_SCHUR_COMPLEX(scomplex, lapack_complex_float, c, C, RowMajor, LAPACK_ROW_MAJOR)
} // end namespace Eigen } // end namespace Eigen

View File

@ -77,15 +77,15 @@ SelfAdjointEigenSolver<Matrix<EIGTYPE, Dynamic, Dynamic, EIGCOLROW> >::compute(c
} }
EIGEN_MKL_EIG_SELFADJ(double, double, double, dsyev, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_EIG_SELFADJ(double, double, double, dsyev, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_EIG_SELFADJ(float, float, float, ssyev, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_EIG_SELFADJ(float, float, float, ssyev, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_EIG_SELFADJ(dcomplex, MKL_Complex16, double, zheev, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_EIG_SELFADJ(dcomplex, lapack_complex_double, double, zheev, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_EIG_SELFADJ(scomplex, MKL_Complex8, float, cheev, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_EIG_SELFADJ(scomplex, lapack_complex_float, float, cheev, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_EIG_SELFADJ(double, double, double, dsyev, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_EIG_SELFADJ(double, double, double, dsyev, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_EIG_SELFADJ(float, float, float, ssyev, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_EIG_SELFADJ(float, float, float, ssyev, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_EIG_SELFADJ(dcomplex, MKL_Complex16, double, zheev, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_EIG_SELFADJ(dcomplex, lapack_complex_double, double, zheev, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_EIG_SELFADJ(scomplex, MKL_Complex8, float, cheev, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_EIG_SELFADJ(scomplex, lapack_complex_float, float, cheev, RowMajor, LAPACK_ROW_MAJOR)
} // end namespace Eigen } // end namespace Eigen

View File

@ -75,8 +75,8 @@ struct partial_lu_impl<EIGTYPE, StorageOrder, lapack_int> \
EIGEN_MKL_LU_PARTPIV(double, double, d) EIGEN_MKL_LU_PARTPIV(double, double, d)
EIGEN_MKL_LU_PARTPIV(float, float, s) EIGEN_MKL_LU_PARTPIV(float, float, s)
EIGEN_MKL_LU_PARTPIV(dcomplex, MKL_Complex16, z) EIGEN_MKL_LU_PARTPIV(dcomplex, lapack_complex_double, z)
EIGEN_MKL_LU_PARTPIV(scomplex, MKL_Complex8, c) EIGEN_MKL_LU_PARTPIV(scomplex, lapack_complex_float, c)
} // end namespace internal } // end namespace internal

View File

@ -86,13 +86,13 @@ ColPivHouseholderQR<Matrix<EIGTYPE, Dynamic, Dynamic, EIGCOLROW, Dynamic, Dynami
EIGEN_MKL_QR_COLPIV(double, double, d, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_QR_COLPIV(double, double, d, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_QR_COLPIV(float, float, s, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_QR_COLPIV(float, float, s, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_QR_COLPIV(dcomplex, MKL_Complex16, z, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_QR_COLPIV(dcomplex, lapack_complex_double, z, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_QR_COLPIV(scomplex, MKL_Complex8, c, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_QR_COLPIV(scomplex, lapack_complex_float, c, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_QR_COLPIV(double, double, d, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_QR_COLPIV(double, double, d, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_QR_COLPIV(float, float, s, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_QR_COLPIV(float, float, s, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_QR_COLPIV(dcomplex, MKL_Complex16, z, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_QR_COLPIV(dcomplex, lapack_complex_double, z, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_QR_COLPIV(scomplex, MKL_Complex8, c, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_QR_COLPIV(scomplex, lapack_complex_float, c, RowMajor, LAPACK_ROW_MAJOR)
} // end namespace Eigen } // end namespace Eigen

View File

@ -60,8 +60,8 @@ struct householder_qr_inplace_blocked<MatrixQR, HCoeffs, EIGTYPE, true> \
EIGEN_MKL_QR_NOPIV(double, double, d) EIGEN_MKL_QR_NOPIV(double, double, d)
EIGEN_MKL_QR_NOPIV(float, float, s) EIGEN_MKL_QR_NOPIV(float, float, s)
EIGEN_MKL_QR_NOPIV(dcomplex, MKL_Complex16, z) EIGEN_MKL_QR_NOPIV(dcomplex, lapack_complex_double, z)
EIGEN_MKL_QR_NOPIV(scomplex, MKL_Complex8, c) EIGEN_MKL_QR_NOPIV(scomplex, lapack_complex_float, c)
} // end namespace internal } // end namespace internal

View File

@ -52,40 +52,40 @@ JacobiSVD<Matrix<EIGTYPE, Dynamic, Dynamic, EIGCOLROW, Dynamic, Dynamic>, ColPiv
/*const RealScalar precision = RealScalar(2) * NumTraits<Scalar>::epsilon();*/ \ /*const RealScalar precision = RealScalar(2) * NumTraits<Scalar>::epsilon();*/ \
m_nonzeroSingularValues = m_diagSize; \ m_nonzeroSingularValues = m_diagSize; \
\ \
lapack_int lda = matrix.outerStride(), ldu, ldvt; \ lapack_int lda = internal::convert_index<lapack_int>(matrix.outerStride()), ldu, ldvt; \
lapack_int matrix_order = MKLCOLROW; \ lapack_int matrix_order = MKLCOLROW; \
char jobu, jobvt; \ char jobu, jobvt; \
MKLTYPE *u, *vt, dummy; \ MKLTYPE *u, *vt, dummy; \
jobu = (m_computeFullU) ? 'A' : (m_computeThinU) ? 'S' : 'N'; \ jobu = (m_computeFullU) ? 'A' : (m_computeThinU) ? 'S' : 'N'; \
jobvt = (m_computeFullV) ? 'A' : (m_computeThinV) ? 'S' : 'N'; \ jobvt = (m_computeFullV) ? 'A' : (m_computeThinV) ? 'S' : 'N'; \
if (computeU()) { \ if (computeU()) { \
ldu = m_matrixU.outerStride(); \ ldu = internal::convert_index<lapack_int>(m_matrixU.outerStride()); \
u = (MKLTYPE*)m_matrixU.data(); \ u = (MKLTYPE*)m_matrixU.data(); \
} else { ldu=1; u=&dummy; }\ } else { ldu=1; u=&dummy; }\
MatrixType localV; \ MatrixType localV; \
ldvt = (m_computeFullV) ? m_cols : (m_computeThinV) ? m_diagSize : 1; \ ldvt = (m_computeFullV) ? internal::convert_index<lapack_int>(m_cols) : (m_computeThinV) ? internal::convert_index<lapack_int>(m_diagSize) : 1; \
if (computeV()) { \ if (computeV()) { \
localV.resize(ldvt, m_cols); \ localV.resize(ldvt, m_cols); \
vt = (MKLTYPE*)localV.data(); \ vt = (MKLTYPE*)localV.data(); \
} else { ldvt=1; vt=&dummy; }\ } else { ldvt=1; vt=&dummy; }\
Matrix<MKLRTYPE, Dynamic, Dynamic> superb; superb.resize(m_diagSize, 1); \ Matrix<MKLRTYPE, Dynamic, Dynamic> superb; superb.resize(m_diagSize, 1); \
MatrixType m_temp; m_temp = matrix; \ MatrixType m_temp; m_temp = matrix; \
LAPACKE_##MKLPREFIX##gesvd( matrix_order, jobu, jobvt, m_rows, m_cols, (MKLTYPE*)m_temp.data(), lda, (MKLRTYPE*)m_singularValues.data(), u, ldu, vt, ldvt, superb.data()); \ LAPACKE_##MKLPREFIX##gesvd( matrix_order, jobu, jobvt, internal::convert_index<lapack_int>(m_rows), internal::convert_index<lapack_int>(m_cols), (MKLTYPE*)m_temp.data(), lda, (MKLRTYPE*)m_singularValues.data(), u, ldu, vt, ldvt, superb.data()); \
if (computeV()) m_matrixV = localV.adjoint(); \ if (computeV()) m_matrixV = localV.adjoint(); \
/* for(int i=0;i<m_diagSize;i++) if (m_singularValues.coeffRef(i) < precision) { m_nonzeroSingularValues--; m_singularValues.coeffRef(i)=RealScalar(0);}*/ \ /* for(int i=0;i<m_diagSize;i++) if (m_singularValues.coeffRef(i) < precision) { m_nonzeroSingularValues--; m_singularValues.coeffRef(i)=RealScalar(0);}*/ \
m_isInitialized = true; \ m_isInitialized = true; \
return *this; \ return *this; \
} }
EIGEN_MKL_SVD(double, double, double, d, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_SVD(double, double, double, d, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_SVD(float, float, float , s, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_SVD(float, float, float , s, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_SVD(dcomplex, MKL_Complex16, double, z, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_SVD(dcomplex, lapack_complex_double, double, z, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_SVD(scomplex, MKL_Complex8, float , c, ColMajor, LAPACK_COL_MAJOR) EIGEN_MKL_SVD(scomplex, lapack_complex_float, float , c, ColMajor, LAPACK_COL_MAJOR)
EIGEN_MKL_SVD(double, double, double, d, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_SVD(double, double, double, d, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_SVD(float, float, float , s, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_SVD(float, float, float , s, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_SVD(dcomplex, MKL_Complex16, double, z, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_SVD(dcomplex, lapack_complex_double, double, z, RowMajor, LAPACK_ROW_MAJOR)
EIGEN_MKL_SVD(scomplex, MKL_Complex8, float , c, RowMajor, LAPACK_ROW_MAJOR) EIGEN_MKL_SVD(scomplex, lapack_complex_float, float , c, RowMajor, LAPACK_ROW_MAJOR)
} // end namespace Eigen } // end namespace Eigen

16287
Eigen/src/misc/lapacke.h Normal file

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,17 @@
#ifndef LAPACK_HEADER_INCLUDED
#define LAPACK_HEADER_INCLUDED
#ifndef LAPACK_GLOBAL
#if defined(LAPACK_GLOBAL_PATTERN_LC) || defined(ADD_)
#define LAPACK_GLOBAL(lcname,UCNAME) lcname##_
#elif defined(LAPACK_GLOBAL_PATTERN_UC) || defined(UPPER)
#define LAPACK_GLOBAL(lcname,UCNAME) UCNAME
#elif defined(LAPACK_GLOBAL_PATTERN_MC) || defined(NOCHANGE)
#define LAPACK_GLOBAL(lcname,UCNAME) lcname
#else
#define LAPACK_GLOBAL(lcname,UCNAME) lcname##_
#endif
#endif
#endif