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

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