use explicit method_id for pinv

This commit is contained in:
conrad
2021-03-11 12:20:32 +10:00
parent 823761e043
commit 3efe33ccc1
3 changed files with 36 additions and 21 deletions
+24 -10
View File
@@ -27,18 +27,26 @@ pinv
(
const Base<typename T1::elem_type,T1>& X,
const typename T1::pod_type tol = 0.0,
const char* method = "dc"
const char* method = nullptr
)
{
arma_extra_debug_sigprint();
typedef typename T1::elem_type eT;
const char sig = (method != nullptr) ? method[0] : char(0);
uword method_id = 0; // default setting
arma_debug_check( ((sig != 's') && (sig != 'd')), "pinv(): unknown method specified" );
if(method != nullptr)
{
const char sig = method[0];
arma_debug_check( ((sig != 's') && (sig != 'd')), "pinv(): unknown method specified" );
if(sig == 's') { method_id = 1; }
if(sig == 'd') { method_id = 2; }
}
return (sig == 'd') ? Op<T1, op_pinv>(X.get_ref(), eT(tol), 1, 0) : Op<T1, op_pinv>(X.get_ref(), eT(tol), 0, 0);
return Op<T1, op_pinv>(X.get_ref(), eT(tol), method_id, uword(0));
}
@@ -51,18 +59,24 @@ pinv
Mat<typename T1::elem_type>& out,
const Base<typename T1::elem_type,T1>& X,
const typename T1::pod_type tol = 0.0,
const char* method = "dc"
const char* method = nullptr
)
{
arma_extra_debug_sigprint();
const char sig = (method != nullptr) ? method[0] : char(0);
uword method_id = 0; // default setting
arma_debug_check( ((sig != 's') && (sig != 'd')), "pinv(): unknown method specified" );
if(method != nullptr)
{
const char sig = method[0];
arma_debug_check( ((sig != 's') && (sig != 'd')), "pinv(): unknown method specified" );
if(sig == 's') { method_id = 1; }
if(sig == 'd') { method_id = 2; }
}
const bool use_divide_and_conquer = (sig == 'd');
const bool status = op_pinv::apply_direct(out, X.get_ref(), tol, use_divide_and_conquer);
const bool status = op_pinv::apply_direct(out, X.get_ref(), tol, method_id);
if(status == false)
{
+2 -2
View File
@@ -26,8 +26,8 @@ class op_pinv
public:
template<typename T1> inline static void apply(Mat<typename T1::elem_type>& out, const Op<T1,op_pinv>& in);
template<typename T1> inline static bool apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::elem_type,T1>& expr, typename T1::pod_type tol, const bool use_divide_and_conquer);
template<typename T1> inline static bool apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::elem_type,T1>& expr, typename T1::pod_type tol, const uword method_id);
};
+10 -9
View File
@@ -29,11 +29,10 @@ op_pinv::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_pinv>& in)
typedef typename T1::pod_type T;
const T tol = access::tmp_real(in.aux);
const T tol = access::tmp_real(in.aux);
const uword method_id = in.aux_uword_a;
const bool use_divide_and_conquer = (in.aux_uword_a == 1);
const bool status = op_pinv::apply_direct(out, in.m, tol, use_divide_and_conquer);
const bool status = op_pinv::apply_direct(out, in.m, tol, method_id);
if(status == false)
{
@@ -46,7 +45,7 @@ op_pinv::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_pinv>& in)
template<typename T1>
inline
bool
op_pinv::apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::elem_type,T1>& expr, typename T1::pod_type tol, const bool use_divide_and_conquer)
op_pinv::apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::elem_type,T1>& expr, typename T1::pod_type tol, const uword method_id)
{
arma_extra_debug_sigprint();
@@ -55,6 +54,10 @@ op_pinv::apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::
arma_debug_check((tol < T(0)), "pinv(): tolerance must be >= 0");
// method_id = 0 -> default setting
// method_id = 1 -> use standard algorithm
// method_id = 2 -> use divide and conquer algorithm
Mat<eT> A(expr.get_ref());
const uword n_rows = A.n_rows;
@@ -63,7 +66,7 @@ op_pinv::apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::
if(A.is_empty()) { out.set_size(n_cols,n_rows); return true; }
#if defined(ARMA_OPTIMISE_SYMPD)
const bool try_sympd = (auxlib::crippled_lapack(A) == false) && (tol == T(0)) && sympd_helper::guess_sympd_anysize(A);
const bool try_sympd = (auxlib::crippled_lapack(A) == false) && (tol == T(0)) && (method_id == uword(0)) && sympd_helper::guess_sympd_anysize(A);
#else
const bool try_sympd = false;
#endif
@@ -89,11 +92,9 @@ op_pinv::apply_direct(Mat<typename T1::elem_type>& out, const Base<typename T1::
Col< T> s;
Mat<eT> V;
bool status = false;
if(n_cols > n_rows) { A = trans(A); }
status = (use_divide_and_conquer) ? auxlib::svd_dc_econ(U, s, V, A) : auxlib::svd_econ(U, s, V, A, 'b');
const bool status = (method_id == uword(2)) ? auxlib::svd_dc_econ(U, s, V, A) : auxlib::svd_econ(U, s, V, A, 'b');
if(status == false) { out.soft_reset(); return false; }