SelfAdjointView.h
Go to the documentation of this file.
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 {
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;
58  typedef MatrixTypeNestedCleaned NestedExpression;
59 
62  typedef typename MatrixType::StorageIndex StorageIndex;
65 
66  enum {
69  TransposeMode = ((int(Mode) & int(Upper)) ? Lower : 0) | ((int(Mode) & int(Lower)) ? Upper : 0)
70  };
71  typedef typename MatrixType::PlainObject PlainObject;
72 
74  explicit inline SelfAdjointView(MatrixType& matrix) : m_matrix(matrix)
75  {
76  EIGEN_STATIC_ASSERT(UpLo==Lower || UpLo==Upper,SELFADJOINTVIEW_ACCEPTS_UPPER_AND_LOWER_MODE_ONLY);
77  }
78 
80  inline Index rows() const EIGEN_NOEXCEPT { return m_matrix.rows(); }
82  inline Index cols() const EIGEN_NOEXCEPT { return m_matrix.cols(); }
84  inline Index outerStride() const EIGEN_NOEXCEPT { return m_matrix.outerStride(); }
86  inline Index innerStride() const EIGEN_NOEXCEPT { return m_matrix.innerStride(); }
87 
92  inline Scalar coeff(Index row, Index col) const
93  {
94  Base::check_coordinates_internal(row, col);
95  return m_matrix.coeff(row, col);
96  }
97 
102  inline Scalar& coeffRef(Index row, Index col)
103  {
105  Base::check_coordinates_internal(row, col);
106  return m_matrix.coeffRef(row, col);
107  }
108 
111  const MatrixTypeNestedCleaned& _expression() const { return m_matrix; }
112 
114  const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
116  MatrixTypeNestedCleaned& nestedExpression() { return m_matrix; }
117 
119  template<typename OtherDerived>
123  {
124  return Product<SelfAdjointView,OtherDerived>(*this, rhs.derived());
125  }
126 
128  template<typename OtherDerived> friend
132  {
133  return Product<OtherDerived,SelfAdjointView>(lhs.derived(),rhs);
134  }
135 
136  friend EIGEN_DEVICE_FUNC
138  operator*(const Scalar& s, const SelfAdjointView& mat)
139  {
140  return (s*mat.nestedExpression()).template selfadjointView<UpLo>();
141  }
142 
153  template<typename DerivedU, typename DerivedV>
155  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, const Scalar& alpha = Scalar(1));
156 
167  template<typename DerivedU>
169  SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const Scalar& alpha = Scalar(1));
170 
181  template<unsigned int TriMode>
183  typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
187  {
190  return typename internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)),
193  }
194 
198  inline const ConjugateReturnType conjugate() const
199  { return ConjugateReturnType(m_matrix.conjugate()); }
200 
204  template<bool Cond>
207  conjugateIf() const
208  {
210  return ReturnType(m_matrix.template conjugateIf<Cond>());
211  }
212 
216  inline const AdjointReturnType adjoint() const
217  { return AdjointReturnType(m_matrix.adjoint()); }
218 
222  inline TransposeReturnType transpose()
223  {
224  EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
225  typename MatrixType::TransposeReturnType tmp(m_matrix);
226  return TransposeReturnType(tmp);
227  }
228 
232  inline const ConstTransposeReturnType transpose() const
233  {
234  return ConstTransposeReturnType(m_matrix.transpose());
235  }
236 
243  typename MatrixType::ConstDiagonalReturnType diagonal() const
244  {
245  return typename MatrixType::ConstDiagonalReturnType(m_matrix);
246  }
247 
249 
250  const LLT<PlainObject, UpLo> llt() const;
251  const LDLT<PlainObject, UpLo> ldlt() const;
252 
254 
259 
261  EigenvaluesReturnType eigenvalues() const;
263  RealScalar operatorNorm() const;
264 
265  protected:
266  MatrixTypeNested m_matrix;
267 };
268 
269 
270 // template<typename OtherDerived, typename MatrixType, unsigned int UpLo>
271 // internal::selfadjoint_matrix_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >
272 // operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView<MatrixType,UpLo>& rhs)
273 // {
274 // return internal::matrix_selfadjoint_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >(lhs.derived(),rhs);
275 // }
276 
277 // selfadjoint to dense matrix
278 
279 namespace internal {
280 
281 // TODO currently a selfadjoint expression has the form SelfAdjointView<.,.>
282 // in the future selfadjoint-ness should be defined by the expression traits
283 // such that Transpose<SelfAdjointView<.,.> > is valid. (currently TriangularBase::transpose() is overloaded to make it work)
284 template<typename MatrixType, unsigned int Mode>
286 {
289 };
290 
291 template<int UpLo, int SetOpposite, typename DstEvaluatorTypeT, typename SrcEvaluatorTypeT, typename Functor, int Version>
292 class triangular_dense_assignment_kernel<UpLo,SelfAdjoint,SetOpposite,DstEvaluatorTypeT,SrcEvaluatorTypeT,Functor,Version>
293  : public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version>
294 {
295 protected:
297  typedef typename Base::DstXprType DstXprType;
298  typedef typename Base::SrcXprType SrcXprType;
299  using Base::m_dst;
300  using Base::m_src;
301  using Base::m_functor;
302 public:
303 
306  typedef typename Base::Scalar Scalar;
308 
309 
310  EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst, const SrcEvaluatorType &src, const Functor &func, DstXprType& dstExpr)
311  : Base(dst, src, func, dstExpr)
312  {}
313 
315  {
316  eigen_internal_assert(row!=col);
317  Scalar tmp = m_src.coeff(row,col);
318  m_functor.assignCoeff(m_dst.coeffRef(row,col), tmp);
319  m_functor.assignCoeff(m_dst.coeffRef(col,row), numext::conj(tmp));
320  }
321 
323  {
324  Base::assignCoeff(id,id);
325  }
326 
328  { eigen_internal_assert(false && "should never be called"); }
329 };
330 
331 } // end namespace internal
332 
333 /***************************************************************************
334 * Implementation of MatrixBase methods
335 ***************************************************************************/
336 
338 template<typename Derived>
339 template<unsigned int UpLo>
342 {
343  return typename ConstSelfAdjointViewReturnType<UpLo>::Type(derived());
344 }
345 
355 template<typename Derived>
356 template<unsigned int UpLo>
359 {
360  return typename SelfAdjointViewReturnType<UpLo>::Type(derived());
361 }
362 
363 } // end namespace Eigen
364 
365 #endif // EIGEN_SELFADJOINTMATRIX_H
Robust Cholesky decomposition of a matrix with pivoting.
Definition: LDLT.h:59
Matrix< RealScalar, internal::traits< MatrixType >::ColsAtCompileTime, 1 > EigenvaluesReturnType
friend EIGEN_DEVICE_FUNC const SelfAdjointView< const EIGEN_SCALAR_BINARYOP_EXPR_RETURN_TYPE(Scalar, MatrixType, product), UpLo > operator*(const Scalar &s, const SelfAdjointView &mat)
SCALAR Scalar
Definition: bench_gemm.cpp:46
Expression of the product of two arbitrary matrices or vectors.
Definition: Product.h:71
EIGEN_DEVICE_FUNC const ConjugateReturnType conjugate() const
MatrixType::StorageIndex StorageIndex
Base class for triangular part in a matrix.
const unsigned int DirectAccessBit
Definition: Constants.h:155
EIGEN_DEVICE_FUNC MatrixType::ConstDiagonalReturnType diagonal() const
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index innerStride() const EIGEN_NOEXCEPT
const unsigned int LvalueBit
Definition: Constants.h:144
SelfAdjointView< const typename MatrixType::AdjointReturnType, TransposeMode > AdjointReturnType
Namespace containing all symbols from the Eigen library.
Definition: jet.h:637
internal::remove_all< typename MatrixType::ConjugateReturnType >::type MatrixConjugateReturnType
EIGEN_DEVICE_FUNC internal::conditional< Cond, ConjugateReturnType, ConstSelfAdjointView >::type conjugateIf() const
MatrixXf MatrixType
MatrixType::PlainObject PlainObject
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:232
NumTraits< Scalar >::Real RealScalar
EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst, const SrcEvaluatorType &src, const Functor &func, DstXprType &dstExpr)
#define EIGEN_STATIC_ASSERT(CONDITION, MSG)
Definition: StaticAssert.h:127
EIGEN_DEVICE_FUNC TransposeReturnType transpose()
storage_kind_to_evaluator_kind< typename MatrixType::StorageKind >::Kind Kind
SelfAdjointView< typename internal::add_const< MatrixType >::type, UpLo > ConstSelfAdjointView
EIGEN_DEVICE_FUNC const AdjointReturnType adjoint() const
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:39
EIGEN_DEVICE_FUNC SelfAdjointViewReturnType< UpLo >::Type selfadjointView()
const unsigned int PacketAccessBit
Definition: Constants.h:94
EIGEN_DEVICE_FUNC const MatrixTypeNestedCleaned & _expression() const
#define EIGEN_STATIC_ASSERT_LVALUE(Derived)
Definition: StaticAssert.h:202
AnnoyingScalar conj(const AnnoyingScalar &x)
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index outerStride() const EIGEN_NOEXCEPT
internal::traits< SelfAdjointView >::Scalar Scalar
The type of coefficients in this matrix.
internal::traits< SelfAdjointView >::MatrixTypeNestedCleaned MatrixTypeNestedCleaned
EIGEN_DEVICE_FUNC MatrixTypeNestedCleaned & nestedExpression()
internal::traits< SelfAdjointView >::MatrixTypeNested MatrixTypeNested
EIGEN_DEVICE_FUNC const MatrixTypeNestedCleaned & nestedExpression() const
TriangularBase< SelfAdjointView > Base
EIGEN_DONT_INLINE void llt(const Mat &A, const Mat &B, Mat &C)
Definition: llt.cpp:5
const unsigned int HereditaryBits
Definition: Constants.h:195
m row(1)
#define EIGEN_NOEXCEPT
Definition: Macros.h:1418
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:66
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
SelfAdjointView< typename MatrixType::TransposeReturnType, TransposeMode > TransposeReturnType
Array< int, Dynamic, 1 > v
ref_selector< MatrixType >::non_const_type MatrixTypeNested
RealScalar alpha
RealScalar s
Expression of a selfadjoint matrix from a triangular part of a dense matrix.
#define EIGEN_CONSTEXPR
Definition: Macros.h:787
EIGEN_DEVICE_FUNC SelfAdjointView(MatrixType &matrix)
#define EIGEN_DEVICE_FUNC
Definition: Macros.h:976
MatrixTypeNested m_matrix
SelfAdjointView< const MatrixConjugateReturnType, UpLo > ConjugateReturnType
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Expression of a triangular part in a matrix.
SelfAdjointView< const typename MatrixType::ConstTransposeReturnType, TransposeMode > ConstTransposeReturnType
m col(1)
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
remove_all< MatrixTypeNested >::type MatrixTypeNestedCleaned
internal::conditional< NumTraits< Scalar >::IsComplex, const CwiseUnaryOp< internal::scalar_conjugate_op< Scalar >, const Derived >, const Derived &>::type ConjugateReturnType
#define eigen_internal_assert(x)
Definition: Macros.h:1043
EIGEN_DEVICE_FUNC const ConstTransposeReturnType transpose() const
Generic expression where a coefficient-wise unary operator is applied to an expression.
Definition: CwiseUnaryOp.h:55
Map< Matrix< T, Dynamic, Dynamic, ColMajor >, 0, OuterStride<> > matrix(T *data, int rows, int cols, int stride)
The matrix class, also used for vectors and row-vectors.
EIGEN_DEVICE_FUNC internal::conditional<(TriMode &(Upper|Lower))==(UpLo &(Upper|Lower)), TriangularView< MatrixType, TriMode >, TriangularView< typename MatrixType::AdjointReturnType, TriMode > >::type triangularView() const
EIGEN_DEVICE_FUNC const Product< SelfAdjointView, OtherDerived > operator*(const MatrixBase< OtherDerived > &rhs) const
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const unsigned int LinearAccessBit
Definition: Constants.h:130
friend EIGEN_DEVICE_FUNC const Product< OtherDerived, SelfAdjointView > operator*(const MatrixBase< OtherDerived > &lhs, const SelfAdjointView &rhs)
MatrixTypeNestedCleaned NestedExpression
Definition: pytypes.h:1370
EIGEN_DEVICE_FUNC Scalar coeff(Index row, Index col) const
EIGEN_DEVICE_FUNC Scalar & coeffRef(Index row, Index col)


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:35:41