simplifications

This commit is contained in:
conrad
2021-03-25 15:40:28 +10:00
parent be701c8e1f
commit 97eefbce85
3 changed files with 40 additions and 52 deletions
+16 -21
View File
@@ -128,7 +128,7 @@ arma_warn_unused
typename enable_if2< is_arma_sparse_type<T1>::value, typename T1::pod_type >::result
norm
(
const T1& X,
const T1& expr,
const uword k = uword(2),
const typename arma_real_or_cx_only<typename T1::elem_type>::result* junk = nullptr
)
@@ -139,19 +139,17 @@ norm
typedef typename T1::elem_type eT;
typedef typename T1::pod_type T;
const SpProxy<T1> P(X);
const unwrap_spmat<T1> U(expr);
const SpMat<eT>& X = U.M;
if(P.get_n_nonzero() == 0) { return T(0); }
if(X.n_nonzero == 0) { return T(0); }
const bool is_vec = (T1::is_xvec) || (T1::is_row) || (T1::is_col) || (P.get_n_rows() == 1) || (P.get_n_cols() == 1);
const bool is_vec = (T1::is_xvec) || (T1::is_row) || (T1::is_col) || (X.n_rows == 1) || (X.n_cols == 1);
if(is_vec)
{
const unwrap_spmat<typename SpProxy<T1>::stored_type> tmp(P.Q);
const SpMat<eT>& A = tmp.M;
// create a fake dense vector to allow reuse of code for dense vectors
Col<eT> fake_vector( access::rwp(A.values), A.n_nonzero, false );
Col<eT> fake_vector( access::rwp(X.values), X.n_nonzero, false );
const Proxy< Col<eT> > P_fake_vector(fake_vector);
@@ -164,8 +162,8 @@ norm
}
else
{
if(k == uword(1)) { return spop_norm::mat_norm_1(P); }
if(k == uword(2)) { return spop_norm::mat_norm_2(P); }
if(k == uword(1)) { return spop_norm::mat_norm_1(X); }
if(k == uword(2)) { return spop_norm::mat_norm_2(X); }
arma_stop_logic_error("norm(): unsupported or unimplemented norm type for sparse matrices");
}
@@ -181,7 +179,7 @@ arma_warn_unused
typename enable_if2< is_arma_sparse_type<T1>::value, typename T1::pod_type >::result
norm
(
const T1& X,
const T1& expr,
const char* method,
const typename arma_real_or_cx_only<typename T1::elem_type>::result* junk = nullptr
)
@@ -192,22 +190,19 @@ norm
typedef typename T1::elem_type eT;
typedef typename T1::pod_type T;
const SpProxy<T1> P(X);
const unwrap_spmat<T1> U(expr);
const SpMat<eT>& X = U.M;
if(P.get_n_nonzero() == 0) { return T(0); }
const unwrap_spmat<typename SpProxy<T1>::stored_type> tmp(P.Q);
const SpMat<eT>& A = tmp.M;
if(X.n_nonzero == 0) { return T(0); }
// create a fake dense vector to allow reuse of code for dense vectors
Col<eT> fake_vector( access::rwp(A.values), A.n_nonzero, false );
Col<eT> fake_vector( access::rwp(X.values), X.n_nonzero, false );
const Proxy< Col<eT> > P_fake_vector(fake_vector);
const char sig = (method != nullptr) ? method[0] : char(0);
const bool is_vec = (T1::is_xvec) || (T1::is_row) || (T1::is_col) || (P.get_n_rows() == 1) || (P.get_n_cols() == 1);
const bool is_vec = (T1::is_xvec) || (T1::is_row) || (T1::is_col) || (X.n_rows == 1) || (X.n_cols == 1);
if(is_vec)
{
@@ -220,7 +215,7 @@ norm
{
const T val = op_norm::vec_norm_min(P_fake_vector);
return (P.get_n_nonzero() < P.get_n_elem()) ? T((std::min)(T(0), val)) : T(val);
return (X.n_nonzero < X.n_elem) ? T((std::min)(T(0), val)) : T(val);
}
else
if( (sig == 'f') || (sig == 'F') )
@@ -234,7 +229,7 @@ norm
{
if( (sig == 'i') || (sig == 'I') || (sig == '+') ) // inf norm
{
return spop_norm::mat_norm_inf(P);
return spop_norm::mat_norm_inf(X);
}
else
if( (sig == 'f') || (sig == 'F') )
+4 -5
View File
@@ -23,14 +23,13 @@ class spop_norm
{
public:
template<typename T1> inline static typename T1::pod_type mat_norm_1(const SpProxy<T1>& P);
template<typename eT> inline static typename get_pod_type<eT>::result mat_norm_1(const SpMat<eT>& X);
template<typename T1> inline static typename T1::pod_type mat_norm_2(const SpProxy<T1>& P, const typename arma_real_only<typename T1::elem_type>::result* junk = nullptr);
template<typename T1> inline static typename T1::pod_type mat_norm_2(const SpProxy<T1>& P, const typename arma_cx_only<typename T1::elem_type>::result* junk = nullptr);
template<typename eT> inline static typename get_pod_type<eT>::result mat_norm_2(const SpMat<eT>& X, const typename arma_real_only<eT>::result* junk = nullptr);
template<typename eT> inline static typename get_pod_type<eT>::result mat_norm_2(const SpMat<eT>& X, const typename arma_cx_only<eT>::result* junk = nullptr);
template<typename T1> inline static typename T1::pod_type mat_norm_inf(const SpProxy<T1>& P);
template<typename eT> inline static typename get_pod_type<eT>::result mat_norm_inf(const SpMat<eT>& X);
};
//! @}
+20 -26
View File
@@ -19,23 +19,23 @@
template<typename T1>
template<typename eT>
inline
typename T1::pod_type
spop_norm::mat_norm_1(const SpProxy<T1>& P)
typename get_pod_type<eT>::result
spop_norm::mat_norm_1(const SpMat<eT>& X)
{
arma_extra_debug_sigprint();
// TODO: this can be sped up with a dedicated implementation
return as_scalar( max( sum(abs(P.Q), 0), 1) );
return as_scalar( max( sum(abs(X), 0), 1) );
}
template<typename T1>
template<typename eT>
inline
typename T1::pod_type
spop_norm::mat_norm_2(const SpProxy<T1>& P, const typename arma_real_only<typename T1::elem_type>::result* junk)
typename get_pod_type<eT>::result
spop_norm::mat_norm_2(const SpMat<eT>& X, const typename arma_real_only<eT>::result* junk)
{
arma_extra_debug_sigprint();
arma_ignore(junk);
@@ -43,12 +43,9 @@ spop_norm::mat_norm_2(const SpProxy<T1>& P, const typename arma_real_only<typena
// norm = sqrt( largest eigenvalue of (A^H)*A ), where ^H is the conjugate transpose
// http://math.stackexchange.com/questions/4368/computing-the-largest-eigenvalue-of-a-very-large-sparse-matrix
typedef typename T1::elem_type eT;
typedef typename T1::pod_type T;
typedef typename get_pod_type<eT>::result T;
const unwrap_spmat<typename SpProxy<T1>::stored_type> tmp(P.Q);
const SpMat<eT>& A = tmp.M;
const SpMat<eT>& A = X;
const SpMat<eT> B = trans(A);
const SpMat<eT> C = (A.n_rows <= A.n_cols) ? (A*B) : (B*A);
@@ -56,21 +53,20 @@ spop_norm::mat_norm_2(const SpProxy<T1>& P, const typename arma_real_only<typena
Col<T> eigval;
eigs_sym(eigval, C, 1);
return (eigval.n_elem > 0) ? std::sqrt(eigval[0]) : T(0);
return (eigval.n_elem > 0) ? T(std::sqrt(eigval[0])) : T(0);
}
template<typename T1>
template<typename eT>
inline
typename T1::pod_type
spop_norm::mat_norm_2(const SpProxy<T1>& P, const typename arma_cx_only<typename T1::elem_type>::result* junk)
typename get_pod_type<eT>::result
spop_norm::mat_norm_2(const SpMat<eT>& X, const typename arma_cx_only<eT>::result* junk)
{
arma_extra_debug_sigprint();
arma_ignore(junk);
typedef typename T1::elem_type eT;
typedef typename T1::pod_type T;
typedef typename get_pod_type<eT>::result T;
// we're calling eigs_gen(), which currently requires ARPACK
#if !defined(ARMA_USE_ARPACK)
@@ -80,9 +76,7 @@ spop_norm::mat_norm_2(const SpProxy<T1>& P, const typename arma_cx_only<typename
}
#endif
const unwrap_spmat<typename SpProxy<T1>::stored_type> tmp(P.Q);
const SpMat<eT>& A = tmp.M;
const SpMat<eT>& A = X;
const SpMat<eT> B = trans(A);
const SpMat<eT> C = (A.n_rows <= A.n_cols) ? (A*B) : (B*A);
@@ -90,20 +84,20 @@ spop_norm::mat_norm_2(const SpProxy<T1>& P, const typename arma_cx_only<typename
Col<eT> eigval;
eigs_gen(eigval, C, 1);
return (eigval.n_elem > 0) ? std::sqrt(std::real(eigval[0])) : T(0);
return (eigval.n_elem > 0) ? T(std::sqrt(std::real(eigval[0]))) : T(0);
}
template<typename T1>
template<typename eT>
inline
typename T1::pod_type
spop_norm::mat_norm_inf(const SpProxy<T1>& P)
typename get_pod_type<eT>::result
spop_norm::mat_norm_inf(const SpMat<eT>& X)
{
arma_extra_debug_sigprint();
// TODO: this can be sped up with a dedicated implementation
return as_scalar( max( sum(abs(P.Q), 1), 0) );
return as_scalar( max( sum(abs(X), 1), 0) );
}