Automatic Differentiation
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
cholesky_factor_constrain.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_PRIM_CONSTRAINT_CHOLESKY_FACTOR_CONSTRAIN_HPP
2#define STAN_MATH_PRIM_CONSTRAINT_CHOLESKY_FACTOR_CONSTRAIN_HPP
3
9#include <cmath>
10#include <stdexcept>
11#include <vector>
12
13namespace stan {
14namespace math {
15
28template <typename T, require_eigen_col_vector_t<T>* = nullptr>
29inline Eigen::Matrix<value_type_t<T>, Eigen::Dynamic, Eigen::Dynamic>
30cholesky_factor_constrain(const T& x, int M, int N) {
31 using std::exp;
32 using T_scalar = value_type_t<T>;
33 check_greater_or_equal("cholesky_factor_constrain",
34 "num rows (must be greater or equal to num cols)", M,
35 N);
36 check_size_match("cholesky_factor_constrain", "constrain size", x.size(),
37 "((N * (N + 1)) / 2 + (M - N) * N)",
38 ((N * (N + 1)) / 2 + (M - N) * N));
39 Eigen::Matrix<T_scalar, Eigen::Dynamic, Eigen::Dynamic> y(M, N);
40 int pos = 0;
41
42 const auto& x_ref = to_ref(x);
43 for (int m = 0; m < N; ++m) {
44 y.row(m).head(m) = x_ref.segment(pos, m);
45 pos += m;
46 y.coeffRef(m, m) = exp(x_ref.coeff(pos++));
47 y.row(m).tail(N - m - 1).setZero();
48 }
49
50 for (int m = N; m < M; ++m) {
51 y.row(m) = x_ref.segment(pos, N);
52 pos += N;
53 }
54 return y;
55}
56
73template <typename T, typename Lp, require_eigen_vector_t<T>* = nullptr,
74 require_convertible_t<return_type_t<T>, Lp>* = nullptr>
75inline Eigen::Matrix<value_type_t<T>, Eigen::Dynamic, Eigen::Dynamic>
76cholesky_factor_constrain(const T& x, int M, int N, Lp& lp) {
77 check_size_match("cholesky_factor_constrain", "x.size()", x.size(),
78 "((N * (N + 1)) / 2 + (M - N) * N)",
79 ((N * (N + 1)) / 2 + (M - N) * N));
80 int pos = 0;
81 const auto& x_ref = to_ref(x);
82 for (int n = 0; n < N; ++n) {
83 pos += n;
84 lp += x_ref.coeff(pos++);
85 }
86 return cholesky_factor_constrain(x_ref, M, N);
87}
88
103template <typename T, require_std_vector_t<T>* = nullptr>
104inline auto cholesky_factor_constrain(const T& x, int M, int N) {
106 x, [M, N](auto&& v) { return cholesky_factor_constrain(v, M, N); });
107}
108
126template <typename T, typename Lp, require_std_vector_t<T>* = nullptr,
127 require_convertible_t<return_type_t<T>, Lp>* = nullptr>
128inline auto cholesky_factor_constrain(const T& x, int M, int N, Lp& lp) {
129 return apply_vector_unary<T>::apply(x, [&lp, M, N](auto&& v) {
130 return cholesky_factor_constrain(v, M, N, lp);
131 });
132}
133
156template <bool Jacobian, typename T, typename Lp,
158inline auto cholesky_factor_constrain(const T& x, int M, int N, Lp& lp) {
159 if constexpr (Jacobian) {
160 return cholesky_factor_constrain(x, M, N, lp);
161 } else {
162 return cholesky_factor_constrain(x, M, N);
163 }
164}
165
166} // namespace math
167} // namespace stan
168
169#endif
require_t< std::is_convertible< std::decay_t< T >, std::decay_t< S > > > require_convertible_t
Require types T and S satisfies std::is_convertible.
typename value_type< T >::type value_type_t
Helper function for accessing underlying type.
void check_greater_or_equal(const char *function, const char *name, const T_y &y, const T_low &low, Idxs... idxs)
Throw an exception if y is not greater or equal than low.
Eigen::Matrix< value_type_t< T >, Eigen::Dynamic, Eigen::Dynamic > cholesky_factor_constrain(const T &x, int M, int N)
Return the Cholesky factor of the specified size read from the specified vector.
ref_type_t< T && > to_ref(T &&a)
This evaluates expensive Eigen expressions.
Definition to_ref.hpp:17
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.
fvar< T > exp(const fvar< T > &x)
Definition exp.hpp:15
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...