1c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath/* 2c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Copyright (c) 2011, Intel Corporation. All rights reserved. 3c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 4c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Redistribution and use in source and binary forms, with or without modification, 5c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath are permitted provided that the following conditions are met: 6c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 7c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Redistributions of source code must retain the above copyright notice, this 8c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath list of conditions and the following disclaimer. 9c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Redistributions in binary form must reproduce the above copyright notice, 10c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath this list of conditions and the following disclaimer in the documentation 11c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath and/or other materials provided with the distribution. 12c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Neither the name of Intel Corporation nor the names of its contributors may 13c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath be used to endorse or promote products derived from this software without 14c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath specific prior written permission. 15c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 16c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND 17c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED 18c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 19c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR 20c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 21c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 22c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON 23c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 24c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 25c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 26c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 27c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ******************************************************************************** 282b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang * Content : Eigen bindings to LAPACKe 29c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * LLt decomposition based on LAPACKE_?potrf function. 30c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ******************************************************************************** 31c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath*/ 32c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 332b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang#ifndef EIGEN_LLT_LAPACKE_H 342b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang#define EIGEN_LLT_LAPACKE_H 35c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 36c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace Eigen { 37c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 38c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace internal { 39c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 402b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wangtemplate<typename Scalar> struct lapacke_llt; 41c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 422b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang#define EIGEN_LAPACKE_LLT(EIGTYPE, BLASTYPE, LAPACKE_PREFIX) \ 432b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wangtemplate<> struct lapacke_llt<EIGTYPE> \ 44c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ \ 45c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename MatrixType> \ 462b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang static inline Index potrf(MatrixType& m, char uplo) \ 47c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { \ 48c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath lapack_int matrix_order; \ 49c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath lapack_int size, lda, info, StorageOrder; \ 50c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath EIGTYPE* a; \ 51c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m.rows()==m.cols()); \ 52c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* Set up parameters for ?potrf */ \ 532b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang size = convert_index<lapack_int>(m.rows()); \ 54c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath StorageOrder = MatrixType::Flags&RowMajorBit?RowMajor:ColMajor; \ 55c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath matrix_order = StorageOrder==RowMajor ? LAPACK_ROW_MAJOR : LAPACK_COL_MAJOR; \ 56c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a = &(m.coeffRef(0,0)); \ 572b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang lda = convert_index<lapack_int>(m.outerStride()); \ 58c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath\ 592b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang info = LAPACKE_##LAPACKE_PREFIX##potrf( matrix_order, uplo, size, (BLASTYPE*)a, lda ); \ 60a829215e078ace896f52702caa0c27608f40e3b0Miao Wang info = (info==0) ? -1 : info>0 ? info-1 : size; \ 61c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return info; \ 62c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } \ 63c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; \ 64c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<> struct llt_inplace<EIGTYPE, Lower> \ 65c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ \ 66c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename MatrixType> \ 672b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang static Index blocked(MatrixType& m) \ 68c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { \ 692b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang return lapacke_llt<EIGTYPE>::potrf(m, 'L'); \ 70c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } \ 71c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename MatrixType, typename VectorType> \ 722b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \ 73c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { return Eigen::internal::llt_rank_update_lower(mat, vec, sigma); } \ 74c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; \ 75c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<> struct llt_inplace<EIGTYPE, Upper> \ 76c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ \ 77c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename MatrixType> \ 782b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang static Index blocked(MatrixType& m) \ 79c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { \ 802b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang return lapacke_llt<EIGTYPE>::potrf(m, 'U'); \ 81c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } \ 82c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename MatrixType, typename VectorType> \ 832b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang static Index rankUpdate(MatrixType& mat, const VectorType& vec, const typename MatrixType::RealScalar& sigma) \ 84c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { \ 85c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Transpose<MatrixType> matt(mat); \ 86c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return llt_inplace<EIGTYPE, Lower>::rankUpdate(matt, vec.conjugate(), sigma); \ 87c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } \ 88c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; 89c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 902b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao WangEIGEN_LAPACKE_LLT(double, double, d) 912b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao WangEIGEN_LAPACKE_LLT(float, float, s) 922b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao WangEIGEN_LAPACKE_LLT(dcomplex, lapack_complex_double, z) 932b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao WangEIGEN_LAPACKE_LLT(scomplex, lapack_complex_float, c) 94c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 95c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace internal 96c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 97c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace Eigen 98c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 992b8756b6f1de65d3f8bffab45be6c44ceb7411fcMiao Wang#endif // EIGEN_LLT_LAPACKE_H 100