Stan  2.10.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 
5 #include <stan/math/prim/mat/fun/Eigen.hpp>
6 #include <stan/model/util.hpp>
7 #include <iostream>
8 #include <limits>
9 #include <stdexcept>
10 #include <vector>
11 
12 namespace stan {
13  namespace mcmc {
14 
15  template <class Model, class Point, class BaseRNG>
17  public:
18  explicit base_hamiltonian(const Model& model)
19  : model_(model) {}
20 
22 
23  typedef Point PointType;
24 
25  virtual double T(Point& z) = 0;
26 
27  double V(Point& z) {
28  return z.V;
29  }
30 
31  virtual double tau(Point& z) = 0;
32 
33  virtual double phi(Point& z) = 0;
34 
35  double H(Point& z) {
36  return T(z) + V(z);
37  }
38 
39  // The time derivative of the virial, G = \sum_{d = 1}^{D} q^{d} p_{d}.
40  virtual double dG_dt(
41  Point& z,
44 
45  // tau = 0.5 p_{i} p_{j} Lambda^{ij} (q)
46  virtual Eigen::VectorXd dtau_dq(
47  Point& z,
50 
51  virtual Eigen::VectorXd dtau_dp(Point& z) = 0;
52 
53  // phi = 0.5 * log | Lambda (q) | + V(q)
54  virtual Eigen::VectorXd dphi_dq(
55  Point& z,
58 
59  virtual void sample_p(Point& z, BaseRNG& rng) = 0;
60 
61  void init(Point& z,
64  this->update_potential_gradient(z, info_writer, error_writer);
65  }
66 
68  Point& z,
71  try {
72  z.V = -stan::model::log_prob_propto<true>(model_, z.q);
73  } catch (const std::exception& e) {
74  this->write_error_msg_(e, error_writer);
75  z.V = std::numeric_limits<double>::infinity();
76  }
77  }
78 
80  Point& z,
83  try {
84  stan::model::gradient(model_, z.q, z.V, z.g, info_writer);
85  z.V = -z.V;
86  } catch (const std::exception& e) {
87  this->write_error_msg_(e, error_writer);
88  z.V = std::numeric_limits<double>::infinity();
89  }
90  z.g = -z.g;
91  }
92 
94  Point& z,
97 
99  Point& z,
102 
104  Point& z,
107  update_potential_gradient(z, info_writer, error_writer);
108  }
109 
110  protected:
111  const Model& model_;
112 
113  void write_error_msg_(const std::exception& e,
115  writer();
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();
123  writer("but if this warning occurs often then your model may be "
124  "either severely ill-conditioned or misspecified.");
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: util.hpp:422
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.