Stan  2.10.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/fun/Eigen.hpp>
5 #include <stan/math/prim/mat/meta/index_type.hpp>
8 #include <boost/random/variate_generator.hpp>
9 #include <boost/random/normal_distribution.hpp>
10 #include <Eigen/Cholesky>
11 
12 namespace stan {
13  namespace mcmc {
14 
15  // Euclidean manifold with dense metric
16  template <class Model, class BaseRNG>
18  : public base_hamiltonian<Model, dense_e_point, BaseRNG> {
19  public:
20  explicit dense_e_metric(const Model& model)
21  : base_hamiltonian<Model, dense_e_point, BaseRNG>(model) {}
22 
23  double T(dense_e_point& z) {
24  return 0.5 * z.p.transpose() * z.mInv * z.p;
25  }
26 
27  double tau(dense_e_point& z) {
28  return T(z);
29  }
30 
31  double phi(dense_e_point& z) {
32  return this->V(z);
33  }
34 
35  double dG_dt(dense_e_point& z,
38  return 2 * T(z) - z.q.dot(z.g);
39  }
40 
41  Eigen::VectorXd dtau_dq(
42  dense_e_point& z,
45  return Eigen::VectorXd::Zero(this->model_.num_params_r());
46  }
47 
48  Eigen::VectorXd dtau_dp(dense_e_point& z) {
49  return z.mInv * z.p;
50  }
51 
52  Eigen::VectorXd dphi_dq(
53  dense_e_point& z,
56  return z.g;
57  }
58 
59  void sample_p(dense_e_point& z, BaseRNG& rng) {
60  typedef typename stan::math::index_type<Eigen::VectorXd>::type idx_t;
61  boost::variate_generator<BaseRNG&, boost::normal_distribution<> >
62  rand_dense_gaus(rng, boost::normal_distribution<>());
63 
64  Eigen::VectorXd u(z.p.size());
65 
66  for (idx_t i = 0; i < u.size(); ++i)
67  u(i) = rand_dense_gaus();
68 
69  z.p = z.mInv.llt().matrixL().solve(u);
70  }
71  };
72 
73  } // mcmc
74 } // stan
75 #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.