Eigen  3.3.3
LLT_LAPACKE.h
00001 /*
00002  Copyright (c) 2011, Intel Corporation. All rights reserved.
00003 
00004  Redistribution and use in source and binary forms, with or without modification,
00005  are permitted provided that the following conditions are met:
00006 
00007  * Redistributions of source code must retain the above copyright notice, this
00008    list of conditions and the following disclaimer.
00009  * Redistributions in binary form must reproduce the above copyright notice,
00010    this list of conditions and the following disclaimer in the documentation
00011    and/or other materials provided with the distribution.
00012  * Neither the name of Intel Corporation nor the names of its contributors may
00013    be used to endorse or promote products derived from this software without
00014    specific prior written permission.
00015 
00016  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
00017  ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
00018  WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
00019  DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
00020  ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
00021  (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
00022  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
00023  ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
00024  (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00025  SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00026 
00027  ********************************************************************************
00028  *   Content : Eigen bindings to LAPACKe
00029  *     LLt decomposition based on LAPACKE_?potrf function.
00030  ********************************************************************************
00031 */
00032 
00033 #ifndef EIGEN_LLT_LAPACKE_H
00034 #define EIGEN_LLT_LAPACKE_H
00035 
00036 namespace Eigen { 
00037 
00038 namespace internal {
00039 
00040 template<typename Scalar> struct lapacke_llt;
00041 
00042 #define EIGEN_LAPACKE_LLT(EIGTYPE, BLASTYPE, LAPACKE_PREFIX) \
00043 template<> struct lapacke_llt<EIGTYPE> \
00044 { \
00045   template<typename MatrixType> \
00046   static inline Index potrf(MatrixType& m, char uplo) \
00047   { \
00048     lapack_int matrix_order; \
00049     lapack_int size, lda, info, StorageOrder; \
00050     EIGTYPE* a; \
00051     eigen_assert(m.rows()==m.cols()); \
00052     /* Set up parameters for ?potrf */ \
00053     size = convert_index<lapack_int>(m.rows()); \
00054     StorageOrder = MatrixType::Flags&RowMajorBit?RowMajor:ColMajor; \
00055     matrix_order = StorageOrder==RowMajor ? LAPACK_ROW_MAJOR : LAPACK_COL_MAJOR; \
00056     a = &(m.coeffRef(0,0)); \
00057     lda = convert_index<lapack_int>(m.outerStride()); \
00058 \
00059     info = LAPACKE_##LAPACKE_PREFIX##potrf( matrix_order, uplo, size, (BLASTYPE*)a, lda ); \
00060     info = (info==0) ? -1 : info>0 ? info-1 : size; \
00061     return info; \
00062   } \
00063 }; \
00064 template<> struct llt_inplace<EIGTYPE, Lower> \
00065 { \
00066   template<typename MatrixType> \
00067   static Index blocked(MatrixType& m) \
00068   { \
00069     return lapacke_llt<EIGTYPE>::potrf(m, 'L'); \
00070   } \
00071   template<typename MatrixType, typename VectorType> \
00072   static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \
00073   { return Eigen::internal::llt_rank_update_lower(mat, vec, sigma); } \
00074 }; \
00075 template<> struct llt_inplace<EIGTYPE, Upper> \
00076 { \
00077   template<typename MatrixType> \
00078   static Index blocked(MatrixType& m) \
00079   { \
00080     return lapacke_llt<EIGTYPE>::potrf(m, 'U'); \
00081   } \
00082   template<typename MatrixType, typename VectorType> \
00083   static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \
00084   { \
00085     Transpose<MatrixType> matt(mat); \
00086     return llt_inplace<EIGTYPE, Lower>::rankUpdate(matt, vec.conjugate(), sigma); \
00087   } \
00088 };
00089 
00090 EIGEN_LAPACKE_LLT(double, double, d)
00091 EIGEN_LAPACKE_LLT(float, float, s)
00092 EIGEN_LAPACKE_LLT(dcomplex, lapack_complex_double, z)
00093 EIGEN_LAPACKE_LLT(scomplex, lapack_complex_float, c)
00094 
00095 } // end namespace internal
00096 
00097 } // end namespace Eigen
00098 
00099 #endif // EIGEN_LLT_LAPACKE_H
 All Classes Functions Variables Typedefs Enumerations Enumerator Friends