Stan  2.14.0
probability, sampling & optimization
dense_e_metric.hpp
Go to the documentation of this file.
1 #ifndef STAN_MCMC_HMC_HAMILTONIANS_DENSE_E_METRIC_HPP
2 #define STAN_MCMC_HMC_HAMILTONIANS_DENSE_E_METRIC_HPP
3 
4 #include <stan/math/prim/mat.hpp>
7 #include <boost/random/variate_generator.hpp>
8 #include <boost/random/normal_distribution.hpp>
9 #include <Eigen/Cholesky>
10 
11 namespace stan {
12  namespace mcmc {
13 
14  // Euclidean manifold with dense metric
15  template <class Model, class BaseRNG>
17  : public base_hamiltonian<Model, dense_e_point, BaseRNG> {
18  public:
19  explicit dense_e_metric(const Model& model)
20  : base_hamiltonian<Model, dense_e_point, BaseRNG>(model) {}
21 
22  double T(dense_e_point& z) {
23  return 0.5 * z.p.transpose() * z.mInv * z.p;
24  }
25 
26  double tau(dense_e_point& z) {
27  return T(z);
28  }
29 
30  double phi(dense_e_point& z) {
31  return this->V(z);
32  }
33 
34  double dG_dt(dense_e_point& z,
37  return 2 * T(z) - z.q.dot(z.g);
38  }
39 
40  Eigen::VectorXd dtau_dq(
41  dense_e_point& z,
44  return Eigen::VectorXd::Zero(this->model_.num_params_r());
45  }
46 
47  Eigen::VectorXd dtau_dp(dense_e_point& z) {
48  return z.mInv * z.p;
49  }
50 
51  Eigen::VectorXd dphi_dq(
52  dense_e_point& z,
55  return z.g;
56  }
57 
58  void sample_p(dense_e_point& z, BaseRNG& rng) {
59  typedef typename stan::math::index_type<Eigen::VectorXd>::type idx_t;
60  boost::variate_generator<BaseRNG&, boost::normal_distribution<> >
61  rand_dense_gaus(rng, boost::normal_distribution<>());
62 
63  Eigen::VectorXd u(z.p.size());
64 
65  for (idx_t i = 0; i < u.size(); ++i)
66  u(i) = rand_dense_gaus();
67 
68  z.p = z.mInv.llt().matrixL().solve(u);
69  }
70  };
71 
72  } // mcmc
73 } // stan
74 #endif
Probability, optimization and sampling library.
double tau(dense_e_point &z)
Eigen::VectorXd dtau_dq(dense_e_point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
double phi(dense_e_point &z)
double dG_dt(dense_e_point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
Eigen::VectorXd q
Definition: ps_point.hpp:45
Eigen::VectorXd p
Definition: ps_point.hpp:46
Eigen::VectorXd dtau_dp(dense_e_point &z)
base_writer is an abstract base class defining the interface for Stan writer callbacks.
Definition: base_writer.hpp:20
Eigen::VectorXd dphi_dq(dense_e_point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
void sample_p(dense_e_point &z, BaseRNG &rng)
Eigen::VectorXd g
Definition: ps_point.hpp:49
double T(dense_e_point &z)
Point in a phase space with a base Euclidean manifold with dense metric.
dense_e_metric(const Model &model)

     [ Stan Home Page ] © 2011–2016, Stan Development Team.