Loading [MathJax]/extensions/TeX/AMSsymbols.js
Automatic Differentiation
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
hessian_times_vector.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_MIX_FUNCTOR_HESSIAN_TIMES_VECTOR_HPP
2#define STAN_MATH_MIX_FUNCTOR_HESSIAN_TIMES_VECTOR_HPP
3
7#include <stdexcept>
8#include <vector>
9
10namespace stan {
11namespace math {
12
13template <typename F>
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) {
18 using Eigen::Matrix;
19
20 // Run nested autodiff in this scope
22
23 Matrix<var, Eigen::Dynamic, 1> x_var(x.size());
24 for (int i = 0; i < x_var.size(); ++i) {
25 x_var(i) = x(i);
26 }
27 var fx_var;
28 var grad_fx_var_dot_v;
29 gradient_dot_vector(f, x_var, v, fx_var, grad_fx_var_dot_v);
30 fx = fx_var.val();
31 grad(grad_fx_var_dot_v.vi_);
32 Hv.resize(x.size());
33 for (int i = 0; i < x.size(); ++i) {
34 Hv(i) = x_var(i).adj();
35 }
36}
37
38template <typename T, typename F, typename EigVec,
40 require_stan_scalar_t<T>* = nullptr>
41void hessian_times_vector(const F& f,
42 const Eigen::Matrix<T, Eigen::Dynamic, 1>& x,
43 const EigVec& v, T& fx,
44 Eigen::Matrix<T, Eigen::Dynamic, 1>& Hv) {
45 using Eigen::Matrix;
46 Matrix<T, Eigen::Dynamic, 1> grad;
47 Matrix<T, Eigen::Dynamic, Eigen::Dynamic> H;
48 hessian(f, x, fx, grad, H);
49 Hv = H * v;
50}
51
57template <typename F, typename XAdj, typename XVec, typename VVec,
58 typename... Args,
60inline void hessian_times_vector(const F& f, XAdj& x_adj, XVec&& x, VVec&& v,
61 Args&&... args) {
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++) {
67 x_fvar(i) = fvar<var>(x_var(i), v(i));
68 }
69 fvar<var> fx_fvar = f(x_fvar, args...);
70 grad(fx_fvar.d_.vi_);
71 x_adj = x_var.adj();
72}
73
74} // namespace math
75} // namespace stan
76
77#endif
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...
Definition hessian.hpp:41
static void grad()
Compute the gradient for all variables starting from the end of the AD tape.
Definition grad.hpp:26
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.
Definition fvar.hpp:61
This template class represents scalars used in forward-mode automatic differentiation,...
Definition fvar.hpp:40