Stan Math Library  2.20.0
reverse mode automatic differentiation
multi_normal_prec_rng.hpp
Go to the documentation of this file.
1 #ifndef STAN_MATH_PRIM_MAT_PROB_MULTI_NORMAL_PREC_RNG_HPP
2 #define STAN_MATH_PRIM_MAT_PROB_MULTI_NORMAL_PREC_RNG_HPP
3 
10 #include <boost/random/normal_distribution.hpp>
11 #include <boost/random/variate_generator.hpp>
12 
13 namespace stan {
14 namespace math {
15 
32 template <typename T_loc, class RNG>
34 multi_normal_prec_rng(const T_loc &mu, const Eigen::MatrixXd &S, RNG &rng) {
35  using boost::normal_distribution;
36  using boost::variate_generator;
37 
38  static const char *function = "multi_normal_prec_rng";
39 
40  check_positive(function, "Precision matrix rows", S.rows());
41  check_finite(function, "Precision matrix", S);
42  check_symmetric(function, "Precision matrix", S);
43 
44  Eigen::LLT<Eigen::MatrixXd> llt_of_S = S.llt();
45  check_pos_definite(function, "precision matrix argument", llt_of_S);
46 
47  vector_seq_view<T_loc> mu_vec(mu);
48  check_positive(function, "number of location parameter vectors",
49  mu_vec.size());
50  size_t size_mu = mu_vec[0].size();
51 
52  size_t N = mu_vec.size();
53 
54  for (size_t i = 1; i < N; i++) {
55  int size_mu_new = mu_vec[i].size();
56  check_size_match(function,
57  "Size of one of the vectors of "
58  "the location variable",
59  size_mu_new,
60  "Size of another vector of the "
61  "location variable",
62  size_mu);
63  }
64 
65  for (size_t i = 0; i < N; i++) {
66  check_finite(function, "Location parameter", mu_vec[i]);
67  }
68 
69  check_size_match(function, "Rows of location parameter", size_mu, "Rows of S",
70  S.rows());
71 
73 
74  variate_generator<RNG &, normal_distribution<>> std_normal_rng(
75  rng, normal_distribution<>(0, 1));
76 
77  for (size_t n = 0; n < N; ++n) {
78  Eigen::VectorXd z(S.cols());
79  for (int i = 0; i < S.cols(); i++)
80  z(i) = std_normal_rng();
81 
82  output[n] = Eigen::VectorXd(mu_vec[n]) + llt_of_S.matrixU().solve(z);
83  }
84 
85  return output.data();
86 }
87 
88 } // namespace math
89 } // namespace stan
90 #endif
void check_finite(const char *function, const char *name, const T_y &y)
Check if y is finite.
void check_size_match(const char *function, const char *name_i, T_size1 i, const char *name_j, T_size2 j)
Check if the provided sizes match.
StdVectorBuilder allocates type T1 values to be used as intermediate values.
StdVectorBuilder< true, Eigen::VectorXd, T_loc >::type multi_normal_prec_rng(const T_loc &mu, const Eigen::MatrixXd &S, RNG &rng)
Return a multivariate normal random variate with the given location and precision using the specified...
This class provides a low-cost wrapper for situations where you either need an Eigen Vector or RowVec...
void check_positive(const char *function, const char *name, const T_y &y)
Check if y is positive.
void check_pos_definite(const char *function, const char *name, const Eigen::Matrix< T_y, -1, -1 > &y)
Check if the specified square, symmetric matrix is positive definite.
void check_symmetric(const char *function, const char *name, const matrix_cl &y)
Check if the matrix_cl is symmetric.

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