diff options
author | Narayan Kamath <narayan@google.com> | 2012-11-02 10:59:05 +0000 |
---|---|---|
committer | Xiaotao Duan <xiaotao@google.com> | 2012-11-07 14:17:48 -0800 |
commit | c981c48f5bc9aefeffc0bcb0cc3934c2fae179dd (patch) | |
tree | 54d1c7d66098154c1d7c5bd414394ef4cf255810 /Eigen/src/Core/products/SelfadjointRank2Update.h | |
parent | 63f67d748682b46d58be31235a0a2d64d81b998c (diff) | |
download | eigen-jb-mr1.1-release.tar.gz |
Initial import of eigen 3.1.1android-cts-4.4_r4android-cts-4.2_r2android-4.4_r1.2.0.1android-4.4_r1.2android-4.4_r1.1.0.1android-4.4_r1.1android-4.4_r1.0.1android-4.4_r1android-4.4_r0.9android-4.4_r0.8android-4.3_r2.3android-4.3_r2.2android-4.3_r2.1android-4.3_r2android-4.3_r1.1android-4.3_r1android-4.3_r0.9.1android-4.3_r0.9android-4.2.2_r1.2android-4.2.2_r1.1android-4.2.2_r1kitkat-releasekitkat-cts-releasejb-mr2.0-releasejb-mr2-releasejb-mr1.1-releasejb-mr1.1-dev
Added a README.android and a MODULE_LICENSE_MPL2 file.
Added empty Android.mk and CleanSpec.mk to optimize Android build.
Non MPL2 license code is disabled in ./Eigen/src/Core/util/NonMPL2.h.
Trying to include such files will lead to an error.
Change-Id: I0e148b7c3e83999bcc4dfaa5809d33bfac2aac32
Diffstat (limited to 'Eigen/src/Core/products/SelfadjointRank2Update.h')
-rw-r--r-- | Eigen/src/Core/products/SelfadjointRank2Update.h | 93 |
1 files changed, 93 insertions, 0 deletions
diff --git a/Eigen/src/Core/products/SelfadjointRank2Update.h b/Eigen/src/Core/products/SelfadjointRank2Update.h new file mode 100644 index 000000000..57a98cc2d --- /dev/null +++ b/Eigen/src/Core/products/SelfadjointRank2Update.h @@ -0,0 +1,93 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr> +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#ifndef EIGEN_SELFADJOINTRANK2UPTADE_H +#define EIGEN_SELFADJOINTRANK2UPTADE_H + +namespace Eigen { + +namespace internal { + +/* Optimized selfadjoint matrix += alpha * uv' + conj(alpha)*vu' + * It corresponds to the Level2 syr2 BLAS routine + */ + +template<typename Scalar, typename Index, typename UType, typename VType, int UpLo> +struct selfadjoint_rank2_update_selector; + +template<typename Scalar, typename Index, typename UType, typename VType> +struct selfadjoint_rank2_update_selector<Scalar,Index,UType,VType,Lower> +{ + static void run(Scalar* mat, Index stride, const UType& u, const VType& v, Scalar alpha) + { + const Index size = u.size(); + for (Index i=0; i<size; ++i) + { + Map<Matrix<Scalar,Dynamic,1> >(mat+stride*i+i, size-i) += + (conj(alpha) * conj(u.coeff(i))) * v.tail(size-i) + + (alpha * conj(v.coeff(i))) * u.tail(size-i); + } + } +}; + +template<typename Scalar, typename Index, typename UType, typename VType> +struct selfadjoint_rank2_update_selector<Scalar,Index,UType,VType,Upper> +{ + static void run(Scalar* mat, Index stride, const UType& u, const VType& v, Scalar alpha) + { + const Index size = u.size(); + for (Index i=0; i<size; ++i) + Map<Matrix<Scalar,Dynamic,1> >(mat+stride*i, i+1) += + (conj(alpha) * conj(u.coeff(i))) * v.head(i+1) + + (alpha * conj(v.coeff(i))) * u.head(i+1); + } +}; + +template<bool Cond, typename T> struct conj_expr_if + : conditional<!Cond, const T&, + CwiseUnaryOp<scalar_conjugate_op<typename traits<T>::Scalar>,T> > {}; + +} // end namespace internal + +template<typename MatrixType, unsigned int UpLo> +template<typename DerivedU, typename DerivedV> +SelfAdjointView<MatrixType,UpLo>& SelfAdjointView<MatrixType,UpLo> +::rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, Scalar alpha) +{ + typedef internal::blas_traits<DerivedU> UBlasTraits; + typedef typename UBlasTraits::DirectLinearAccessType ActualUType; + typedef typename internal::remove_all<ActualUType>::type _ActualUType; + typename internal::add_const_on_value_type<ActualUType>::type actualU = UBlasTraits::extract(u.derived()); + + typedef internal::blas_traits<DerivedV> VBlasTraits; + typedef typename VBlasTraits::DirectLinearAccessType ActualVType; + typedef typename internal::remove_all<ActualVType>::type _ActualVType; + typename internal::add_const_on_value_type<ActualVType>::type actualV = VBlasTraits::extract(v.derived()); + + // If MatrixType is row major, then we use the routine for lower triangular in the upper triangular case and + // vice versa, and take the complex conjugate of all coefficients and vector entries. + + enum { IsRowMajor = (internal::traits<MatrixType>::Flags&RowMajorBit) ? 1 : 0 }; + Scalar actualAlpha = alpha * UBlasTraits::extractScalarFactor(u.derived()) + * internal::conj(VBlasTraits::extractScalarFactor(v.derived())); + if (IsRowMajor) + actualAlpha = internal::conj(actualAlpha); + + internal::selfadjoint_rank2_update_selector<Scalar, Index, + typename internal::remove_all<typename internal::conj_expr_if<IsRowMajor ^ UBlasTraits::NeedToConjugate,_ActualUType>::type>::type, + typename internal::remove_all<typename internal::conj_expr_if<IsRowMajor ^ VBlasTraits::NeedToConjugate,_ActualVType>::type>::type, + (IsRowMajor ? int(UpLo==Upper ? Lower : Upper) : UpLo)> + ::run(_expression().const_cast_derived().data(),_expression().outerStride(),actualU,actualV,actualAlpha); + + return *this; +} + +} // end namespace Eigen + +#endif // EIGEN_SELFADJOINTRANK2UPTADE_H |