10 #ifndef EIGEN_SELFADJOINT_PRODUCT_H
11 #define EIGEN_SELFADJOINT_PRODUCT_H
22 template<
typename Scalar,
typename Index,
int UpLo,
bool ConjLhs,
bool ConjRhs>
23 struct selfadjoint_rank1_update<Scalar,Index,
ColMajor,UpLo,ConjLhs,ConjRhs>
25 static void run(Index size, Scalar* mat, Index stride,
const Scalar* vecX,
const Scalar* vecY,
const Scalar& alpha)
27 internal::conj_if<ConjRhs> cj;
28 typedef Map<const Matrix<Scalar,Dynamic,1> > OtherMap;
29 typedef typename internal::conditional<ConjLhs,typename OtherMap::ConjugateReturnType,const OtherMap&>::type ConjLhsType;
30 for (Index i=0; i<size; ++i)
32 Map<Matrix<Scalar,Dynamic,1> >(mat+stride*i+(UpLo==
Lower ? i : 0), (UpLo==
Lower ? size-i : (i+1)))
33 += (alpha * cj(vecY[i])) * ConjLhsType(OtherMap(vecX+(UpLo==
Lower ? i : 0),UpLo==
Lower ? size-i : (i+1)));
38 template<
typename Scalar,
typename Index,
int UpLo,
bool ConjLhs,
bool ConjRhs>
39 struct selfadjoint_rank1_update<Scalar,Index,
RowMajor,UpLo,ConjLhs,ConjRhs>
41 static void run(Index size, Scalar* mat, Index stride,
const Scalar* vecX,
const Scalar* vecY,
const Scalar& alpha)
43 selfadjoint_rank1_update<Scalar,Index,ColMajor,UpLo==Lower?Upper:Lower,ConjRhs,ConjLhs>::run(size,mat,stride,vecY,vecX,alpha);
47 template<
typename MatrixType,
typename OtherType,
int UpLo,
bool OtherIsVector = OtherType::IsVectorAtCompileTime>
48 struct selfadjoint_product_selector;
50 template<
typename MatrixType,
typename OtherType,
int UpLo>
51 struct selfadjoint_product_selector<MatrixType,OtherType,UpLo,true>
53 static void run(MatrixType& mat,
const OtherType& other,
const typename MatrixType::Scalar& alpha)
55 typedef typename MatrixType::Scalar Scalar;
56 typedef internal::blas_traits<OtherType> OtherBlasTraits;
57 typedef typename OtherBlasTraits::DirectLinearAccessType ActualOtherType;
58 typedef typename internal::remove_all<ActualOtherType>::type _ActualOtherType;
59 typename internal::add_const_on_value_type<ActualOtherType>::type actualOther = OtherBlasTraits::extract(other.derived());
61 Scalar actualAlpha = alpha * OtherBlasTraits::extractScalarFactor(other.derived());
65 UseOtherDirectly = _ActualOtherType::InnerStrideAtCompileTime==1
67 internal::gemv_static_vector_if<Scalar,OtherType::SizeAtCompileTime,OtherType::MaxSizeAtCompileTime,!UseOtherDirectly> static_other;
69 ei_declare_aligned_stack_constructed_variable(Scalar, actualOtherPtr, other.size(),
70 (UseOtherDirectly ?
const_cast<Scalar*
>(actualOther.data()) : static_other.data()));
73 Map<typename _ActualOtherType::PlainObject>(actualOtherPtr, actualOther.size()) = actualOther;
75 selfadjoint_rank1_update<Scalar,Index,StorageOrder,UpLo,
76 OtherBlasTraits::NeedToConjugate && NumTraits<Scalar>::IsComplex,
77 (!OtherBlasTraits::NeedToConjugate) && NumTraits<Scalar>::IsComplex>
78 ::run(other.size(), mat.data(), mat.outerStride(), actualOtherPtr, actualOtherPtr, actualAlpha);
82 template<
typename MatrixType,
typename OtherType,
int UpLo>
83 struct selfadjoint_product_selector<MatrixType,OtherType,UpLo,false>
85 static void run(MatrixType& mat,
const OtherType& other,
const typename MatrixType::Scalar& alpha)
87 typedef typename MatrixType::Scalar Scalar;
88 typedef internal::blas_traits<OtherType> OtherBlasTraits;
89 typedef typename OtherBlasTraits::DirectLinearAccessType ActualOtherType;
90 typedef typename internal::remove_all<ActualOtherType>::type _ActualOtherType;
91 typename internal::add_const_on_value_type<ActualOtherType>::type actualOther = OtherBlasTraits::extract(other.derived());
93 Scalar actualAlpha = alpha * OtherBlasTraits::extractScalarFactor(other.derived());
95 enum { IsRowMajor = (internal::traits<MatrixType>::Flags&
RowMajorBit) ? 1 : 0 };
97 internal::general_matrix_matrix_triangular_product<Index,
98 Scalar, _ActualOtherType::Flags&
RowMajorBit ?
RowMajor :
ColMajor, OtherBlasTraits::NeedToConjugate && NumTraits<Scalar>::IsComplex,
99 Scalar, _ActualOtherType::Flags&
RowMajorBit ? ColMajor :
RowMajor, (!OtherBlasTraits::NeedToConjugate) && NumTraits<Scalar>::IsComplex,
101 ::run(mat.cols(), actualOther.cols(),
102 &actualOther.coeffRef(0,0), actualOther.outerStride(), &actualOther.coeffRef(0,0), actualOther.outerStride(),
103 mat.data(), mat.outerStride(), actualAlpha);
109 template<
typename MatrixType,
unsigned int UpLo>
110 template<
typename DerivedU>
112 ::rankUpdate(
const MatrixBase<DerivedU>& u,
const Scalar& alpha)
114 selfadjoint_product_selector<MatrixType,DerivedU,UpLo>::run(_expression().const_cast_derived(), u.derived(), alpha);
121 #endif // EIGEN_SELFADJOINT_PRODUCT_H
SelfAdjointView & rankUpdate(const MatrixBase< DerivedU > &u, const MatrixBase< DerivedV > &v, const Scalar &alpha=Scalar(1))
Definition: Constants.h:320
const unsigned int RowMajorBit
Definition: Constants.h:61
Definition: Constants.h:322
Definition: Constants.h:204