Eigen  3.3.7
SelfAdjointView.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN_SELFADJOINTMATRIX_H
11 #define EIGEN_SELFADJOINTMATRIX_H
12 
13 namespace Eigen {
14 
31 namespace internal {
32 template<typename MatrixType, unsigned int UpLo>
33 struct traits<SelfAdjointView<MatrixType, UpLo> > : traits<MatrixType>
34 {
35  typedef typename ref_selector<MatrixType>::non_const_type MatrixTypeNested;
36  typedef typename remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
37  typedef MatrixType ExpressionType;
38  typedef typename MatrixType::PlainObject FullMatrixType;
39  enum {
40  Mode = UpLo | SelfAdjoint,
41  FlagsLvalueBit = is_lvalue<MatrixType>::value ? LvalueBit : 0,
42  Flags = MatrixTypeNestedCleaned::Flags & (HereditaryBits|FlagsLvalueBit)
43  & (~(PacketAccessBit | DirectAccessBit | LinearAccessBit)) // FIXME these flags should be preserved
44  };
45 };
46 }
47 
48 
49 template<typename _MatrixType, unsigned int UpLo> class SelfAdjointView
50  : public TriangularBase<SelfAdjointView<_MatrixType, UpLo> >
51 {
52  public:
53 
54  typedef _MatrixType MatrixType;
56  typedef typename internal::traits<SelfAdjointView>::MatrixTypeNested MatrixTypeNested;
57  typedef typename internal::traits<SelfAdjointView>::MatrixTypeNestedCleaned MatrixTypeNestedCleaned;
58  typedef MatrixTypeNestedCleaned NestedExpression;
59 
61  typedef typename internal::traits<SelfAdjointView>::Scalar Scalar;
62  typedef typename MatrixType::StorageIndex StorageIndex;
63  typedef typename internal::remove_all<typename MatrixType::ConjugateReturnType>::type MatrixConjugateReturnType;
64 
65  enum {
66  Mode = internal::traits<SelfAdjointView>::Mode,
67  Flags = internal::traits<SelfAdjointView>::Flags,
68  TransposeMode = ((Mode & Upper) ? Lower : 0) | ((Mode & Lower) ? Upper : 0)
69  };
70  typedef typename MatrixType::PlainObject PlainObject;
71 
72  EIGEN_DEVICE_FUNC
73  explicit inline SelfAdjointView(MatrixType& matrix) : m_matrix(matrix)
74  {
75  EIGEN_STATIC_ASSERT(UpLo==Lower || UpLo==Upper,SELFADJOINTVIEW_ACCEPTS_UPPER_AND_LOWER_MODE_ONLY);
76  }
77 
78  EIGEN_DEVICE_FUNC
79  inline Index rows() const { return m_matrix.rows(); }
80  EIGEN_DEVICE_FUNC
81  inline Index cols() const { return m_matrix.cols(); }
82  EIGEN_DEVICE_FUNC
83  inline Index outerStride() const { return m_matrix.outerStride(); }
84  EIGEN_DEVICE_FUNC
85  inline Index innerStride() const { return m_matrix.innerStride(); }
86 
90  EIGEN_DEVICE_FUNC
91  inline Scalar coeff(Index row, Index col) const
92  {
93  Base::check_coordinates_internal(row, col);
94  return m_matrix.coeff(row, col);
95  }
96 
100  EIGEN_DEVICE_FUNC
101  inline Scalar& coeffRef(Index row, Index col)
102  {
103  EIGEN_STATIC_ASSERT_LVALUE(SelfAdjointView);
104  Base::check_coordinates_internal(row, col);
105  return m_matrix.coeffRef(row, col);
106  }
107 
109  EIGEN_DEVICE_FUNC
110  const MatrixTypeNestedCleaned& _expression() const { return m_matrix; }
111 
112  EIGEN_DEVICE_FUNC
113  const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
114  EIGEN_DEVICE_FUNC
115  MatrixTypeNestedCleaned& nestedExpression() { return m_matrix; }
116 
118  template<typename OtherDerived>
119  EIGEN_DEVICE_FUNC
120  const Product<SelfAdjointView,OtherDerived>
122  {
123  return Product<SelfAdjointView,OtherDerived>(*this, rhs.derived());
124  }
125 
127  template<typename OtherDerived> friend
128  EIGEN_DEVICE_FUNC
131  {
133  }
134 
135  friend EIGEN_DEVICE_FUNC
137  operator*(const Scalar& s, const SelfAdjointView& mat)
138  {
139  return (s*mat.nestedExpression()).template selfadjointView<UpLo>();
140  }
141 
152  template<typename DerivedU, typename DerivedV>
153  EIGEN_DEVICE_FUNC
154  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, const Scalar& alpha = Scalar(1));
155 
166  template<typename DerivedU>
167  EIGEN_DEVICE_FUNC
168  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const Scalar& alpha = Scalar(1));
169 
180  template<unsigned int TriMode>
181  EIGEN_DEVICE_FUNC
182  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
183  TriangularView<MatrixType,TriMode>,
184  TriangularView<typename MatrixType::AdjointReturnType,TriMode> >::type
186  {
187  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), MatrixType&, typename MatrixType::ConstTransposeReturnType>::type tmp1(m_matrix);
188  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), MatrixType&, typename MatrixType::AdjointReturnType>::type tmp2(tmp1);
189  return typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
192  }
193 
194  typedef SelfAdjointView<const MatrixConjugateReturnType,UpLo> ConjugateReturnType;
196  EIGEN_DEVICE_FUNC
197  inline const ConjugateReturnType conjugate() const
198  { return ConjugateReturnType(m_matrix.conjugate()); }
199 
202  EIGEN_DEVICE_FUNC
203  inline const AdjointReturnType adjoint() const
204  { return AdjointReturnType(m_matrix.adjoint()); }
205 
208  EIGEN_DEVICE_FUNC
210  {
211  EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
212  typename MatrixType::TransposeReturnType tmp(m_matrix);
213  return TransposeReturnType(tmp);
214  }
215 
218  EIGEN_DEVICE_FUNC
219  inline const ConstTransposeReturnType transpose() const
220  {
221  return ConstTransposeReturnType(m_matrix.transpose());
222  }
223 
229  EIGEN_DEVICE_FUNC
230  typename MatrixType::ConstDiagonalReturnType diagonal() const
231  {
232  return typename MatrixType::ConstDiagonalReturnType(m_matrix);
233  }
234 
236 
237  const LLT<PlainObject, UpLo> llt() const;
238  const LDLT<PlainObject, UpLo> ldlt() const;
239 
241 
246 
247  EIGEN_DEVICE_FUNC
249  EIGEN_DEVICE_FUNC
250  RealScalar operatorNorm() const;
251 
252  protected:
253  MatrixTypeNested m_matrix;
254 };
255 
256 
257 // template<typename OtherDerived, typename MatrixType, unsigned int UpLo>
258 // internal::selfadjoint_matrix_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >
259 // operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView<MatrixType,UpLo>& rhs)
260 // {
261 // return internal::matrix_selfadjoint_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >(lhs.derived(),rhs);
262 // }
263 
264 // selfadjoint to dense matrix
265 
266 namespace internal {
267 
268 // TODO currently a selfadjoint expression has the form SelfAdjointView<.,.>
269 // in the future selfadjoint-ness should be defined by the expression traits
270 // such that Transpose<SelfAdjointView<.,.> > is valid. (currently TriangularBase::transpose() is overloaded to make it work)
271 template<typename MatrixType, unsigned int Mode>
272 struct evaluator_traits<SelfAdjointView<MatrixType,Mode> >
273 {
274  typedef typename storage_kind_to_evaluator_kind<typename MatrixType::StorageKind>::Kind Kind;
275  typedef SelfAdjointShape Shape;
276 };
277 
278 template<int UpLo, int SetOpposite, typename DstEvaluatorTypeT, typename SrcEvaluatorTypeT, typename Functor, int Version>
279 class triangular_dense_assignment_kernel<UpLo,SelfAdjoint,SetOpposite,DstEvaluatorTypeT,SrcEvaluatorTypeT,Functor,Version>
280  : public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version>
281 {
282 protected:
283  typedef generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version> Base;
284  typedef typename Base::DstXprType DstXprType;
285  typedef typename Base::SrcXprType SrcXprType;
286  using Base::m_dst;
287  using Base::m_src;
288  using Base::m_functor;
289 public:
290 
291  typedef typename Base::DstEvaluatorType DstEvaluatorType;
292  typedef typename Base::SrcEvaluatorType SrcEvaluatorType;
293  typedef typename Base::Scalar Scalar;
294  typedef typename Base::AssignmentTraits AssignmentTraits;
295 
296 
297  EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst, const SrcEvaluatorType &src, const Functor &func, DstXprType& dstExpr)
298  : Base(dst, src, func, dstExpr)
299  {}
300 
301  EIGEN_DEVICE_FUNC void assignCoeff(Index row, Index col)
302  {
303  eigen_internal_assert(row!=col);
304  Scalar tmp = m_src.coeff(row,col);
305  m_functor.assignCoeff(m_dst.coeffRef(row,col), tmp);
306  m_functor.assignCoeff(m_dst.coeffRef(col,row), numext::conj(tmp));
307  }
308 
309  EIGEN_DEVICE_FUNC void assignDiagonalCoeff(Index id)
310  {
311  Base::assignCoeff(id,id);
312  }
313 
314  EIGEN_DEVICE_FUNC void assignOppositeCoeff(Index, Index)
315  { eigen_internal_assert(false && "should never be called"); }
316 };
317 
318 } // end namespace internal
319 
320 /***************************************************************************
321 * Implementation of MatrixBase methods
322 ***************************************************************************/
323 
325 template<typename Derived>
326 template<unsigned int UpLo>
327 typename MatrixBase<Derived>::template ConstSelfAdjointViewReturnType<UpLo>::Type
329 {
330  return typename ConstSelfAdjointViewReturnType<UpLo>::Type(derived());
331 }
332 
342 template<typename Derived>
343 template<unsigned int UpLo>
344 typename MatrixBase<Derived>::template SelfAdjointViewReturnType<UpLo>::Type
346 {
347  return typename SelfAdjointViewReturnType<UpLo>::Type(derived());
348 }
349 
350 } // end namespace Eigen
351 
352 #endif // EIGEN_SELFADJOINTMATRIX_H
Robust Cholesky decomposition of a matrix with pivoting.
Definition: LDLT.h:50
Matrix< RealScalar, internal::traits< MatrixType >::ColsAtCompileTime, 1 > EigenvaluesReturnType
Definition: SelfAdjointView.h:245
const LDLT< PlainObject, UpLo > ldlt() const
Definition: LDLT.h:655
SelfAdjointView & rankUpdate(const MatrixBase< DerivedU > &u, const MatrixBase< DerivedV > &v, const Scalar &alpha=Scalar(1))
Expression of the product of two arbitrary matrices or vectors.
Definition: Product.h:71
Base class for triangular part in a matrix.
Definition: TriangularMatrix.h:27
const unsigned int DirectAccessBit
Definition: Constants.h:150
const unsigned int LvalueBit
Definition: Constants.h:139
TransposeReturnType transpose()
Definition: SelfAdjointView.h:209
Namespace containing all symbols from the Eigen library.
Definition: Core:306
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:150
const Product< SelfAdjointView, OtherDerived > operator*(const MatrixBase< OtherDerived > &rhs) const
Definition: SelfAdjointView.h:121
NumTraits< Scalar >::Real RealScalar
Definition: SelfAdjointView.h:243
Scalar coeff(Index row, Index col) const
Definition: SelfAdjointView.h:91
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:38
const unsigned int PacketAccessBit
Definition: Constants.h:89
internal::traits< SelfAdjointView >::Scalar Scalar
The type of coefficients in this matrix.
Definition: SelfAdjointView.h:61
RealScalar operatorNorm() const
Computes the L2 operator norm.
Definition: MatrixBaseEigenvalues.h:151
const ConstTransposeReturnType transpose() const
Definition: SelfAdjointView.h:219
Definition: Constants.h:204
Derived & derived()
Definition: EigenBase.h:45
const ConjugateReturnType conjugate() const
Definition: SelfAdjointView.h:197
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:56
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:33
Expression of a selfadjoint matrix from a triangular part of a dense matrix.
Definition: SelfAdjointView.h:49
MatrixType::ConstDiagonalReturnType diagonal() const
Definition: SelfAdjointView.h:230
friend const Product< OtherDerived, SelfAdjointView > operator*(const MatrixBase< OtherDerived > &lhs, const SelfAdjointView &rhs)
Definition: SelfAdjointView.h:130
Definition: Constants.h:206
Definition: Eigen_Colamd.h:50
Definition: Constants.h:220
Expression of a triangular part in a matrix.
Definition: TriangularMatrix.h:186
Scalar & coeffRef(Index row, Index col)
Definition: SelfAdjointView.h:101
EigenvaluesReturnType eigenvalues() const
Computes the eigenvalues of a matrix.
Definition: MatrixBaseEigenvalues.h:88
internal::conditional<(TriMode &(Upper|Lower))==(UpLo &(Upper|Lower)), TriangularView< MatrixType, TriMode >, TriangularView< typename MatrixType::AdjointReturnType, TriMode > >::type triangularView() const
Definition: SelfAdjointView.h:185
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:178
const AdjointReturnType adjoint() const
Definition: SelfAdjointView.h:203
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const unsigned int LinearAccessBit
Definition: Constants.h:125
const LLT< PlainObject, UpLo > llt() const
Definition: LLT.h:535