BasicPreconditioners.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) 2011 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_BASIC_PRECONDITIONERS_H
11 #define EIGEN_BASIC_PRECONDITIONERS_H
12 
13 namespace Eigen {
14 
32 template <typename _Scalar>
34 {
35  typedef _Scalar Scalar;
37  typedef typename Vector::Index Index;
38 
39  public:
40  // this typedef is only to export the scalar type and compile-time dimensions to solve_retval
42 
44 
45  template<typename MatType>
46  DiagonalPreconditioner(const MatType& mat) : m_invdiag(mat.cols())
47  {
48  compute(mat);
49  }
50 
51  Index rows() const { return m_invdiag.size(); }
52  Index cols() const { return m_invdiag.size(); }
53 
54  template<typename MatType>
56  {
57  return *this;
58  }
59 
60  template<typename MatType>
61  DiagonalPreconditioner& factorize(const MatType& mat)
62  {
63  m_invdiag.resize(mat.cols());
64  for(int j=0; j<mat.outerSize(); ++j)
65  {
66  typename MatType::InnerIterator it(mat,j);
67  while(it && it.index()!=j) ++it;
68  if(it && it.index()==j)
69  m_invdiag(j) = Scalar(1)/it.value();
70  else
71  m_invdiag(j) = 0;
72  }
73  m_isInitialized = true;
74  return *this;
75  }
76 
77  template<typename MatType>
78  DiagonalPreconditioner& compute(const MatType& mat)
79  {
80  return factorize(mat);
81  }
82 
83  template<typename Rhs, typename Dest>
84  void _solve(const Rhs& b, Dest& x) const
85  {
86  x = m_invdiag.array() * b.array() ;
87  }
88 
89  template<typename Rhs> inline const internal::solve_retval<DiagonalPreconditioner, Rhs>
90  solve(const MatrixBase<Rhs>& b) const
91  {
92  eigen_assert(m_isInitialized && "DiagonalPreconditioner is not initialized.");
93  eigen_assert(m_invdiag.size()==b.rows()
94  && "DiagonalPreconditioner::solve(): invalid number of rows of the right hand side matrix b");
95  return internal::solve_retval<DiagonalPreconditioner, Rhs>(*this, b.derived());
96  }
97 
98  protected:
99  Vector m_invdiag;
101 };
102 
103 namespace internal {
104 
105 template<typename _MatrixType, typename Rhs>
106 struct solve_retval<DiagonalPreconditioner<_MatrixType>, Rhs>
107  : solve_retval_base<DiagonalPreconditioner<_MatrixType>, Rhs>
108 {
110  EIGEN_MAKE_SOLVE_HELPERS(Dec,Rhs)
111 
112  template<typename Dest> void evalTo(Dest& dst) const
113  {
114  dec()._solve(rhs(),dst);
115  }
116 };
117 
118 }
119 
126 {
127  public:
128 
130 
131  template<typename MatrixType>
133 
134  template<typename MatrixType>
136 
137  template<typename MatrixType>
138  IdentityPreconditioner& factorize(const MatrixType& ) { return *this; }
139 
140  template<typename MatrixType>
141  IdentityPreconditioner& compute(const MatrixType& ) { return *this; }
142 
143  template<typename Rhs>
144  inline const Rhs& solve(const Rhs& b) const { return b; }
145 };
146 
147 } // end namespace Eigen
148 
149 #endif // EIGEN_BASIC_PRECONDITIONERS_H
void _solve(const Rhs &b, Dest &x) const
A preconditioner based on the digonal entries.
const internal::solve_retval< DiagonalPreconditioner, Rhs > solve(const MatrixBase< Rhs > &b) const
internal::traits< Derived >::Index Index
iterative scaling algorithm to equilibrate rows and column norms in matrices
Definition: matrix.hpp:471
IdentityPreconditioner & analyzePattern(const MatrixType &)
DiagonalPreconditioner & factorize(const MatType &mat)
DiagonalPreconditioner & compute(const MatType &mat)
EIGEN_STRONG_INLINE void resize(Index nbRows, Index nbCols)
DiagonalPreconditioner & analyzePattern(const MatType &)
void rhs(const real_t *x, real_t *f)
IdentityPreconditioner & compute(const MatrixType &)
Matrix< Scalar, Dynamic, 1 > Vector
const Rhs & solve(const Rhs &b) const
Matrix< Scalar, Dynamic, Dynamic > MatrixType
DiagonalPreconditioner(const MatType &mat)
A naive preconditioner which approximates any matrix as the identity matrix.
#define EIGEN_MAKE_SOLVE_HELPERS(DecompositionType, Rhs)
Definition: Solve.h:61
#define eigen_assert(x)
IdentityPreconditioner & factorize(const MatrixType &)
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
IdentityPreconditioner(const MatrixType &)


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Mon Jun 10 2019 12:34:28