Automatic Differentiation
 
Loading...
Searching...
No Matches
read_cov_matrix.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_PRIM_FUN_READ_COV_MATRIX_HPP
2#define STAN_MATH_PRIM_FUN_READ_COV_MATRIX_HPP
3
8
9namespace stan {
10namespace math {
11
25template <typename T_CPCs, typename T_sds,
26 require_all_eigen_vector_t<T_CPCs, T_sds>* = nullptr,
27 require_vt_same<T_CPCs, T_sds>* = nullptr>
28Eigen::Matrix<value_type_t<T_CPCs>, Eigen::Dynamic, Eigen::Dynamic>
29read_cov_matrix(const T_CPCs& CPCs, const T_sds& sds,
30 value_type_t<T_CPCs>& log_prob) {
31 Eigen::Matrix<value_type_t<T_CPCs>, Eigen::Dynamic, Eigen::Dynamic> L
32 = read_cov_L(CPCs, sds, log_prob);
34}
35
46template <typename T_CPCs, typename T_sds,
49Eigen::Matrix<value_type_t<T_CPCs>, Eigen::Dynamic, Eigen::Dynamic>
50read_cov_matrix(const T_CPCs& CPCs, const T_sds& sds) {
51 size_t K = sds.rows();
52 Eigen::DiagonalMatrix<value_type_t<T_CPCs>, Eigen::Dynamic> D(K);
53 D.diagonal() = sds;
54 Eigen::Matrix<value_type_t<T_CPCs>, Eigen::Dynamic, Eigen::Dynamic> L
55 = D * read_corr_L(CPCs, K);
57}
58
59} // namespace math
60} // namespace stan
61
62#endif
require_all_t< is_eigen_vector< std::decay_t< Types > >... > require_all_eigen_vector_t
Require all of the types satisfy is_eigen_vector.
require_t< std::is_same< value_type_t< std::decay_t< T > >, value_type_t< std::decay_t< S > > > > require_vt_same
Value types of T and S satisfies std::is_same.
typename value_type< T >::type value_type_t
Helper function for accessing underlying type.
Eigen::Matrix< value_type_t< EigMat >, EigMat::RowsAtCompileTime, EigMat::RowsAtCompileTime > multiply_lower_tri_self_transpose(const EigMat &m)
Eigen::Matrix< value_type_t< T_CPCs >, Eigen::Dynamic, Eigen::Dynamic > read_cov_matrix(const T_CPCs &CPCs, const T_sds &sds, value_type_t< T_CPCs > &log_prob)
A generally worse alternative to call prior to evaluating the density of an elliptical distribution.
Eigen::Matrix< value_type_t< T_CPCs >, Eigen::Dynamic, Eigen::Dynamic > read_cov_L(const T_CPCs &CPCs, const T_sds &sds, value_type_t< T_CPCs > &log_prob)
This is the function that should be called prior to evaluating the density of any elliptical distribu...
Eigen::Matrix< value_type_t< T >, Eigen::Dynamic, Eigen::Dynamic > read_corr_L(const T &CPCs, size_t K)
Return the Cholesky factor of the correlation matrix of the specified dimensionality corresponding to...
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...