Stan Math Library  2.14.0
reverse mode automatic differentiation
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 
7 #include <boost/shared_ptr.hpp>
8 
9 namespace stan {
10  namespace math {
11 
12  template <typename T, int R, int C>
13  class LDLT_factor;
14 
55  template <int R, int C, typename T>
56  class LDLT_factor<T, R, C> {
57  public:
59  : N_(0), ldltP_(new Eigen::LDLT<Eigen::Matrix<T, R, C> >()) {}
60 
61  explicit LDLT_factor(const Eigen::Matrix<T, R, C>& A)
62  : N_(0), ldltP_(new Eigen::LDLT<Eigen::Matrix<T, R, C> >()) {
63  compute(A);
64  }
65 
66  inline void compute(const Eigen::Matrix<T, R, C>& A) {
67  check_square("LDLT_factor", "A", A);
68  N_ = A.rows();
69  ldltP_->compute(A);
70  }
71 
72  inline bool success() const {
73  if (ldltP_->info() != Eigen::Success)
74  return false;
75  if (!(ldltP_->isPositive()))
76  return false;
77  Eigen::Matrix<T, Eigen::Dynamic, 1> ldltP_diag(ldltP_->vectorD());
78  for (int i = 0; i < ldltP_diag.size(); ++i)
79  if (ldltP_diag(i) <= 0 || is_nan(ldltP_diag(i)))
80  return false;
81  return true;
82  }
83 
84  inline T log_abs_det() const {
85  return ldltP_->vectorD().array().log().sum();
86  }
87 
88  inline void inverse(Eigen::Matrix<T, R, C>& invA) const {
89  invA.setIdentity(N_);
90  ldltP_->solveInPlace(invA);
91  }
92 
93 #if EIGEN_VERSION_AT_LEAST(3, 3, 0)
94  template <typename Rhs>
95  inline const
96  Eigen::Solve<Eigen::LDLT<Eigen::Matrix<T, R, C> >, Rhs>
97  solve(const Eigen::MatrixBase<Rhs>& b) const {
98  return ldltP_->solve(b);
99  }
100 #else
101  template <typename Rhs>
102  inline const
103  Eigen::internal::solve_retval<Eigen::LDLT< Eigen::Matrix<T, R, C> >, Rhs>
104  solve(const Eigen::MatrixBase<Rhs>& b) const {
105  return ldltP_->solve(b);
106  }
107 #endif
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 }
134 #endif
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>.
LDLT_factor(const Eigen::Matrix< T, R, C > &A)
Definition: LDLT_factor.hpp:61
(Expert) Numerical traits for algorithmic differentiation variables.
const Eigen::internal::solve_retval< Eigen::LDLT< Eigen::Matrix< T, R, C > >, Rhs > solve(const Eigen::MatrixBase< Rhs > &b) const
void inverse(Eigen::Matrix< T, R, C > &invA) const
Definition: LDLT_factor.hpp:88
Eigen::LDLT< Eigen::Matrix< T, R, C > > matrixLDLT() const
Eigen::Matrix< T, R, C > solveRight(const Eigen::Matrix< T, R, C > &B) const
void check_square(const char *function, const char *name, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y)
Check if the specified matrix is square.
boost::shared_ptr< Eigen::LDLT< Eigen::Matrix< T, R, C > > > ldltP_
Eigen::Matrix< T, Eigen::Dynamic, 1 > vectorD() const
void compute(const Eigen::Matrix< T, R, C > &A)
Definition: LDLT_factor.hpp:66
int is_nan(const fvar< T > &x)
Returns 1 if the input&#39;s value is NaN and 0 otherwise.
Definition: is_nan.hpp:21
Eigen::Matrix< T, C, R > transpose(const Eigen::Matrix< T, R, C > &m)
Definition: transpose.hpp:12
int N_

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