check for nonfinite matrices

include/armadillo_bits/glue_solve_meat.hpp

@@ -130,6 +130,8 @@ glue_solve_gen_full::apply(Mat<eT>& out, const Base<eT,T1>& A_expr, const Base<e
     if(refine)        { arma_debug_warn_level(2,  "solve(): option 'refine' ignored for forced approximate solution"       ); }
     if(likely_sympd)  { arma_debug_warn_level(2,  "solve(): option 'likely_sympd' ignored for forced approximate solution" ); }
     
+    if(arma_config::check_nonfinite && A.has_nonfinite())  { arma_debug_warn_level(3, "solve(): given matrix has non-finite elements"); return false; }
+    
     return auxlib::solve_approx_svd(out, A, B_expr.get_ref());  // A is overwritten
     }
   
@@ -150,6 +152,18 @@ glue_solve_gen_full::apply(Mat<eT>& out, const Base<eT,T1>& A_expr, const Base<e
     
     const bool try_sympd = arma_config::optimise_sympd && ((no_sympd || auxlib::crippled_lapack(A) || is_band || is_triu || is_tril) ? false : (likely_sympd ? true : sympd_helper::guess_sympd(A, uword(16))));
     
+    if(arma_config::check_nonfinite)
+      {
+      bool has_nonfinite = false;
+      
+           if(is_triu  )  { has_nonfinite = trimat_helper::has_nonfinite_triu(A); }
+      else if(is_tril  )  { has_nonfinite = trimat_helper::has_nonfinite_tril(A); }
+      else if(try_sympd)  { has_nonfinite = trimat_helper::has_nonfinite_tril(A); }
+      else                { has_nonfinite = A.has_nonfinite();                    }
+            
+      if(has_nonfinite)  { arma_debug_warn_level(3, "solve(): given matrix has non-finite elements"); return false; }
+      }
+    
     if(fast)
       {
       // fast mode: solvers without refinement and without rcond estimate
@@ -322,6 +336,8 @@ glue_solve_gen_full::apply(Mat<eT>& out, const Base<eT,T1>& A_expr, const Base<e
     if(refine)        { arma_debug_warn_level(2,  "solve(): option 'refine' ignored for non-square matrix"       ); }
     if(likely_sympd)  { arma_debug_warn_level(2,  "solve(): option 'likely_sympd' ignored for non-square matrix" ); }
     
+    if(arma_config::check_nonfinite && A.has_nonfinite())  { arma_debug_warn_level(3, "solve(): given matrix has non-finite elements"); return false; }
+    
     if(fast)
       {
       status = auxlib::solve_rect_fast(out, A, B_expr.get_ref());  // A is overwritten
@@ -406,6 +422,16 @@ glue_solve_tri_default::apply(Mat<eT>& actual_out, const Base<eT,T1>& A_expr, co
   
   arma_debug_check( (A.is_square() == false), "solve(): matrix marked as triangular must be square sized" );
   
+  if(arma_config::check_nonfinite)
+    {
+    bool has_nonfinite = false;
+    
+    if(triu)  { has_nonfinite = trimat_helper::has_nonfinite_triu(A); }
+    if(tril)  { has_nonfinite = trimat_helper::has_nonfinite_tril(A); }
+          
+    if(has_nonfinite)  { arma_debug_warn_level(3, "solve(): given matrix has non-finite elements"); return false; }
+    }
+  
   const uword layout   = (triu) ? uword(0) : uword(1);
   const bool  is_alias = U.is_alias(actual_out);
   
@@ -518,6 +544,16 @@ glue_solve_tri_full::apply(Mat<eT>& actual_out, const Base<eT,T1>& A_expr, const
   
   arma_debug_check( (A.is_square() == false), "solve(): matrix marked as triangular must be square sized" );
   
+  if(arma_config::check_nonfinite)
+    {
+    bool has_nonfinite = false;
+    
+    if(triu)  { has_nonfinite = trimat_helper::has_nonfinite_triu(A); }
+    if(tril)  { has_nonfinite = trimat_helper::has_nonfinite_tril(A); }
+          
+    if(has_nonfinite)  { arma_debug_warn_level(3, "solve(): given matrix has non-finite elements"); return false; }
+    }
+  
   const uword layout   = (triu) ? uword(0) : uword(1);
   const bool  is_alias = U.is_alias(actual_out);
   

include/armadillo_bits/trimat_helper.hpp

@@ -107,6 +107,58 @@ is_tril(const Mat<eT>& A)
 
 
 
+template<typename eT>
+inline
+bool
+has_nonfinite_tril(const Mat<eT>& A)
+  {
+  arma_extra_debug_sigprint();
+  
+  // NOTE: assuming that A has a square size
+  
+  const eT*   colptr = A.memptr();
+  const uword N      = A.n_rows;
+  
+  for(uword i=0; i<N; ++i)
+    {
+    const uword len = N-i;
+    
+    if(arrayops::is_finite(&(colptr[i]), len) == false)  { return true; }
+    
+    colptr += N;
+    }
+  
+  return false;
+  }
+
+
+
+template<typename eT>
+inline
+bool
+has_nonfinite_triu(const Mat<eT>& A)
+  {
+  arma_extra_debug_sigprint();
+  
+  // NOTE: assuming that A has a square size
+  
+  const eT*   colptr = A.memptr();
+  const uword N      = A.n_rows;
+  
+  for(uword i=0; i<N; ++i)
+    {
+    const uword len = i+1;
+    
+    if(arrayops::is_finite(colptr, len) == false)  { return true; }
+    
+    colptr += N;
+    }
+  
+  return false;
+  }
+
+
+
 }  // end of namespace trimat_helper