fml  0.1-0
Fused Matrix Library
qr.hh
1 // This file is part of fml which is released under the Boost Software
2 // License, Version 1.0. See accompanying file LICENSE or copy at
3 // https://www.boost.org/LICENSE_1_0.txt
4 
5 #ifndef FML_PAR_CPU_LINALG_QR_H
6 #define FML_PAR_CPU_LINALG_QR_H
7 #pragma once
8 
9 
10 #include "../parmat.hh"
11 
12 #include "../../../cpu/linalg/linalg_blas.hh"
13 #include "../../../cpu/linalg/linalg_invert.hh"
14 #include "../../../cpu/linalg/linalg_qr.hh"
15 
16 #include "qr_allreduce.hh"
17 
18 
19 namespace fml
20 {
21 namespace linalg
22 {
23  namespace internals
24  {
25  template <typename REAL>
26  void qr_R(const int root, parmat_cpu<REAL> &x, cpumat<REAL> &R,
27  cpumat<REAL> &R_local, cpuvec<REAL> &qraux)
28  {
29  const len_t n = x.ncols();
30 
31  linalg::qr(false, x.data_obj(), qraux);
32  linalg::qr_R(x.data_obj(), R_local);
33 
34  R.resize(n, n);
35  tsqr::qr_allreduce(root, n, n, R_local.data_ptr(), R.data_ptr(),
36  x.get_comm().get_comm());
37  }
38  }
39 
40  template <typename REAL>
41  void qr_R(const int root, parmat_cpu<REAL> &x, cpumat<REAL> &R)
42  {
43  if (x.nrows() < (len_global_t) x.ncols())
44  throw std::runtime_error("impossible dimensions");
45 
46  cpumat<REAL> R_local;
47  cpuvec<REAL> qraux;
48 
49  internals::qr_R(root, x, R, R_local, qraux);
50  }
51 
52 
53 
54  namespace internals
55  {
56  template <typename REAL>
57  void qr_Q(const parmat_cpu<REAL> &x, parmat_cpu<REAL> &x_cpy,
58  cpumat<REAL> &R, cpumat<REAL> &R_local,
59  cpuvec<REAL> &qraux, parmat_cpu<REAL> &Q)
60  {
61  copy::cpu2cpu(x.data_obj(), x_cpy.data_obj());
62  qr_R(mpi::REDUCE_TO_ALL, x_cpy, R, R_local, qraux);
63  linalg::trinv(true, false, R);
64  matmult(x, R, Q);
65  }
66  }
67 
68  template <typename REAL>
69  void qr_Q(parmat_cpu<REAL> &x, cpuvec<REAL> &qraux, parmat_cpu<REAL> &Q)
70  {
71  cpumat<REAL> R, R_local;
72 
73  qr_R(mpi::REDUCE_TO_ALL, x, R, R_local, qraux);
74  linalg::trinv(true, false, R);
75  matmult(x, R, Q);
76  }
77 }
78 }
79 
80 
81 #endif
fml::copy::cpu2cpu
void cpu2cpu(const cpuvec< REAL_IN > &cpu_in, cpuvec< REAL_OUT > &cpu_out)
Copy data from a CPU object to another.
Definition: copy.hh:40
fml::linalg::qr
void qr(const bool pivot, cpumat< REAL > &x, cpuvec< REAL > &qraux)
Computes the QR decomposition.
Definition: qr.hh:94
fml
Core namespace.
Definition: dimops.hh:10
fml::linalg::qr_Q
void qr_Q(const cpumat< REAL > &QR, const cpuvec< REAL > &qraux, cpumat< REAL > &Q, cpuvec< REAL > &work)
Recover the Q matrix from a QR decomposition.
Definition: qr.hh:120
fml::linalg::trinv
void trinv(const bool upper, const bool unit_diag, cpumat< REAL > &x)
Compute the matrix inverse of a triangular matrix.
Definition: invert.hh:87
fml::linalg::qr_R
void qr_R(const cpumat< REAL > &QR, cpumat< REAL > &R)
Recover the R matrix from a QR decomposition.
Definition: qr.hh:162
fml::linalg::matmult
void matmult(const bool transx, const bool transy, const REAL alpha, const cpumat< REAL > &x, const cpumat< REAL > &y, cpumat< REAL > &ret)
Computes ret = alpha*op(x)*op(y) where op(A) is A or A^T.
Definition: matmult.hh:43