Stan Math Library  2.8.0
reverse mode automatic differentiation
 All Classes Namespaces Files Functions Variables Typedefs Enumerator Friends Macros Groups
LDLT_factor.hpp
Go to the documentation of this file.
1 #ifndef STAN_MATH_PRIM_MAT_FUN_LDLT_FACTOR_HPP
2 #define STAN_MATH_PRIM_MAT_FUN_LDLT_FACTOR_HPP
3 
5 #include <boost/shared_ptr.hpp>
8 
9 namespace stan {
10 
11  namespace math {
12 
13  // This class is conceptually similar to the corresponding Eigen class
14  // Any spd matrix A can be decomposed as LDL' where L is unit
15  // lower-triangular and D is diagonal with positive diagonal elements
16 
17  template<typename T, int R, int C>
18  class LDLT_factor;
19 
57  template<int R, int C, typename T>
58  class LDLT_factor<T, R, C> {
59  public:
61  : N_(0), _ldltP(new Eigen::LDLT< Eigen::Matrix<T, R, C> >()) {}
62 
63  explicit LDLT_factor(const Eigen::Matrix<T, R, C> &A)
64  : N_(0), _ldltP(new Eigen::LDLT< Eigen::Matrix<T, R, C> >()) {
65  compute(A);
66  }
67 
68  inline void compute(const Eigen::Matrix<T, R, C> &A) {
69  stan::math::check_square("LDLT_factor", "A", A);
70  N_ = A.rows();
71  _ldltP->compute(A);
72  }
73 
74  inline bool success() const {
75  using stan::math::is_nan;
76  // bool ret;
77  // ret = _ldltP->info() == Eigen::Success;
78  // ret = ret && _ldltP->isPositive();
79  // ret = ret && (_ldltP->vectorD().array() > 0).all();
80  // return ret;
81 
82  if (_ldltP->info() != Eigen::Success)
83  return false;
84  if (!(_ldltP->isPositive()))
85  return false;
86  Eigen::Matrix<T, Eigen::Dynamic, 1> ldltP_diag(_ldltP->vectorD());
87  for (int i = 0; i < ldltP_diag.size(); ++i)
88  if (ldltP_diag(i) <= 0 || is_nan(ldltP_diag(i)))
89  return false;
90  return true;
91  }
92 
93  inline T log_abs_det() const {
94  return _ldltP->vectorD().array().log().sum();
95  }
96 
97  inline void inverse(Eigen::Matrix<T, R, C> &invA) const {
98  invA.setIdentity(N_);
99  _ldltP->solveInPlace(invA);
100  }
101 
102  template<typename Rhs>
103  inline const
104  Eigen::internal::solve_retval<Eigen::LDLT< Eigen::Matrix<T, R, C> >, Rhs>
105  solve(const Eigen::MatrixBase<Rhs>& b) const {
106  return _ldltP->solve(b);
107  }
108 
109  inline Eigen::Matrix<T, R, C>
110  solveRight(const Eigen::Matrix<T, R, C> &B) const {
111  return _ldltP->solve(B.transpose()).transpose();
112  }
113 
114  inline Eigen::Matrix<T, Eigen::Dynamic, 1> vectorD() const {
115  return _ldltP->vectorD();
116  }
117 
118  inline Eigen::LDLT<Eigen::Matrix<T, R, C> > matrixLDLT() const {
119  return _ldltP->matrixLDLT();
120  }
121 
122  inline size_t rows() const { return N_; }
123  inline size_t cols() const { return N_; }
124 
125  typedef size_t size_type;
126  typedef double value_type;
127 
128  size_t N_;
129  boost::shared_ptr< Eigen::LDLT< Eigen::Matrix<T, R, C> > > _ldltP;
130  };
131  }
132 }
133 #endif
void inverse(Eigen::Matrix< T, R, C > &invA) const
Definition: LDLT_factor.hpp:97
const Eigen::internal::solve_retval< Eigen::LDLT< Eigen::Matrix< T, R, C > >, Rhs > solve(const Eigen::MatrixBase< Rhs > &b) const
boost::shared_ptr< Eigen::LDLT< Eigen::Matrix< T, R, C > > > _ldltP
LDLT_factor(const Eigen::Matrix< T, R, C > &A)
Definition: LDLT_factor.hpp:63
Eigen::Matrix< T, R, C > solveRight(const Eigen::Matrix< T, R, C > &B) const
Eigen::LDLT< Eigen::Matrix< T, R, C > > matrixLDLT() const
boost::shared_ptr< Eigen::LDLT< Eigen::Matrix< double, R1, C1 > > > _ldltP
This share_ptr is used to prevent copying the LDLT factorizations for mdivide_left_ldlt(ldltA, b) when ldltA is a LDLT_factor<double>.
Eigen::Matrix< T, Eigen::Dynamic, 1 > vectorD() const
void compute(const Eigen::Matrix< T, R, C > &A)
Definition: LDLT_factor.hpp:68
int is_nan(const fvar< T > &x)
Returns 1 if the input's value is NaN and 0 otherwise.
Definition: is_nan.hpp:22
Eigen::Matrix< T, C, R > transpose(const Eigen::Matrix< T, R, C > &m)
Definition: transpose.hpp:12
bool check_square(const char *function, const char *name, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y)
Return true if the specified matrix is square.
int N_

     [ Stan Home Page ] © 2011–2015, Stan Development Team.