 |
fml
0.1-0
Fused Matrix Library
|
5 #ifndef FML_CPU_LINALG_LINALG_NORM_H
6 #define FML_CPU_LINALG_LINALG_NORM_H
13 #include "../../_internals/linalgutils.hh"
14 #include "../../_internals/omp.hh"
16 #include "../cpumat.hh"
17 #include "../cpuvec.hh"
20 #include "linalg_qr.hh"
21 #include "linalg_svd.hh"
37 template <
typename REAL>
40 const len_t m = x.
nrows();
41 const len_t n = x.
ncols();
46 #pragma omp parallel for if(m*n > fml::omp::OMP_MIN_SIZE) reduction(max:norm)
47 for (len_t j=0; j<n; j++)
51 #pragma omp simd reduction(+:tmp)
52 for (len_t i=0; i<m; i++)
53 tmp += fabs(x_d[i + m*j]);
78 template <
typename REAL>
81 const len_t m = x.
nrows();
82 const len_t n = x.
ncols();
91 for (len_t j=0; j<n; j++)
93 for (len_t i=0; i<m; i++)
94 mars_d[i] += fabs(x_d[i + m*j]);
97 for (len_t i=0; i<m; i++)
117 template <
typename REAL>
120 const len_t m = x.
nrows();
121 const len_t n = x.
ncols();
127 for (len_t j=0; j<n; j++)
128 lapack::lassq(m, x_d + m*j, 1, &scale, &sumsq);
130 return scale * sqrtf(sumsq);
144 template <
typename REAL>
147 const len_t m = x.
nrows();
148 const len_t n = x.
ncols();
153 #pragma omp parallel for if(m*n > fml::omp::OMP_MIN_SIZE) reduction(max:norm)
154 for (len_t j=0; j<n; j++)
156 for (len_t i=0; i<m; i++)
158 REAL tmp = fabs(x_d[i + m*j]);
187 template <
typename REAL>
202 template <
typename REAL>
203 REAL cond_square_internals(
const char norm,
cpumat<REAL> &x)
205 const len_t n = x.
nrows();
215 lapack::gecon(norm, n, x.
data_ptr(), n, xnorm, &ret, work.data_ptr(),
216 work2.data_ptr(), &info);
218 fml::linalgutils::check_info(info,
"gecon");
220 return ((REAL)1)/ret;
223 template <
typename REAL>
224 REAL cond_nonsquare_internals(
const char norm, cpumat<REAL> &x)
226 const len_t m = x.nrows();
227 const len_t n = x.ncols();
240 aux.resize(R.nrows());
241 lapack::trcon(norm,
'U',
'N', n, R.data_ptr(), n, &ret,
242 x.data_ptr(), aux.data_ptr(), &info);
250 aux.resize(L.nrows());
251 lapack::trcon(norm,
'L',
'N', m, L.data_ptr(), m, &ret,
252 x.data_ptr(), aux.data_ptr(), &info);
255 fml::linalgutils::check_info(info,
"trcon");
257 return ((REAL)1)/ret;
280 template <
typename REAL>
284 return cond_square_internals(
'1', x);
286 return cond_nonsquare_internals(
'1', x);
310 template <
typename REAL>
314 return cond_square_internals(
'I', x);
316 return cond_nonsquare_internals(
'I', x);
338 template <
typename REAL>
349 #pragma omp parallel for if(s.size() > fml::omp::OMP_MIN_SIZE) reduction(max:max) reduction(min:min)
350 for (len_t i=1; i<s.
size(); i++)
354 if (s_d[i] > 0 && s_d[i] < min)
Matrix class for data held on a single CPU.
Definition: cpumat.hh:36
REAL norm_I(const cpumat< REAL > &x)
Computes the infinity matrix norm, the maximum absolute row sum.
Definition: linalg_norm.hh:79
REAL cond_2(cpumat< REAL > &x)
Estimates the condition number under the 2 norm.
Definition: linalg_norm.hh:339
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 cond_1(cpumat< REAL > &x)
Estimates the condition number under the 1-norm.
Definition: linalg_norm.hh:281
void qr(const bool pivot, cpumat< REAL > &x, cpuvec< REAL > &qraux)
Computes the QR decomposition.
Definition: linalg_qr.hh:94
REAL norm_1(const cpumat< REAL > &x)
Computes the 1 matrix norm, the maximum absolute column sum.
Definition: linalg_norm.hh:38
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
T get(const len_t i) const
Get the specified value.
Definition: cpuvec.hh:514
REAL cond_I(cpumat< REAL > &x)
Estimates the condition number under the infinity norm.
Definition: linalg_norm.hh:311
void svd(cpumat< REAL > &x, cpuvec< REAL > &s)
Computes the singular value decomposition.
Definition: linalg_svd.hh:97
REAL norm_2(cpumat< REAL > &x)
Computes the 2/spectral matrix norm.
Definition: linalg_norm.hh:188
REAL * data_ptr()
Pointer to the internal array.
Definition: unimat.hh:40
void lq_L(const cpumat< REAL > &LQ, cpumat< REAL > &L)
Recover the L matrix from an LQ decomposition.
Definition: linalg_qr.hh:247
Core namespace.
Definition: dimops.hh:10
len_t size() const
Number of elements in the vector.
Definition: univec.hh:26
void qr_R(const cpumat< REAL > &QR, cpumat< REAL > &R)
Recover the R matrix from a QR decomposition.
Definition: linalg_qr.hh:162
void lq(cpumat< REAL > &x, cpuvec< REAL > &lqaux)
Computes the LQ decomposition.
Definition: linalg_qr.hh:223
void fill_zero()
Set all values to zero.
Definition: cpuvec.hh:316
REAL norm_M(const cpumat< REAL > &x)
Computes the maximum modulus matrix norm.
Definition: linalg_norm.hh:145
REAL norm_F(const cpumat< REAL > &x)
Computes the Frobenius/Euclidean matrix norm.
Definition: linalg_norm.hh:118