1*bf2c3715SXin Li // This file is part of Eigen, a lightweight C++ template library 2*bf2c3715SXin Li // for linear algebra. 3*bf2c3715SXin Li // 4*bf2c3715SXin Li // Copyright (C) 2012 Chen-Pang He <[email protected]> 5*bf2c3715SXin Li // 6*bf2c3715SXin Li // This Source Code Form is subject to the terms of the Mozilla 7*bf2c3715SXin Li // Public License v. 2.0. If a copy of the MPL was not distributed 8*bf2c3715SXin Li // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9*bf2c3715SXin Li 10*bf2c3715SXin Li #ifndef EIGEN_SELFADJOINT_PACKED_PRODUCT_H 11*bf2c3715SXin Li #define EIGEN_SELFADJOINT_PACKED_PRODUCT_H 12*bf2c3715SXin Li 13*bf2c3715SXin Li namespace internal { 14*bf2c3715SXin Li 15*bf2c3715SXin Li /* Optimized matrix += alpha * uv' 16*bf2c3715SXin Li * The matrix is in packed form. 17*bf2c3715SXin Li */ 18*bf2c3715SXin Li template<typename Scalar, typename Index, int StorageOrder, int UpLo, bool ConjLhs, bool ConjRhs> 19*bf2c3715SXin Li struct selfadjoint_packed_rank1_update; 20*bf2c3715SXin Li 21*bf2c3715SXin Li template<typename Scalar, typename Index, int UpLo, bool ConjLhs, bool ConjRhs> 22*bf2c3715SXin Li struct selfadjoint_packed_rank1_update<Scalar,Index,ColMajor,UpLo,ConjLhs,ConjRhs> 23*bf2c3715SXin Li { 24*bf2c3715SXin Li typedef typename NumTraits<Scalar>::Real RealScalar; 25*bf2c3715SXin Li static void run(Index size, Scalar* mat, const Scalar* vec, RealScalar alpha) 26*bf2c3715SXin Li { 27*bf2c3715SXin Li typedef Map<const Matrix<Scalar,Dynamic,1> > OtherMap; 28*bf2c3715SXin Li typedef typename conj_expr_if<ConjLhs,OtherMap>::type ConjRhsType; 29*bf2c3715SXin Li conj_if<ConjRhs> cj; 30*bf2c3715SXin Li 31*bf2c3715SXin Li for (Index i=0; i<size; ++i) 32*bf2c3715SXin Li { 33*bf2c3715SXin Li Map<Matrix<Scalar,Dynamic,1> >(mat, UpLo==Lower ? size-i : (i+1)) += alpha * cj(vec[i]) * ConjRhsType(OtherMap(vec+(UpLo==Lower ? i : 0), UpLo==Lower ? size-i : (i+1))); 34*bf2c3715SXin Li //FIXME This should be handled outside. 35*bf2c3715SXin Li mat[UpLo==Lower ? 0 : i] = numext::real(mat[UpLo==Lower ? 0 : i]); 36*bf2c3715SXin Li mat += UpLo==Lower ? size-i : (i+1); 37*bf2c3715SXin Li } 38*bf2c3715SXin Li } 39*bf2c3715SXin Li }; 40*bf2c3715SXin Li 41*bf2c3715SXin Li template<typename Scalar, typename Index, int UpLo, bool ConjLhs, bool ConjRhs> 42*bf2c3715SXin Li struct selfadjoint_packed_rank1_update<Scalar,Index,RowMajor,UpLo,ConjLhs,ConjRhs> 43*bf2c3715SXin Li { 44*bf2c3715SXin Li typedef typename NumTraits<Scalar>::Real RealScalar; 45*bf2c3715SXin Li static void run(Index size, Scalar* mat, const Scalar* vec, RealScalar alpha) 46*bf2c3715SXin Li { 47*bf2c3715SXin Li selfadjoint_packed_rank1_update<Scalar,Index,ColMajor,UpLo==Lower?Upper:Lower,ConjRhs,ConjLhs>::run(size,mat,vec,alpha); 48*bf2c3715SXin Li } 49*bf2c3715SXin Li }; 50*bf2c3715SXin Li 51*bf2c3715SXin Li } // end namespace internal 52*bf2c3715SXin Li 53*bf2c3715SXin Li #endif // EIGEN_SELFADJOINT_PACKED_PRODUCT_H 54