1#ifndef STAN_MATH_MIX_FUNCTOR_HESSIAN_TIMES_VECTOR_HPP
2#define STAN_MATH_MIX_FUNCTOR_HESSIAN_TIMES_VECTOR_HPP
15 const F& f,
const Eigen::Matrix<double, Eigen::Dynamic, 1>& x,
16 const Eigen::Matrix<double, Eigen::Dynamic, 1>& v,
double& fx,
17 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) {
28 var grad_fx_var_dot_v;
31 grad(grad_fx_var_dot_v.vi_);
33 for (
int i = 0; i < x.size(); ++i) {
34 Hv(i) = x_var(i).adj();
38template <
typename T,
typename F,
typename EigVec,
42 const Eigen::Matrix<T, Eigen::Dynamic, 1>& x,
43 const EigVec& v, T& fx,
44 Eigen::Matrix<T, Eigen::Dynamic, 1>& Hv) {
46 Matrix<T, Eigen::Dynamic, 1>
grad;
47 Matrix<T, Eigen::Dynamic, Eigen::Dynamic> H;
57template <
typename F,
typename XAdj,
typename XVec,
typename VVec,
63 const Eigen::Index x_size = x.size();
64 Eigen::Matrix<var, Eigen::Dynamic, 1> x_var = std::forward<XVec>(x);
65 Eigen::Matrix<fvar<var>, Eigen::Dynamic, 1> x_fvar(x_size);
66 for (Eigen::Index i = 0; i < x_size; i++) {
A class following the RAII idiom to start and recover nested autodiff scopes.
require_t< is_eigen_vector< std::decay_t< T > > > require_eigen_vector_t
Require type satisfies is_eigen_vector.
require_all_t< is_eigen_vector< std::decay_t< Types > >... > require_all_eigen_vector_t
Require all of the types satisfy is_eigen_vector.
require_t< is_stan_scalar< std::decay_t< T > > > require_stan_scalar_t
Require type satisfies is_stan_scalar.
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)
void hessian(const F &f, const Eigen::Matrix< T, Eigen::Dynamic, 1 > &x, T &fx, Eigen::Matrix< T, Eigen::Dynamic, 1 > &grad, Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > &H)
Calculate the value, the gradient, and the Hessian, of the specified function at the specified argume...
static void grad()
Compute the gradient for all variables starting from the end of the AD tape.
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)
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...
Scalar d_
The tangent (derivative) of this variable.
This template class represents scalars used in forward-mode automatic differentiation,...