1 #ifndef STAN_MATH_MIX_MAT_FUNCTOR_HESSIAN_TIMES_VECTOR_HPP
2 #define STAN_MATH_MIX_MAT_FUNCTOR_HESSIAN_TIMES_VECTOR_HPP
16 const Eigen::Matrix<double, Eigen::Dynamic, 1>& x,
17 const Eigen::Matrix<double, Eigen::Dynamic, 1>& v,
19 Eigen::Matrix<double, Eigen::Dynamic, 1>& Hv) {
23 Matrix<var, Eigen::Dynamic, 1> x_var(x.size());
24 for (
int i = 0; i < x_var.size(); ++i)
27 var grad_fx_var_dot_v;
32 for (
int i = 0; i < x.size(); ++i)
33 Hv(i) = x_var(i).adj();
34 }
catch (
const std::exception&
e) {
40 template <
typename T,
typename F>
43 const Eigen::Matrix<T, Eigen::Dynamic, 1>& x,
44 const Eigen::Matrix<T, Eigen::Dynamic, 1>& v,
46 Eigen::Matrix<T, Eigen::Dynamic, 1>& Hv) {
48 Matrix<T, Eigen::Dynamic, 1>
grad;
49 Matrix<T, Eigen::Dynamic, Eigen::Dynamic> H;
void hessian(const F &f, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &x, double &fx, Eigen::Matrix< double, Eigen::Dynamic, 1 > &grad, Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > &H)
Calculate the value, the gradient, and the Hessian, of the specified function at the specified argume...
Independent (input) and dependent (output) variables for gradients.
void hessian_times_vector(const F &f, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &x, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &v, double &fx, Eigen::Matrix< double, Eigen::Dynamic, 1 > &Hv)
static void grad(vari *vi)
Compute the gradient for all variables starting from the specified root variable implementation.
vari * vi_
Pointer to the implementation of this variable.
double e()
Return the base of the natural logarithm.
static void recover_memory_nested()
Recover only the memory used for the top nested call.
double val() const
Return the value of this variable.
void gradient_dot_vector(const F &f, const Eigen::Matrix< T1, Eigen::Dynamic, 1 > &x, const Eigen::Matrix< T2, Eigen::Dynamic, 1 > &v, T1 &fx, T1 &grad_fx_dot_v)
static void start_nested()
Record the current position so that recover_memory_nested() can find it.