1 #ifndef STAN_MATH_FWD_MAT_FUN_UNIT_VECTOR_CONSTRAIN_HPP
2 #define STAN_MATH_FWD_MAT_FUN_UNIT_VECTOR_CONSTRAIN_HPP
18 template <
typename T,
int R,
int C>
19 inline Eigen::Matrix<fvar<T>, R, C>
24 Matrix<T, R, C> y_t(y.size());
25 for (
int k = 0; k < y.size(); ++k)
26 y_t.coeffRef(k) = y.coeff(k).val_;
28 Matrix<T, R, C> unit_vector_y_t
30 Matrix<fvar<T>, R, C> unit_vector_y(y.size());
31 for (
int k = 0; k < y.size(); ++k)
32 unit_vector_y.coeffRef(k).val_ = unit_vector_y_t.coeff(k);
34 const T squared_norm =
dot_self(y_t);
35 const T norm =
sqrt(squared_norm);
36 const T inv_norm =
inv(norm);
37 Matrix<T, Eigen::Dynamic, Eigen::Dynamic> J
41 for (
int m = 0; m < y.size(); ++m) {
42 J.coeffRef(m, m) += inv_norm;
44 for (
int k = 0; k < y.size(); ++k) {
46 unit_vector_y.coeffRef(k).d_ = J.coeff(k, m);
52 template <
typename T,
int R,
int C>
53 inline Eigen::Matrix<fvar<T>, R, C>
56 lp -= 0.5 * squared_norm;
fvar< T > sqrt(const fvar< T > &x)
Eigen::Matrix< fvar< T >, R, R > tcrossprod(const Eigen::Matrix< fvar< T >, R, C > &m)
stan::return_type< T1, T2 >::type divide(const T1 &v, const T2 &c)
fvar< T > dot_self(const Eigen::Matrix< fvar< T >, R, C > &v)
Eigen::Matrix< fvar< T >, R, C > unit_vector_constrain(const Eigen::Matrix< fvar< T >, R, C > &y)
fvar< T > inv(const fvar< T > &x)