1 #ifndef STAN_MATH_REV_MAT_FUN_MULTIPLY_HPP
2 #define STAN_MATH_REV_MAT_FUN_MULTIPLY_HPP
12 #include <boost/utility/enable_if.hpp>
13 #include <boost/type_traits.hpp>
14 #include <boost/math/tools/promotion.hpp>
26 template <
typename T1,
typename T2>
29 (boost::is_scalar<T1>::value || boost::is_same<T1, var>::value)
30 && (boost::is_scalar<T2>::value || boost::is_same<T2, var>::value),
31 typename boost::math::tools::promote_args<T1, T2>::type>::type
42 template<
typename T1,
typename T2,
int R2,
int C2>
43 inline Eigen::Matrix<var, R2, C2>
44 multiply(
const T1& c,
const Eigen::Matrix<T2, R2, C2>& m) {
56 template<
typename T1,
int R1,
int C1,
typename T2>
57 inline Eigen::Matrix<var, R1, C1>
58 multiply(
const Eigen::Matrix<T1, R1, C1>& m,
const T2& c) {
72 template<
typename T1,
int R1,
int C1,
typename T2,
int R2,
int C2>
74 boost::enable_if_c< boost::is_same<T1, var>::value ||
75 boost::is_same<T2, var>::value,
76 Eigen::Matrix<var, R1, C2> >::type
78 const Eigen::Matrix<T2, R2, C2>& m2) {
82 Eigen::Matrix<var, R1, C2> result(m1.rows(), m2.cols());
83 for (
int i = 0; i < m1.rows(); i++) {
84 typename Eigen::Matrix<T1, R1, C1>::ConstRowXpr crow(m1.row(i));
85 for (
int j = 0; j < m2.cols(); j++) {
86 typename Eigen::Matrix<T2, R2, C2>::ConstColXpr ccol(m2.col(j));
89 result(i, j) =
var(
new dot_product_vari<T1, T2>(crow, ccol));
91 dot_product_vari<T1, T2> *v2
92 =
static_cast<dot_product_vari<T1, T2>*
>(result(0, j).vi_);
94 =
var(
new dot_product_vari<T1, T2>(crow, ccol, NULL, v2));
98 dot_product_vari<T1, T2> *v1
99 =
static_cast<dot_product_vari<T1, T2>*
>(result(i, 0).vi_);
101 =
var(
new dot_product_vari<T1, T2>(crow, ccol, v1, NULL));
103 dot_product_vari<T1, T2> *v1
104 =
static_cast<dot_product_vari<T1, T2>*
>(result(i, 0).vi_);
105 dot_product_vari<T1, T2> *v2
106 =
static_cast<dot_product_vari<T1, T2>*
>(result(0, j).vi_);
108 =
var(
new dot_product_vari<T1, T2>(crow, ccol, v1, v2));
125 template <
typename T1,
int C1,
typename T2,
int R2>
127 boost::enable_if_c< boost::is_same<T1, var>::value ||
128 boost::is_same<T2, var>::value, var >::type
130 const Eigen::Matrix<T2, R2, 1>& v) {
131 if (rv.size() != v.size())
Eigen::Matrix< fvar< T >, R1, C1 > multiply(const Eigen::Matrix< fvar< T >, R1, C1 > &m, const fvar< T > &c)
Independent (input) and dependent (output) variables for gradients.
bool check_multiplicable(const char *function, const char *name1, const T1 &y1, const char *name2, const T2 &y2)
Return true if the matrices can be multiplied.
fvar< T > dot_product(const Eigen::Matrix< fvar< T >, R1, C1 > &v1, const Eigen::Matrix< fvar< T >, R2, C2 > &v2)
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
Throw a domain error with a consistently formatted message.
std::vector< var > to_var(const std::vector< double > &v)
Converts argument to an automatic differentiation variable.