1 #ifndef STAN_MATH_REV_MAT_FUN_MDIVIDE_LEFT_LDLT_HPP 2 #define STAN_MATH_REV_MAT_FUN_MDIVIDE_LEFT_LDLT_HPP 14 template <
int R1,
int C1,
int R2,
int C2>
24 boost::shared_ptr<Eigen::LDLT<Eigen::Matrix<double, R1, C1> > >
ldltP_;
25 Eigen::Matrix<double, R2, C2>
C_;
38 template <
int R1,
int C1,
int R2,
int C2>
49 const Eigen::Matrix<var, R2, C2> &B)
53 variRefB_(reinterpret_cast<
vari **>(
56 variRefC_(reinterpret_cast<
vari **>(
60 alloc_ldlt_(A.alloc_) {
62 alloc_->
C_.resize(M_, N_);
63 for (
int j = 0; j < N_; j++) {
64 for (
int i = 0; i < M_; i++) {
65 variRefB_[pos] = B(i, j).vi_;
66 alloc_->
C_(i, j) = B(i, j).val();
71 alloc_ldlt_->
ldlt_.solveInPlace(alloc_->
C_);
74 for (
int j = 0; j < N_; j++) {
75 for (
int i = 0; i < M_; i++) {
76 variRefC_[pos] =
new vari(alloc_->
C_(i, j),
false);
83 Eigen::Matrix<double, R1, C1> adjA(M_, M_);
84 Eigen::Matrix<double, R2, C2> adjB(M_, N_);
87 for (
int j = 0; j < N_; j++)
88 for (
int i = 0; i < M_; i++)
89 adjB(i, j) = variRefC_[pos++]->
adj_;
91 alloc_ldlt_->
ldlt_.solveInPlace(adjB);
92 adjA.noalias() = -adjB * alloc_->
C_.transpose();
94 for (
int j = 0; j < M_; j++)
95 for (
int i = 0; i < M_; i++)
96 alloc_ldlt_->
variA_(i, j)->adj_ += adjA(i, j);
99 for (
int j = 0; j < N_; j++)
100 for (
int i = 0; i < M_; i++)
101 variRefB_[pos++]->adj_ += adjB(i, j);
115 template <
int R1,
int C1,
int R2,
int C2>
125 const Eigen::Matrix<var, R2, C2> &B)
129 variRefB_(reinterpret_cast<
vari **>(
132 variRefC_(reinterpret_cast<
vari **>(
140 alloc_->
C_.resize(M_, N_);
141 for (
int j = 0; j < N_; j++) {
142 for (
int i = 0; i < M_; i++) {
143 variRefB_[pos] = B(i, j).vi_;
144 alloc_->
C_(i, j) = B(i, j).val();
150 alloc_->
ldltP_->solveInPlace(alloc_->
C_);
153 for (
int j = 0; j < N_; j++) {
154 for (
int i = 0; i < M_; i++) {
155 variRefC_[pos] =
new vari(alloc_->
C_(i, j),
false);
162 Eigen::Matrix<double, R2, C2> adjB(M_, N_);
165 for (
int j = 0; j < adjB.cols(); j++)
166 for (
int i = 0; i < adjB.rows(); i++)
167 adjB(i, j) = variRefC_[pos++]->
adj_;
169 alloc_->
ldltP_->solveInPlace(adjB);
172 for (
int j = 0; j < adjB.cols(); j++)
173 for (
int i = 0; i < adjB.rows(); i++)
174 variRefB_[pos++]->adj_ += adjB(i, j);
188 template <
int R1,
int C1,
int R2,
int C2>
198 const Eigen::Matrix<double, R2, C2> &B)
202 variRefC_(reinterpret_cast<
vari **>(
206 alloc_ldlt_(A.alloc_) {
208 alloc_ldlt_->
ldlt_.solveInPlace(alloc_->
C_);
211 for (
int j = 0; j < N_; j++) {
212 for (
int i = 0; i < M_; i++) {
213 variRefC_[pos] =
new vari(alloc_->
C_(i, j),
false);
220 Eigen::Matrix<double, R1, C1> adjA(M_, M_);
221 Eigen::Matrix<double, R1, C2> adjC(M_, N_);
224 for (
int j = 0; j < adjC.cols(); j++)
225 for (
int i = 0; i < adjC.rows(); i++)
226 adjC(i, j) = variRefC_[pos++]->
adj_;
228 adjA = -alloc_ldlt_->
ldlt_.solve(adjC * alloc_->
C_.transpose());
230 for (
int j = 0; j < adjA.cols(); j++)
231 for (
int i = 0; i < adjA.rows(); i++)
232 alloc_ldlt_->
variA_(i, j)->adj_ += adjA(i, j);
244 template <
int R1,
int C1,
int R2,
int C2>
247 Eigen::Matrix<var, R1, C2> res(b.rows(), b.cols());
255 for (
int j = 0; j < res.cols(); j++)
256 for (
int i = 0; i < res.rows(); i++)
257 res(i, j).vi_ = baseVari->
variRefC_[pos++];
269 template <
int R1,
int C1,
int R2,
int C2>
272 Eigen::Matrix<var, R1, C2> res(b.rows(), b.cols());
280 for (
int j = 0; j < res.cols(); j++)
281 for (
int i = 0; i < res.rows(); i++)
282 res(i, j).vi_ = baseVari->
variRefC_[pos++];
294 template <
int R1,
int C1,
int R2,
int C2>
297 Eigen::Matrix<var, R1, C2> res(b.rows(), b.cols());
305 for (
int j = 0; j < res.cols(); j++)
306 for (
int i = 0; i < res.rows(); i++)
307 res(i, j).vi_ = baseVari->
variRefC_[pos++];
const LDLT_alloc< R1, C1 > * alloc_ldlt_
int rows(const Eigen::Matrix< T, R, C > &m)
Return the number of rows in the specified matrix, vector, or row vector.
Eigen::LDLT< Eigen::Matrix< double, R, C > > ldlt_
Eigen::Matrix< vari *, R, C > variA_
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
mdivide_left_ldlt_alloc< R1, C1, R2, C2 > * alloc_
Eigen::Matrix< fvar< T2 >, R1, C2 > mdivide_left_ldlt(const LDLT_factor< double, R1, C1 > &A, const Eigen::Matrix< fvar< T2 >, R2, C2 > &b)
Returns the solution of the system Ax=b given an LDLT_factor of A.
The variable implementation base class.
const LDLT_alloc< R1, C1 > * alloc_ldlt_
mdivide_left_ldlt_alloc< R1, C1, R2, C2 > * alloc_
LDLT_factor is a thin wrapper on Eigen::LDLT to allow for reusing factorizations and efficient autodi...
boost::shared_ptr< ldlt_t > ldltP_
mdivide_left_ldlt_alloc< R1, C1, R2, C2 > * alloc_
Eigen::Matrix< double, R2, C2 > C_
virtual ~mdivide_left_ldlt_alloc()
int cols(const Eigen::Matrix< T, R, C > &m)
Return the number of columns in the specified matrix, vector, or row vector.
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>.
void check_multiplicable(const char *function, const char *name1, const T1 &y1, const char *name2, const T2 &y2)
Check if the matrices can be multiplied.
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
mdivide_left_ldlt_vv_vari(const LDLT_factor< var, R1, C1 > &A, const Eigen::Matrix< var, R2, C2 > &B)
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
mdivide_left_ldlt_dv_vari(const LDLT_factor< double, R1, C1 > &A, const Eigen::Matrix< var, R2, C2 > &B)
A chainable_alloc is an object which is constructed and destructed normally but the memory lifespan i...
double adj_
The adjoint of this variable, which is the partial derivative of this variable with respect to the ro...
The vari for mdivide_left_ldlt(A, b) which handles the chain() call for all elements of the result...
This struct always provides access to the autodiff stack using the singleton pattern.
mdivide_left_ldlt_vd_vari(const LDLT_factor< var, R1, C1 > &A, const Eigen::Matrix< double, R2, C2 > &B)
The vari for mdivide_left_ldlt(A, b) which handles the chain() call for all elements of the result...
The vari for mdivide_left_ldlt(A, b) which handles the chain() call for all elements of the result...