Stan  2.10.0
probability, sampling & optimization
ps_point.hpp
Go to the documentation of this file.
1 #ifndef STAN_MCMC_HMC_HAMILTONIANS_PS_POINT_HPP
2 #define STAN_MCMC_HMC_HAMILTONIANS_PS_POINT_HPP
3 
5 #include <stan/math/prim/mat/fun/Eigen.hpp>
6 #include <boost/lexical_cast.hpp>
7 #include <string>
8 #include <vector>
9 
10 namespace stan {
11  namespace mcmc {
12  using Eigen::Dynamic;
13 
17  class ps_point {
18  friend class ps_point_test;
19 
20  public:
21  explicit ps_point(int n)
22  : q(n), p(n), V(0), g(n) {}
23 
24  ps_point(const ps_point& z)
25  : q(z.q.size()), p(z.p.size()), V(z.V), g(z.g.size()) {
26  fast_vector_copy_<double>(q, z.q);
27  fast_vector_copy_<double>(p, z.p);
28  fast_vector_copy_<double>(g, z.g);
29  }
30 
32  if (this == &z)
33  return *this;
34 
35  fast_vector_copy_<double>(q, z.q);
36 
37  V = z.V;
38 
39  fast_vector_copy_<double>(p, z.p);
40  fast_vector_copy_<double>(g, z.g);
41 
42  return *this;
43  }
44 
45  Eigen::VectorXd q;
46  Eigen::VectorXd p;
47 
48  double V;
49  Eigen::VectorXd g;
50 
51  virtual void get_param_names(std::vector<std::string>& model_names,
52  std::vector<std::string>& names) {
53  for (int i = 0; i < q.size(); ++i)
54  names.push_back(model_names.at(i));
55  for (int i = 0; i < q.size(); ++i)
56  names.push_back(std::string("p_") + model_names.at(i));
57  for (int i = 0; i < q.size(); ++i)
58  names.push_back(std::string("g_") + model_names.at(i));
59  }
60 
61  virtual void get_params(std::vector<double>& values) {
62  for (int i = 0; i < q.size(); ++i)
63  values.push_back(q(i));
64  for (int i = 0; i < q.size(); ++i)
65  values.push_back(p(i));
66  for (int i = 0; i < q.size(); ++i)
67  values.push_back(g(i));
68  }
69 
75  virtual void
77 
78  protected:
79  template <typename T>
80  static inline void
81  fast_vector_copy_(Eigen::Matrix<T, Dynamic, 1>& v_to,
82  const Eigen::Matrix<T, Dynamic, 1>& v_from) {
83  int sz = v_from.size();
84  v_to.resize(sz);
85  if (sz > 0)
86  std::memcpy(&v_to(0), &v_from(0), v_from.size() * sizeof(double));
87  }
88 
89  template <typename T>
90  static inline void
91  fast_matrix_copy_(Eigen::Matrix<T, Dynamic, Dynamic>& v_to,
92  const Eigen::Matrix<T, Dynamic, Dynamic>& v_from) {
93  int nr = v_from.rows();
94  int nc = v_from.cols();
95  v_to.resize(nr, nc);
96  if (nr > 0 && nc > 0)
97  std::memcpy(&v_to(0), &v_from(0), v_from.size() * sizeof(double));
98  }
99  };
100 
101  } // mcmc
102 } // stan
103 #endif
virtual void get_params(std::vector< double > &values)
Definition: ps_point.hpp:61
ps_point & operator=(const ps_point &z)
Definition: ps_point.hpp:31
Probability, optimization and sampling library.
Point in a generic phase space.
Definition: ps_point.hpp:17
virtual void get_param_names(std::vector< std::string > &model_names, std::vector< std::string > &names)
Definition: ps_point.hpp:51
Eigen::VectorXd q
Definition: ps_point.hpp:45
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
virtual void write_metric(stan::interface_callbacks::writer::base_writer &writer)
Writes the metric.
Definition: ps_point.hpp:76
Eigen::VectorXd g
Definition: ps_point.hpp:49
ps_point(const ps_point &z)
Definition: ps_point.hpp:24
friend class ps_point_test
Definition: ps_point.hpp:18
static void fast_matrix_copy_(Eigen::Matrix< T, Dynamic, Dynamic > &v_to, const Eigen::Matrix< T, Dynamic, Dynamic > &v_from)
Definition: ps_point.hpp:91
static void fast_vector_copy_(Eigen::Matrix< T, Dynamic, 1 > &v_to, const Eigen::Matrix< T, Dynamic, 1 > &v_from)
Definition: ps_point.hpp:81

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