Automatic Differentiation
 
Loading...
Searching...
No Matches
qr_thin_R.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_OPENCL_PRIM_QR_THIN_R_HPP
2#define STAN_MATH_OPENCL_PRIM_QR_THIN_R_HPP
3#ifdef STAN_OPENCL
7
8namespace stan {
9namespace math {
10
18template <typename T_m,
19 require_all_kernel_expressions_and_none_scalar_t<T_m>* = nullptr>
21 check_nonzero_size("qr_thin_R(OpenCL)", "m", m);
22
23 matrix_cl<double> mat_eval = std::forward<T_m>(m);
25
26 qr_decomposition_cl(m, Q, R);
27
28 auto R_block
29 = block_zero_based(R, 0, 0, std::min(m.rows(), m.cols()), m.cols());
30 R = select(rowwise_broadcast(diagonal(R) < 0.0), -R_block, R_block).eval();
31 return R;
32}
33} // namespace math
34} // namespace stan
35
36#endif
37#endif
Represents an arithmetic matrix on the OpenCL device.
Definition matrix_cl.hpp:47
select_< as_operation_cl_t< T_condition >, as_operation_cl_t< T_then >, as_operation_cl_t< T_else > > select(T_condition &&condition, T_then &&then, T_else &&els)
Selection operation on kernel generator expressions.
Definition select.hpp:148
auto block_zero_based(T &&a, int start_row, int start_col, int rows, int cols)
Block of a kernel generator expression.
auto rowwise_broadcast(T &&a)
Broadcast an expression in rowwise dimmension.
auto diagonal(T &&a)
Diagonal of a kernel generator expression.
Definition diagonal.hpp:136
void qr_decomposition_cl(const matrix_cl< double > &A, matrix_cl< double > &Q, matrix_cl< double > &R, int r=100)
Calculates QR decomposition of A using the block Householder algorithm.
matrix_cl< double > qr_thin_R(T_m &&m)
Returns the orthogonal factor of the thin QR decomposition.
Definition qr_thin_R.hpp:20
void check_nonzero_size(const char *function, const char *name, const T_y &y)
Check if the specified matrix/vector is of non-zero size.
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...