gradient.hpp

Go to the documentation of this file.

#ifndef STAN_MATH_REV_FUNCTOR_GRADIENT_HPP
#define STAN_MATH_REV_FUNCTOR_GRADIENT_HPP
 
#include <stan/math/rev/meta.hpp>
#include <stan/math/rev/core.hpp>
#include <stan/math/prim/fun/Eigen.hpp>
#include <stdexcept>
#include <vector>
 
namespace stan {
namespace math {
 
template <typename F>
inline void gradient(const F& f,
                     const Eigen::Matrix<double, Eigen::Dynamic, 1>& x,
                     double& fx,
                     Eigen::Matrix<double, Eigen::Dynamic, 1>& grad_fx) {
  nested_rev_autodiff nested;
 
  Eigen::Matrix<var, Eigen::Dynamic, 1> x_var(x);
  var fx_var = f(x_var);
  fx = fx_var.val();
  grad_fx.resize(x.size());
  grad(fx_var.vi_);
  grad_fx = x_var.adj();
}
 
template <typename F, typename EigVec, typename InputIt,
          require_eigen_vector_vt<std::is_arithmetic, EigVec>* = nullptr>
inline void gradient(const F& f, const EigVec& x, double& fx,
                     InputIt first_grad_fx, InputIt last_grad_fx) {
  nested_rev_autodiff nested;
 
  if (last_grad_fx - first_grad_fx != x.size()) {
    std::stringstream s;
    s << "gradient(): iterator and gradient different sizes; iterator size = "
      << last_grad_fx - first_grad_fx << "; grad size = " << x.size()
      << std::endl;
    throw std::invalid_argument(s.str());
  }
 
  Eigen::Matrix<var, Eigen::Dynamic, 1> x_var(x);
  var fx_var = f(x_var);
  fx = fx_var.val();
  grad(fx_var.vi_);
  for (Eigen::VectorXd::Index i = 0; i < x_var.size(); ++i) {
    *first_grad_fx++ = x_var.coeff(i).adj();
  }
}
 
}  // namespace math
}  // namespace stan
#endif