Eigen  3.3.3
SelfAdjointView.h
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
00005 //
00006 // This Source Code Form is subject to the terms of the Mozilla
00007 // Public License v. 2.0. If a copy of the MPL was not distributed
00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
00009 
00010 #ifndef EIGEN_SELFADJOINTMATRIX_H
00011 #define EIGEN_SELFADJOINTMATRIX_H
00012 
00013 namespace Eigen { 
00014 
00031 namespace internal {
00032 template<typename MatrixType, unsigned int UpLo>
00033 struct traits<SelfAdjointView<MatrixType, UpLo> > : traits<MatrixType>
00034 {
00035   typedef typename ref_selector<MatrixType>::non_const_type MatrixTypeNested;
00036   typedef typename remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
00037   typedef MatrixType ExpressionType;
00038   typedef typename MatrixType::PlainObject FullMatrixType;
00039   enum {
00040     Mode = UpLo | SelfAdjoint,
00041     FlagsLvalueBit = is_lvalue<MatrixType>::value ? LvalueBit : 0,
00042     Flags =  MatrixTypeNestedCleaned::Flags & (HereditaryBits|FlagsLvalueBit)
00043            & (~(PacketAccessBit | DirectAccessBit | LinearAccessBit)) // FIXME these flags should be preserved
00044   };
00045 };
00046 }
00047 
00048 
00049 template<typename _MatrixType, unsigned int UpLo> class SelfAdjointView
00050   : public TriangularBase<SelfAdjointView<_MatrixType, UpLo> >
00051 {
00052   public:
00053 
00054     typedef _MatrixType MatrixType;
00055     typedef TriangularBase<SelfAdjointView> Base;
00056     typedef typename internal::traits<SelfAdjointView>::MatrixTypeNested MatrixTypeNested;
00057     typedef typename internal::traits<SelfAdjointView>::MatrixTypeNestedCleaned MatrixTypeNestedCleaned;
00058     typedef MatrixTypeNestedCleaned NestedExpression;
00059 
00061     typedef typename internal::traits<SelfAdjointView>::Scalar Scalar; 
00062     typedef typename MatrixType::StorageIndex StorageIndex;
00063     typedef typename internal::remove_all<typename MatrixType::ConjugateReturnType>::type MatrixConjugateReturnType;
00064 
00065     enum {
00066       Mode = internal::traits<SelfAdjointView>::Mode,
00067       Flags = internal::traits<SelfAdjointView>::Flags,
00068       TransposeMode = ((Mode & Upper) ? Lower : 0) | ((Mode & Lower) ? Upper : 0)
00069     };
00070     typedef typename MatrixType::PlainObject PlainObject;
00071 
00072     EIGEN_DEVICE_FUNC
00073     explicit inline SelfAdjointView(MatrixType& matrix) : m_matrix(matrix)
00074     {}
00075 
00076     EIGEN_DEVICE_FUNC
00077     inline Index rows() const { return m_matrix.rows(); }
00078     EIGEN_DEVICE_FUNC
00079     inline Index cols() const { return m_matrix.cols(); }
00080     EIGEN_DEVICE_FUNC
00081     inline Index outerStride() const { return m_matrix.outerStride(); }
00082     EIGEN_DEVICE_FUNC
00083     inline Index innerStride() const { return m_matrix.innerStride(); }
00084 
00088     EIGEN_DEVICE_FUNC
00089     inline Scalar coeff(Index row, Index col) const
00090     {
00091       Base::check_coordinates_internal(row, col);
00092       return m_matrix.coeff(row, col);
00093     }
00094 
00098     EIGEN_DEVICE_FUNC
00099     inline Scalar& coeffRef(Index row, Index col)
00100     {
00101       EIGEN_STATIC_ASSERT_LVALUE(SelfAdjointView);
00102       Base::check_coordinates_internal(row, col);
00103       return m_matrix.coeffRef(row, col);
00104     }
00105 
00107     EIGEN_DEVICE_FUNC
00108     const MatrixTypeNestedCleaned& _expression() const { return m_matrix; }
00109 
00110     EIGEN_DEVICE_FUNC
00111     const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
00112     EIGEN_DEVICE_FUNC
00113     MatrixTypeNestedCleaned& nestedExpression() { return m_matrix; }
00114 
00116     template<typename OtherDerived>
00117     EIGEN_DEVICE_FUNC
00118     const Product<SelfAdjointView,OtherDerived>
00119     operator*(const MatrixBase<OtherDerived>& rhs) const
00120     {
00121       return Product<SelfAdjointView,OtherDerived>(*this, rhs.derived());
00122     }
00123 
00125     template<typename OtherDerived> friend
00126     EIGEN_DEVICE_FUNC
00127     const Product<OtherDerived,SelfAdjointView>
00128     operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView& rhs)
00129     {
00130       return Product<OtherDerived,SelfAdjointView>(lhs.derived(),rhs);
00131     }
00132     
00133     friend EIGEN_DEVICE_FUNC
00134     const SelfAdjointView<const EIGEN_SCALAR_BINARYOP_EXPR_RETURN_TYPE(Scalar,MatrixType,product),UpLo>
00135     operator*(const Scalar& s, const SelfAdjointView& mat)
00136     {
00137       return (s*mat.nestedExpression()).template selfadjointView<UpLo>();
00138     }
00139 
00150     template<typename DerivedU, typename DerivedV>
00151     EIGEN_DEVICE_FUNC
00152     SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, const Scalar& alpha = Scalar(1));
00153 
00164     template<typename DerivedU>
00165     EIGEN_DEVICE_FUNC
00166     SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const Scalar& alpha = Scalar(1));
00167 
00178     template<unsigned int TriMode>
00179     EIGEN_DEVICE_FUNC
00180     typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
00181                                    TriangularView<MatrixType,TriMode>,
00182                                    TriangularView<typename MatrixType::AdjointReturnType,TriMode> >::type
00183     triangularView() const
00184     {
00185       typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), MatrixType&, typename MatrixType::ConstTransposeReturnType>::type tmp1(m_matrix);
00186       typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), MatrixType&, typename MatrixType::AdjointReturnType>::type tmp2(tmp1);
00187       return typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
00188                                    TriangularView<MatrixType,TriMode>,
00189                                    TriangularView<typename MatrixType::AdjointReturnType,TriMode> >::type(tmp2);
00190     }
00191 
00192     typedef SelfAdjointView<const MatrixConjugateReturnType,Mode> ConjugateReturnType;
00194     EIGEN_DEVICE_FUNC
00195     inline const ConjugateReturnType conjugate() const
00196     { return ConjugateReturnType(m_matrix.conjugate()); }
00197 
00198     typedef SelfAdjointView<const typename MatrixType::AdjointReturnType,TransposeMode> AdjointReturnType;
00200     EIGEN_DEVICE_FUNC
00201     inline const AdjointReturnType adjoint() const
00202     { return AdjointReturnType(m_matrix.adjoint()); }
00203 
00204     typedef SelfAdjointView<typename MatrixType::TransposeReturnType,TransposeMode> TransposeReturnType;
00206     EIGEN_DEVICE_FUNC
00207     inline TransposeReturnType transpose()
00208     {
00209       EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
00210       typename MatrixType::TransposeReturnType tmp(m_matrix);
00211       return TransposeReturnType(tmp);
00212     }
00213 
00214     typedef SelfAdjointView<const typename MatrixType::ConstTransposeReturnType,TransposeMode> ConstTransposeReturnType;
00216     EIGEN_DEVICE_FUNC
00217     inline const ConstTransposeReturnType transpose() const
00218     {
00219       return ConstTransposeReturnType(m_matrix.transpose());
00220     }
00221 
00227     EIGEN_DEVICE_FUNC
00228     typename MatrixType::ConstDiagonalReturnType diagonal() const
00229     {
00230       return typename MatrixType::ConstDiagonalReturnType(m_matrix);
00231     }
00232 
00234 
00235     const LLT<PlainObject, UpLo> llt() const;
00236     const LDLT<PlainObject, UpLo> ldlt() const;
00237 
00239 
00241     typedef typename NumTraits<Scalar>::Real RealScalar;
00243     typedef Matrix<RealScalar, internal::traits<MatrixType>::ColsAtCompileTime, 1> EigenvaluesReturnType;
00244 
00245     EIGEN_DEVICE_FUNC
00246     EigenvaluesReturnType eigenvalues() const;
00247     EIGEN_DEVICE_FUNC
00248     RealScalar operatorNorm() const;
00249 
00250   protected:
00251     MatrixTypeNested m_matrix;
00252 };
00253 
00254 
00255 // template<typename OtherDerived, typename MatrixType, unsigned int UpLo>
00256 // internal::selfadjoint_matrix_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >
00257 // operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView<MatrixType,UpLo>& rhs)
00258 // {
00259 //   return internal::matrix_selfadjoint_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >(lhs.derived(),rhs);
00260 // }
00261 
00262 // selfadjoint to dense matrix
00263 
00264 namespace internal {
00265 
00266 // TODO currently a selfadjoint expression has the form SelfAdjointView<.,.>
00267 //      in the future selfadjoint-ness should be defined by the expression traits
00268 //      such that Transpose<SelfAdjointView<.,.> > is valid. (currently TriangularBase::transpose() is overloaded to make it work)
00269 template<typename MatrixType, unsigned int Mode>
00270 struct evaluator_traits<SelfAdjointView<MatrixType,Mode> >
00271 {
00272   typedef typename storage_kind_to_evaluator_kind<typename MatrixType::StorageKind>::Kind Kind;
00273   typedef SelfAdjointShape Shape;
00274 };
00275 
00276 template<int UpLo, int SetOpposite, typename DstEvaluatorTypeT, typename SrcEvaluatorTypeT, typename Functor, int Version>
00277 class triangular_dense_assignment_kernel<UpLo,SelfAdjoint,SetOpposite,DstEvaluatorTypeT,SrcEvaluatorTypeT,Functor,Version>
00278   : public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version>
00279 {
00280 protected:
00281   typedef generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version> Base;
00282   typedef typename Base::DstXprType DstXprType;
00283   typedef typename Base::SrcXprType SrcXprType;
00284   using Base::m_dst;
00285   using Base::m_src;
00286   using Base::m_functor;
00287 public:
00288   
00289   typedef typename Base::DstEvaluatorType DstEvaluatorType;
00290   typedef typename Base::SrcEvaluatorType SrcEvaluatorType;
00291   typedef typename Base::Scalar Scalar;
00292   typedef typename Base::AssignmentTraits AssignmentTraits;
00293   
00294   
00295   EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst, const SrcEvaluatorType &src, const Functor &func, DstXprType& dstExpr)
00296     : Base(dst, src, func, dstExpr)
00297   {}
00298   
00299   EIGEN_DEVICE_FUNC void assignCoeff(Index row, Index col)
00300   {
00301     eigen_internal_assert(row!=col);
00302     Scalar tmp = m_src.coeff(row,col);
00303     m_functor.assignCoeff(m_dst.coeffRef(row,col), tmp);
00304     m_functor.assignCoeff(m_dst.coeffRef(col,row), numext::conj(tmp));
00305   }
00306   
00307   EIGEN_DEVICE_FUNC void assignDiagonalCoeff(Index id)
00308   {
00309     Base::assignCoeff(id,id);
00310   }
00311   
00312   EIGEN_DEVICE_FUNC void assignOppositeCoeff(Index, Index)
00313   { eigen_internal_assert(false && "should never be called"); }
00314 };
00315 
00316 } // end namespace internal
00317 
00318 /***************************************************************************
00319 * Implementation of MatrixBase methods
00320 ***************************************************************************/
00321 
00323 template<typename Derived>
00324 template<unsigned int UpLo>
00325 typename MatrixBase<Derived>::template ConstSelfAdjointViewReturnType<UpLo>::Type
00326 MatrixBase<Derived>::selfadjointView() const
00327 {
00328   return typename ConstSelfAdjointViewReturnType<UpLo>::Type(derived());
00329 }
00330 
00340 template<typename Derived>
00341 template<unsigned int UpLo>
00342 typename MatrixBase<Derived>::template SelfAdjointViewReturnType<UpLo>::Type
00343 MatrixBase<Derived>::selfadjointView()
00344 {
00345   return typename SelfAdjointViewReturnType<UpLo>::Type(derived());
00346 }
00347 
00348 } // end namespace Eigen
00349 
00350 #endif // EIGEN_SELFADJOINTMATRIX_H
 All Classes Functions Variables Typedefs Enumerations Enumerator Friends