1 #ifndef STAN_MATH_REV_MAT_FUN_MDIVIDE_LEFT_HPP 2 #define STAN_MATH_REV_MAT_FUN_MDIVIDE_LEFT_HPP 16 template <
int R1,
int C1,
int R2,
int C2>
28 const Eigen::Matrix<var, R2, C2> &B)
32 A_(reinterpret_cast<double *>(
35 C_(reinterpret_cast<double *>(
38 variRefA_(reinterpret_cast<
vari **>(
41 variRefB_(reinterpret_cast<
vari **>(
44 variRefC_(reinterpret_cast<
vari **>(
53 variRefA_[pos] = A(i, j).vi_;
54 A_[pos++] = A(i, j).val();
61 variRefB_[pos] = B(i, j).vi_;
62 C_[pos++] = B(i, j).val();
66 Matrix<double, R1, C2> C(M_, N_);
67 C = Map<Matrix<double, R1, C2> >(
C_,
M_,
N_);
69 C = Map<Matrix<double, R1, C1> >(
A_,
M_,
M_).colPivHouseholderQr().solve(C);
75 variRefC_[pos] =
new vari(C_[pos],
false);
84 Eigen::Matrix<double, R1, C1> adjA(M_, M_);
85 Eigen::Matrix<double, R2, C2> adjB(M_, N_);
86 Eigen::Matrix<double, R1, C2> adjC(M_, N_);
89 for (
size_type j = 0; j < adjC.cols(); j++)
90 for (
size_type i = 0; i < adjC.rows(); i++)
91 adjC(i, j) = variRefC_[pos++]->
adj_;
93 adjB = Map<Matrix<double, R1, C1> >(
A_,
M_,
M_)
95 .colPivHouseholderQr()
101 for (
size_type j = 0; j < adjA.cols(); j++)
102 for (
size_type i = 0; i < adjA.rows(); i++)
103 variRefA_[pos++]->
adj_ += adjA(i, j);
106 for (
size_type j = 0; j < adjB.cols(); j++)
107 for (
size_type i = 0; i < adjB.rows(); i++)
108 variRefB_[pos++]->
adj_ += adjB(i, j);
112 template <
int R1,
int C1,
int R2,
int C2>
123 const Eigen::Matrix<var, R2, C2> &B)
127 A_(reinterpret_cast<double *>(
130 C_(reinterpret_cast<double *>(
133 variRefB_(reinterpret_cast<
vari **>(
136 variRefC_(reinterpret_cast<
vari **>(
152 variRefB_[pos] = B(i, j).vi_;
153 C_[pos++] = B(i, j).val();
157 Matrix<double, R1, C2> C(M_, N_);
158 C = Map<Matrix<double, R1, C2> >(
C_,
M_,
N_);
160 C = Map<Matrix<double, R1, C1> >(
A_,
M_,
M_).colPivHouseholderQr().solve(C);
166 variRefC_[pos] =
new vari(C_[pos],
false);
175 Eigen::Matrix<double, R2, C2> adjB(M_, N_);
176 Eigen::Matrix<double, R1, C2> adjC(M_, N_);
179 for (
size_type j = 0; j < adjC.cols(); j++)
180 for (
size_type i = 0; i < adjC.rows(); i++)
181 adjC(i, j) = variRefC_[pos++]->
adj_;
183 adjB = Map<Matrix<double, R1, C1> >(
A_,
M_,
M_)
185 .colPivHouseholderQr()
189 for (
size_type j = 0; j < adjB.cols(); j++)
190 for (
size_type i = 0; i < adjB.rows(); i++)
191 variRefB_[pos++]->
adj_ += adjB(i, j);
195 template <
int R1,
int C1,
int R2,
int C2>
206 const Eigen::Matrix<double, R2, C2> &B)
210 A_(reinterpret_cast<double *>(
213 C_(reinterpret_cast<double *>(
216 variRefA_(reinterpret_cast<
vari **>(
219 variRefC_(reinterpret_cast<
vari **>(
228 variRefA_[pos] = A(i, j).vi_;
229 A_[pos++] = A(i, j).val();
233 Matrix<double, R1, C2> C(M_, N_);
234 C = Map<Matrix<double, R1, C1> >(
A_,
M_,
M_).colPivHouseholderQr().solve(B);
240 variRefC_[pos] =
new vari(C_[pos],
false);
249 Eigen::Matrix<double, R1, C1> adjA(M_, M_);
250 Eigen::Matrix<double, R1, C2> adjC(M_, N_);
253 for (
size_type j = 0; j < adjC.cols(); j++)
254 for (
size_type i = 0; i < adjC.rows(); i++)
255 adjC(i, j) = variRefC_[pos++]->
adj_;
258 adjA = -Map<Matrix<double, R1, C1> >(
A_,
M_,
M_)
260 .colPivHouseholderQr()
262 * Map<Matrix<double, R1, C2> >(C_, M_, N_).
transpose());
265 for (
size_type j = 0; j < adjA.cols(); j++)
266 for (
size_type i = 0; i < adjA.rows(); i++)
267 variRefA_[pos++]->
adj_ += adjA(i, j);
272 template <
int R1,
int C1,
int R2,
int C2>
274 const Eigen::Matrix<var, R1, C1> &A,
const Eigen::Matrix<var, R2, C2> &b) {
275 Eigen::Matrix<var, R1, C2> res(b.rows(), b.cols());
288 for (
size_type j = 0; j < res.cols(); j++)
289 for (
size_type i = 0; i < res.rows(); i++)
290 res(i, j).vi_ = baseVari->
variRefC_[pos++];
295 template <
int R1,
int C1,
int R2,
int C2>
297 const Eigen::Matrix<var, R1, C1> &A,
298 const Eigen::Matrix<double, R2, C2> &b) {
299 Eigen::Matrix<var, R1, C2> res(b.rows(), b.cols());
312 for (
size_type j = 0; j < res.cols(); j++)
313 for (
size_type i = 0; i < res.rows(); i++)
314 res(i, j).vi_ = baseVari->
variRefC_[pos++];
319 template <
int R1,
int C1,
int R2,
int C2>
321 const Eigen::Matrix<double, R1, C1> &A,
322 const Eigen::Matrix<var, R2, C2> &b) {
323 Eigen::Matrix<var, R1, C2> res(b.rows(), b.cols());
336 for (
size_type j = 0; j < res.cols(); j++)
337 for (
size_type i = 0; i < res.rows(); i++)
338 res(i, j).vi_ = baseVari->
variRefC_[pos++];
Eigen::Matrix< fvar< T >, R1, C2 > mdivide_left(const Eigen::Matrix< fvar< T >, R1, C1 > &A, const Eigen::Matrix< fvar< T >, R2, C2 > &b)
int rows(const Eigen::Matrix< T, R, C > &m)
Return the number of rows in the specified matrix, vector, or row vector.
mdivide_left_vv_vari(const Eigen::Matrix< var, R1, C1 > &A, const Eigen::Matrix< var, R2, C2 > &B)
void check_square(const char *function, const char *name, const matrix_cl &y)
Check if the matrix_cl is square.
The variable implementation base class.
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
mdivide_left_vd_vari(const Eigen::Matrix< var, R1, C1 > &A, const Eigen::Matrix< double, R2, C2 > &B)
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic >::Index size_type
Type for sizes and indexes in an Eigen matrix with double e.
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
int cols(const Eigen::Matrix< T, R, C > &m)
Return the number of columns in the specified matrix, vector, or row vector.
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.
matrix_cl transpose(const matrix_cl &src)
Takes the transpose of the matrix on the OpenCL device.
vari(double x)
Construct a variable implementation from a value.
virtual void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
double adj_
The adjoint of this variable, which is the partial derivative of this variable with respect to the ro...
mdivide_left_dv_vari(const Eigen::Matrix< double, R1, C1 > &A, const Eigen::Matrix< var, R2, C2 > &B)
This struct always provides access to the autodiff stack using the singleton pattern.