Eigen  3.3.4
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 
76  EIGEN_DEVICE_FUNC
77  inline Index rows() const { return m_matrix.rows(); }
78  EIGEN_DEVICE_FUNC
79  inline Index cols() const { return m_matrix.cols(); }
80  EIGEN_DEVICE_FUNC
81  inline Index outerStride() const { return m_matrix.outerStride(); }
82  EIGEN_DEVICE_FUNC
83  inline Index innerStride() const { return m_matrix.innerStride(); }
84 
88  EIGEN_DEVICE_FUNC
89  inline Scalar coeff(Index row, Index col) const
90  {
91  Base::check_coordinates_internal(row, col);
92  return m_matrix.coeff(row, col);
93  }
94 
98  EIGEN_DEVICE_FUNC
99  inline Scalar& coeffRef(Index row, Index col)
100  {
101  EIGEN_STATIC_ASSERT_LVALUE(SelfAdjointView);
102  Base::check_coordinates_internal(row, col);
103  return m_matrix.coeffRef(row, col);
104  }
105 
107  EIGEN_DEVICE_FUNC
108  const MatrixTypeNestedCleaned& _expression() const { return m_matrix; }
109 
110  EIGEN_DEVICE_FUNC
111  const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
112  EIGEN_DEVICE_FUNC
113  MatrixTypeNestedCleaned& nestedExpression() { return m_matrix; }
114 
116  template<typename OtherDerived>
117  EIGEN_DEVICE_FUNC
120  {
121  return Product<SelfAdjointView,OtherDerived>(*this, rhs.derived());
122  }
123 
125  template<typename OtherDerived> friend
126  EIGEN_DEVICE_FUNC
128  operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView& rhs)
129  {
131  }
132 
133  friend EIGEN_DEVICE_FUNC
135  operator*(const Scalar& s, const SelfAdjointView& mat)
136  {
137  return (s*mat.nestedExpression()).template selfadjointView<UpLo>();
138  }
139 
150  template<typename DerivedU, typename DerivedV>
151  EIGEN_DEVICE_FUNC
152  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, const Scalar& alpha = Scalar(1));
153 
164  template<typename DerivedU>
165  EIGEN_DEVICE_FUNC
166  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const Scalar& alpha = Scalar(1));
167 
178  template<unsigned int TriMode>
179  EIGEN_DEVICE_FUNC
180  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
184  {
185  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), MatrixType&, typename MatrixType::ConstTransposeReturnType>::type tmp1(m_matrix);
186  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), MatrixType&, typename MatrixType::AdjointReturnType>::type tmp2(tmp1);
187  return typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
190  }
191 
194  EIGEN_DEVICE_FUNC
195  inline const ConjugateReturnType conjugate() const
196  { return ConjugateReturnType(m_matrix.conjugate()); }
197 
200  EIGEN_DEVICE_FUNC
201  inline const AdjointReturnType adjoint() const
202  { return AdjointReturnType(m_matrix.adjoint()); }
203 
206  EIGEN_DEVICE_FUNC
207  inline TransposeReturnType transpose()
208  {
209  EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
210  typename MatrixType::TransposeReturnType tmp(m_matrix);
211  return TransposeReturnType(tmp);
212  }
213 
216  EIGEN_DEVICE_FUNC
217  inline const ConstTransposeReturnType transpose() const
218  {
219  return ConstTransposeReturnType(m_matrix.transpose());
220  }
221 
227  EIGEN_DEVICE_FUNC
228  typename MatrixType::ConstDiagonalReturnType diagonal() const
229  {
230  return typename MatrixType::ConstDiagonalReturnType(m_matrix);
231  }
232 
234 
235  const LLT<PlainObject, UpLo> llt() const;
236  const LDLT<PlainObject, UpLo> ldlt() const;
237 
239 
244 
245  EIGEN_DEVICE_FUNC
246  EigenvaluesReturnType eigenvalues() const;
247  EIGEN_DEVICE_FUNC
248  RealScalar operatorNorm() const;
249 
250  protected:
251  MatrixTypeNested m_matrix;
252 };
253 
254 
255 // template<typename OtherDerived, typename MatrixType, unsigned int UpLo>
256 // internal::selfadjoint_matrix_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >
257 // operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView<MatrixType,UpLo>& rhs)
258 // {
259 // return internal::matrix_selfadjoint_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >(lhs.derived(),rhs);
260 // }
261 
262 // selfadjoint to dense matrix
263 
264 namespace internal {
265 
266 // TODO currently a selfadjoint expression has the form SelfAdjointView<.,.>
267 // in the future selfadjoint-ness should be defined by the expression traits
268 // such that Transpose<SelfAdjointView<.,.> > is valid. (currently TriangularBase::transpose() is overloaded to make it work)
269 template<typename MatrixType, unsigned int Mode>
270 struct evaluator_traits<SelfAdjointView<MatrixType,Mode> >
271 {
272  typedef typename storage_kind_to_evaluator_kind<typename MatrixType::StorageKind>::Kind Kind;
273  typedef SelfAdjointShape Shape;
274 };
275 
276 template<int UpLo, int SetOpposite, typename DstEvaluatorTypeT, typename SrcEvaluatorTypeT, typename Functor, int Version>
277 class triangular_dense_assignment_kernel<UpLo,SelfAdjoint,SetOpposite,DstEvaluatorTypeT,SrcEvaluatorTypeT,Functor,Version>
278  : public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version>
279 {
280 protected:
281  typedef generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version> Base;
282  typedef typename Base::DstXprType DstXprType;
283  typedef typename Base::SrcXprType SrcXprType;
284  using Base::m_dst;
285  using Base::m_src;
286  using Base::m_functor;
287 public:
288 
289  typedef typename Base::DstEvaluatorType DstEvaluatorType;
290  typedef typename Base::SrcEvaluatorType SrcEvaluatorType;
291  typedef typename Base::Scalar Scalar;
292  typedef typename Base::AssignmentTraits AssignmentTraits;
293 
294 
295  EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst, const SrcEvaluatorType &src, const Functor &func, DstXprType& dstExpr)
296  : Base(dst, src, func, dstExpr)
297  {}
298 
299  EIGEN_DEVICE_FUNC void assignCoeff(Index row, Index col)
300  {
301  eigen_internal_assert(row!=col);
302  Scalar tmp = m_src.coeff(row,col);
303  m_functor.assignCoeff(m_dst.coeffRef(row,col), tmp);
304  m_functor.assignCoeff(m_dst.coeffRef(col,row), numext::conj(tmp));
305  }
306 
307  EIGEN_DEVICE_FUNC void assignDiagonalCoeff(Index id)
308  {
309  Base::assignCoeff(id,id);
310  }
311 
312  EIGEN_DEVICE_FUNC void assignOppositeCoeff(Index, Index)
313  { eigen_internal_assert(false && "should never be called"); }
314 };
315 
316 } // end namespace internal
317 
318 /***************************************************************************
319 * Implementation of MatrixBase methods
320 ***************************************************************************/
321 
323 template<typename Derived>
324 template<unsigned int UpLo>
325 typename MatrixBase<Derived>::template ConstSelfAdjointViewReturnType<UpLo>::Type
327 {
328  return typename ConstSelfAdjointViewReturnType<UpLo>::Type(derived());
329 }
330 
340 template<typename Derived>
341 template<unsigned int UpLo>
342 typename MatrixBase<Derived>::template SelfAdjointViewReturnType<UpLo>::Type
344 {
345  return typename SelfAdjointViewReturnType<UpLo>::Type(derived());
346 }
347 
348 } // end namespace Eigen
349 
350 #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:243
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:207
Namespace containing all symbols from the Eigen library.
Definition: Core:287
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:119
NumTraits< Scalar >::Real RealScalar
Definition: SelfAdjointView.h:241
Scalar coeff(Index row, Index col) const
Definition: SelfAdjointView.h:89
Derived & derived()
Definition: EigenBase.h:45
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
const ConstTransposeReturnType transpose() const
Definition: SelfAdjointView.h:217
Definition: Constants.h:204
const ConjugateReturnType conjugate() const
Definition: SelfAdjointView.h:195
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:52
const Product< MatrixDerived, PermutationDerived, AliasFreeProduct > operator*(const MatrixBase< MatrixDerived > &matrix, const PermutationBase< PermutationDerived > &permutation)
Definition: PermutationMatrix.h:543
Expression of a selfadjoint matrix from a triangular part of a dense matrix.
Definition: SelfAdjointView.h:49
MatrixType::ConstDiagonalReturnType diagonal() const
Definition: SelfAdjointView.h:228
friend const Product< OtherDerived, SelfAdjointView > operator*(const MatrixBase< OtherDerived > &lhs, const SelfAdjointView &rhs)
Definition: SelfAdjointView.h:128
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:99
internal::conditional<(TriMode &(Upper|Lower))==(UpLo &(Upper|Lower)), TriangularView< MatrixType, TriMode >, TriangularView< typename MatrixType::AdjointReturnType, TriMode > >::type triangularView() const
Definition: SelfAdjointView.h:183
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:178
const AdjointReturnType adjoint() const
Definition: SelfAdjointView.h:201
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const unsigned int LinearAccessBit
Definition: Constants.h:125