Stan  2.14.0
probability, sampling & optimization
base_hamiltonian.hpp
Go to the documentation of this file.
1 #ifndef STAN_MCMC_HMC_HAMILTONIANS_BASE_HAMILTONIAN_HPP
2 #define STAN_MCMC_HMC_HAMILTONIANS_BASE_HAMILTONIAN_HPP
3 
7 #include <Eigen/Dense>
8 #include <iostream>
9 #include <limits>
10 #include <stdexcept>
11 #include <vector>
12 
13 namespace stan {
14  namespace mcmc {
15 
16  template <class Model, class Point, class BaseRNG>
18  public:
19  explicit base_hamiltonian(const Model& model)
20  : model_(model) {}
21 
23 
24  typedef Point PointType;
25 
26  virtual double T(Point& z) = 0;
27 
28  double V(Point& z) {
29  return z.V;
30  }
31 
32  virtual double tau(Point& z) = 0;
33 
34  virtual double phi(Point& z) = 0;
35 
36  double H(Point& z) {
37  return T(z) + V(z);
38  }
39 
40  // The time derivative of the virial, G = \sum_{d = 1}^{D} q^{d} p_{d}.
41  virtual double dG_dt(
42  Point& z,
45 
46  // tau = 0.5 p_{i} p_{j} Lambda^{ij} (q)
47  virtual Eigen::VectorXd dtau_dq(
48  Point& z,
51 
52  virtual Eigen::VectorXd dtau_dp(Point& z) = 0;
53 
54  // phi = 0.5 * log | Lambda (q) | + V(q)
55  virtual Eigen::VectorXd dphi_dq(
56  Point& z,
59 
60  virtual void sample_p(Point& z, BaseRNG& rng) = 0;
61 
62  void init(Point& z,
65  this->update_potential_gradient(z, info_writer, error_writer);
66  }
67 
69  Point& z,
72  try {
73  z.V = -stan::model::log_prob_propto<true>(model_, z.q);
74  } catch (const std::exception& e) {
75  this->write_error_msg_(e, error_writer);
76  z.V = std::numeric_limits<double>::infinity();
77  }
78  }
79 
81  Point& z,
84  try {
85  stan::model::gradient(model_, z.q, z.V, z.g, info_writer);
86  z.V = -z.V;
87  } catch (const std::exception& e) {
88  this->write_error_msg_(e, error_writer);
89  z.V = std::numeric_limits<double>::infinity();
90  }
91  z.g = -z.g;
92  }
93 
95  Point& z,
98 
100  Point& z,
103 
105  Point& z,
108  update_potential_gradient(z, info_writer, error_writer);
109  }
110 
111  protected:
112  const Model& model_;
113 
114  void write_error_msg_(const std::exception& e,
116  writer("Informational Message: The current Metropolis proposal "
117  "is about to be rejected because of the following issue:");
118  writer(e.what());
119  writer("If this warning occurs sporadically, such as for highly "
120  "constrained variable types like covariance matrices, then "
121  "the sampler is fine,");
122  writer("but if this warning occurs often then your model may be "
123  "either severely ill-conditioned or misspecified.");
124  writer();
125  }
126  };
127 
128  } // mcmc
129 } // stan
130 
131 #endif
virtual double T(Point &z)=0
Probability, optimization and sampling library.
virtual Eigen::VectorXd dtau_dp(Point &z)=0
void write_error_msg_(const std::exception &e, interface_callbacks::writer::base_writer &writer)
virtual double tau(Point &z)=0
void update_metric_gradient(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
void gradient(const M &model, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &x, double &f, Eigen::Matrix< double, Eigen::Dynamic, 1 > &grad_f, std::ostream *msgs=0)
Definition: gradient.hpp:14
virtual double phi(Point &z)=0
void update_gradients(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
void update_potential_gradient(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
virtual Eigen::VectorXd dphi_dq(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)=0
virtual void sample_p(Point &z, BaseRNG &rng)=0
base_writer is an abstract base class defining the interface for Stan writer callbacks.
Definition: base_writer.hpp:20
virtual double dG_dt(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)=0
void update_metric(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
void update_potential(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
base_hamiltonian(const Model &model)
void init(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)
virtual Eigen::VectorXd dtau_dq(Point &z, interface_callbacks::writer::base_writer &info_writer, interface_callbacks::writer::base_writer &error_writer)=0

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