1 #ifndef STAN_MATH_REV_MAT_FUN_CHOLESKY_DECOMPOSE_HPP
2 #define STAN_MATH_REV_MAT_FUN_CHOLESKY_DECOMPOSE_HPP
42 const Eigen::Matrix<double, -1, -1>& L_A)
50 size_t accum_i = accum;
54 size_t pos = j + accum_i;
80 using Eigen::RowMajor;
81 Matrix<double, -1, -1, RowMajor> adjL(
M_,
M_);
82 Matrix<double, -1, -1, RowMajor> LA(
M_,
M_);
83 Matrix<double, -1, -1, RowMajor> adjA(
M_,
M_);
94 for (
int i = M_ - 1; i >= 0; --i) {
95 for (
int j = i; j >= 0; --j) {
97 adjA.coeffRef(i, j) = 0.5 * adjL.coeff(i, j)
100 adjA.coeffRef(i, j) = adjL.coeff(i, j)
102 adjL.coeffRef(j, j) -= adjL.coeff(i, j)
103 * LA.coeff(i, j) / LA.coeff(j, j);
105 for (
int k = j - 1; k >=0; --k) {
106 adjL.coeffRef(i, k) -= adjA.coeff(i, j)
108 adjL.coeffRef(j, k) -= adjA.coeff(i, j)
130 Eigen::Matrix<var, -1, -1>
136 Eigen::LLT<Eigen::MatrixXd> L_factor
137 = L_A.selfadjointView<Eigen::Lower>().llt();
139 L_A = L_factor.matrixL();
148 Eigen::Matrix<
var, -1, -1> L(A.rows(), A.cols());
150 size_t accum_i = accum;
151 for (
size_type j = 0; j < L.cols(); ++j) {
152 for (
size_type i = j; i < L.cols(); ++i) {
154 size_t pos = j + accum_i;
158 L.coeffRef(k, j).vi_ = &dummy;
vari(const double x)
Construct a variable implementation from a value.
int rows(const Eigen::Matrix< T, R, C > &m)
Return the number of rows in the specified matrix, vector, or row vector.
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
cholesky_decompose_v_vari(const Eigen::Matrix< var,-1,-1 > &A, const Eigen::Matrix< double,-1,-1 > &L_A)
double value_of_rec(const fvar< T > &v)
Return the value of the specified variable.
The variable implementation base class.
bool check_symmetric(const char *function, const char *name, const Eigen::Matrix< T_y, Dynamic, Dynamic > &y)
Return true if the specified matrix is symmetric.
Independent (input) and dependent (output) variables for gradients.
const double val_
The value of this variable.
bool check_pos_definite(const char *function, const char *name, const Eigen::Matrix< T_y, Dynamic, Dynamic > &y)
Return true if the specified square, symmetric matrix is positive definite.
Empty struct for use in boost::condtional<is_constant_struct<T1>::value, T1, dummy>::type as false co...
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic >::Index size_type
Type for sizes and indexes in an Eigen matrix with double e.
vari * vi_
Pointer to the implementation of this variable.
Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > cholesky_decompose(const Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > &m)
Return the lower-triangular Cholesky factor (i.e., matrix square root) of the specified square...
double adj_
The adjoint of this variable, which is the partial derivative of this variable with respect to the ro...
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.