Stan  2.10.0
probability, sampling & optimization
unit_e_metric.hpp
Go to the documentation of this file.
1 #ifndef STAN_MCMC_HMC_HAMILTONIANS_UNIT_E_METRIC_HPP
2 #define STAN_MCMC_HMC_HAMILTONIANS_UNIT_E_METRIC_HPP
3 
6 #include <boost/random/variate_generator.hpp>
7 #include <boost/random/normal_distribution.hpp>
8 
9 namespace stan {
10  namespace mcmc {
11 
12  // Euclidean manifold with unit metric
13  template <class Model, class BaseRNG>
15  : public base_hamiltonian<Model, unit_e_point, BaseRNG> {
16  public:
17  explicit unit_e_metric(const Model& model)
18  : base_hamiltonian<Model, unit_e_point, BaseRNG>(model) {}
19 
20  double T(unit_e_point& z) {
21  return 0.5 * z.p.squaredNorm();
22  }
23 
24  double tau(unit_e_point& z) {
25  return T(z);
26  }
27 
28  double phi(unit_e_point& z) {
29  return this->V(z);
30  }
31 
32  double dG_dt(unit_e_point& z,
35  return 2 * T(z) - z.q.dot(z.g);
36  }
37 
38  Eigen::VectorXd dtau_dq(
39  unit_e_point& z,
42  return Eigen::VectorXd::Zero(this->model_.num_params_r());
43  }
44 
45  Eigen::VectorXd dtau_dp(unit_e_point& z) {
46  return z.p;
47  }
48 
49  Eigen::VectorXd dphi_dq(
50  unit_e_point& z,
53  return z.g;
54  }
55 
56  void sample_p(unit_e_point& z, BaseRNG& rng) {
57  boost::variate_generator<BaseRNG&, boost::normal_distribution<> >
58  rand_unit_gaus(rng, boost::normal_distribution<>());
59 
60  for (int i = 0; i < z.p.size(); ++i)
61  z.p(i) = rand_unit_gaus();
62  }
63  };
64 
65  } // mcmc
66 } // stan
67 #endif
Probability, optimization and sampling library.
void sample_p(unit_e_point &z, BaseRNG &rng)
double tau(unit_e_point &z)
double phi(unit_e_point &z)
Eigen::VectorXd dtau_dq(unit_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
Point in a phase space with a base Euclidean manifold with unit metric.
Eigen::VectorXd p
Definition: ps_point.hpp:46
base_writer is an abstract base class defining the interface for Stan writer callbacks.
Definition: base_writer.hpp:20
double dG_dt(unit_e_point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
Eigen::VectorXd g
Definition: ps_point.hpp:49
double T(unit_e_point &z)
Eigen::VectorXd dphi_dq(unit_e_point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
unit_e_metric(const Model &model)
Eigen::VectorXd dtau_dp(unit_e_point &z)

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