 |
fml
0.1-0
Fused Matrix Library
|
5 #ifndef FML_CPU_LINALG_LINALG_MISC_H
6 #define FML_CPU_LINALG_LINALG_MISC_H
13 #include "../../_internals/linalgutils.hh"
14 #include "../../_internals/omp.hh"
16 #include "../cpumat.hh"
17 #include "../cpuvec.hh"
20 #include "linalg_lu.hh"
45 template <
typename REAL>
48 const len_t m = x.
nrows();
50 throw std::runtime_error(
"'x' must be a square matrix");
74 for (
int i=0; i<m; i++)
76 if (ipiv[i] != (i + 1))
82 #pragma omp parallel for reduction(+:mod) reduction(*:sgn)
83 for (
int i=0; i<m; i++)
85 const REAL d = a[i + m*i];
108 template <
typename REAL>
112 const len_t m = x.
nrows();
113 const len_t minmn = std::min(m, x.
ncols());
116 #pragma omp parallel for simd if(minmn > fml::omp::OMP_MIN_SIZE) reduction(+:tr)
117 for (len_t i=0; i<minmn; i++)
Matrix class for data held on a single CPU.
Definition: cpumat.hh:36
bool is_square() const
Is the matrix square?
Definition: unimat.hh:34
T * data_ptr()
Pointer to the internal array.
Definition: univec.hh:28
len_t nrows() const
Number of rows.
Definition: unimat.hh:36
REAL trace(const cpumat< REAL > &x)
Computes the trace, i.e. the sum of the diagonal.
Definition: linalg_misc.hh:109
void lu(cpumat< REAL > &x, cpuvec< int > &p, int &info)
Computes the PLU factorization with partial pivoting.
Definition: linalg_lu.hh:48
Vector class for data held on a single CPU.
Definition: cpuvec.hh:31
len_t ncols() const
Number of columns.
Definition: unimat.hh:38
void det(cpumat< REAL > &x, int &sign, REAL &modulus)
Computes the determinant in logarithmic form.
Definition: linalg_misc.hh:46
REAL * data_ptr()
Pointer to the internal array.
Definition: unimat.hh:40
Core namespace.
Definition: dimops.hh:10