Eigen  3.2.92
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>::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 // FIXME could also be called SelfAdjointWrapper to be consistent with DiagonalWrapper ??
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 
60  typedef typename internal::traits<SelfAdjointView>::Scalar Scalar;
61  typedef typename MatrixType::StorageIndex StorageIndex;
62 
63  enum {
64  Mode = internal::traits<SelfAdjointView>::Mode,
65  Flags = internal::traits<SelfAdjointView>::Flags
66  };
67  typedef typename MatrixType::PlainObject PlainObject;
68 
69  EIGEN_DEVICE_FUNC
70  explicit inline SelfAdjointView(MatrixType& matrix) : m_matrix(matrix)
71  {}
72 
73  EIGEN_DEVICE_FUNC
74  inline Index rows() const { return m_matrix.rows(); }
75  EIGEN_DEVICE_FUNC
76  inline Index cols() const { return m_matrix.cols(); }
77  EIGEN_DEVICE_FUNC
78  inline Index outerStride() const { return m_matrix.outerStride(); }
79  EIGEN_DEVICE_FUNC
80  inline Index innerStride() const { return m_matrix.innerStride(); }
81 
85  EIGEN_DEVICE_FUNC
86  inline Scalar coeff(Index row, Index col) const
87  {
88  Base::check_coordinates_internal(row, col);
89  return m_matrix.coeff(row, col);
90  }
91 
95  EIGEN_DEVICE_FUNC
96  inline Scalar& coeffRef(Index row, Index col)
97  {
98  EIGEN_STATIC_ASSERT_LVALUE(SelfAdjointView);
99  Base::check_coordinates_internal(row, col);
100  return m_matrix.const_cast_derived().coeffRef(row, col);
101  }
102 
104  EIGEN_DEVICE_FUNC
105  const MatrixTypeNestedCleaned& _expression() const { return m_matrix; }
106 
107  EIGEN_DEVICE_FUNC
108  const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
109  EIGEN_DEVICE_FUNC
110  MatrixTypeNestedCleaned& nestedExpression() { return *const_cast<MatrixTypeNestedCleaned*>(&m_matrix); }
111 
113  template<typename OtherDerived>
114  EIGEN_DEVICE_FUNC
115  const Product<SelfAdjointView,OtherDerived>
117  {
118  return Product<SelfAdjointView,OtherDerived>(*this, rhs.derived());
119  }
120 
122  template<typename OtherDerived> friend
123  EIGEN_DEVICE_FUNC
126  {
127  return Product<OtherDerived,SelfAdjointView>(lhs.derived(),rhs);
128  }
129 
130  friend EIGEN_DEVICE_FUNC
132  operator*(const Scalar& s, const SelfAdjointView& mat)
133  {
134  return (s*mat.nestedExpression()).template selfadjointView<UpLo>();
135  }
136 
147  template<typename DerivedU, typename DerivedV>
148  EIGEN_DEVICE_FUNC
149  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, const Scalar& alpha = Scalar(1));
150 
161  template<typename DerivedU>
162  EIGEN_DEVICE_FUNC
163  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const Scalar& alpha = Scalar(1));
164 
166 
167  const LLT<PlainObject, UpLo> llt() const;
168  const LDLT<PlainObject, UpLo> ldlt() const;
169 
171 
176 
177  EIGEN_DEVICE_FUNC
178  EigenvaluesReturnType eigenvalues() const;
179  EIGEN_DEVICE_FUNC
180  RealScalar operatorNorm() const;
181 
182  protected:
183  MatrixTypeNested m_matrix;
184 };
185 
186 
187 // template<typename OtherDerived, typename MatrixType, unsigned int UpLo>
188 // internal::selfadjoint_matrix_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >
189 // operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView<MatrixType,UpLo>& rhs)
190 // {
191 // return internal::matrix_selfadjoint_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >(lhs.derived(),rhs);
192 // }
193 
194 // selfadjoint to dense matrix
195 
196 namespace internal {
197 
198 // TODO currently a selfadjoint expression has the form SelfAdjointView<.,.>
199 // in the future selfadjoint-ness should be defined by the expression traits
200 // such that Transpose<SelfAdjointView<.,.> > is valid. (currently TriangularBase::transpose() is overloaded to make it work)
201 template<typename MatrixType, unsigned int Mode>
202 struct evaluator_traits<SelfAdjointView<MatrixType,Mode> >
203 {
204  typedef typename storage_kind_to_evaluator_kind<typename MatrixType::StorageKind>::Kind Kind;
205  typedef SelfAdjointShape Shape;
206 
207  static const int AssumeAliasing = 0;
208 };
209 
210 template<int UpLo, int SetOpposite, typename DstEvaluatorTypeT, typename SrcEvaluatorTypeT, typename Functor, int Version>
211 class triangular_dense_assignment_kernel<UpLo,SelfAdjoint,SetOpposite,DstEvaluatorTypeT,SrcEvaluatorTypeT,Functor,Version>
212  : public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version>
213 {
214 protected:
215  typedef generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version> Base;
216  typedef typename Base::DstXprType DstXprType;
217  typedef typename Base::SrcXprType SrcXprType;
218  using Base::m_dst;
219  using Base::m_src;
220  using Base::m_functor;
221 public:
222 
223  typedef typename Base::DstEvaluatorType DstEvaluatorType;
224  typedef typename Base::SrcEvaluatorType SrcEvaluatorType;
225  typedef typename Base::Scalar Scalar;
226  typedef typename Base::AssignmentTraits AssignmentTraits;
227 
228 
229  EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst, const SrcEvaluatorType &src, const Functor &func, DstXprType& dstExpr)
230  : Base(dst, src, func, dstExpr)
231  {}
232 
233  EIGEN_DEVICE_FUNC void assignCoeff(Index row, Index col)
234  {
235  eigen_internal_assert(row!=col);
236  Scalar tmp = m_src.coeff(row,col);
237  m_functor.assignCoeff(m_dst.coeffRef(row,col), tmp);
238  m_functor.assignCoeff(m_dst.coeffRef(col,row), numext::conj(tmp));
239  }
240 
241  EIGEN_DEVICE_FUNC void assignDiagonalCoeff(Index id)
242  {
243  Base::assignCoeff(id,id);
244  }
245 
246  EIGEN_DEVICE_FUNC void assignOppositeCoeff(Index, Index)
247  { eigen_internal_assert(false && "should never be called"); }
248 };
249 
250 } // end namespace internal
251 
252 /***************************************************************************
253 * Implementation of MatrixBase methods
254 ***************************************************************************/
255 
256 template<typename Derived>
257 template<unsigned int UpLo>
258 typename MatrixBase<Derived>::template ConstSelfAdjointViewReturnType<UpLo>::Type
259 MatrixBase<Derived>::selfadjointView() const
260 {
261  return typename ConstSelfAdjointViewReturnType<UpLo>::Type(derived());
262 }
263 
264 template<typename Derived>
265 template<unsigned int UpLo>
266 typename MatrixBase<Derived>::template SelfAdjointViewReturnType<UpLo>::Type
267 MatrixBase<Derived>::selfadjointView()
268 {
269  return typename SelfAdjointViewReturnType<UpLo>::Type(derived());
270 }
271 
272 } // end namespace Eigen
273 
274 #endif // EIGEN_SELFADJOINTMATRIX_H
Matrix< RealScalar, internal::traits< MatrixType >::ColsAtCompileTime, 1 > EigenvaluesReturnType
Definition: SelfAdjointView.h:175
const Product< SelfAdjointView, OtherDerived > operator*(const MatrixBase< OtherDerived > &rhs) const
Definition: SelfAdjointView.h:116
Expression of the product of two arbitrary matrices or vectors.
Definition: Product.h:107
Base class for triangular part in a matrix.
Definition: TriangularMatrix.h:27
const unsigned int DirectAccessBit
Definition: Constants.h:149
const unsigned int LvalueBit
Definition: Constants.h:138
Definition: LDLT.h:16
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:107
NumTraits< Scalar >::Real RealScalar
Definition: SelfAdjointView.h:173
Scalar coeff(Index row, Index col) const
Definition: SelfAdjointView.h:86
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:37
const unsigned int PacketAccessBit
Definition: Constants.h:88
internal::traits< SelfAdjointView >::Scalar Scalar
The type of coefficients in this matrix.
Definition: SelfAdjointView.h:60
Expression of a selfadjoint matrix from a triangular part of a dense matrix.
Definition: SelfAdjointView.h:49
friend const Product< OtherDerived, SelfAdjointView > operator*(const MatrixBase< OtherDerived > &lhs, const SelfAdjointView &rhs)
Definition: SelfAdjointView.h:125
Definition: Constants.h:220
Definition: Eigen_Colamd.h:54
Scalar & coeffRef(Index row, Index col)
Definition: SelfAdjointView.h:96
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:178
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const unsigned int LinearAccessBit
Definition: Constants.h:124