@@ -50,19 +50,26 @@ guess_sympd_worker(const Mat<eT>& A)
5050 const eT* A_mem = A.memptr ();
5151 const eT* A_col = A_mem;
5252
53+ bool diag_below_tol = true ;
54+
5355 eT max_diag = eT (0 );
5456
5557 for (uword j=0 ; j < N; ++j)
5658 {
5759 const eT A_jj = A_col[j];
5860
59- if (A_jj <= eT (0 )) { return false ; }
61+ if ( A_jj <= eT (0 )) { return false ; }
62+ if (arma_isfinite (A_jj) == false ) { return false ; }
63+
64+ if (A_jj >= tol) { diag_below_tol = false ; }
6065
6166 max_diag = (A_jj > max_diag) ? A_jj : max_diag;
6267
6368 A_col += N;
6469 }
6570
71+ if (diag_below_tol) { return false ; } // assume matrix is suspect if all diagonal elements are close to zero
72+
6673 A_col = A_mem;
6774
6875 const uword Nm1 = N-1 ;
@@ -128,6 +135,8 @@ guess_sympd_worker(const Mat<eT>& A)
128135 const eT* A_mem = A.memptr ();
129136 const eT* A_col = A_mem;
130137
138+ bool diag_below_tol = true ;
139+
131140 T max_diag = T (0 );
132141
133142 for (uword j=0 ; j < N; ++j)
@@ -138,15 +147,21 @@ guess_sympd_worker(const Mat<eT>& A)
138147 const T A_jj_rabs = std::abs (A_jj_r);
139148 const T A_jj_iabs = std::abs (A_jj_i);
140149
141- if (A_jj_r <= T (0 ) ) { return false ; } // real should be positive
142- if (A_jj_iabs > tol ) { return false ; } // imag should be approx zero
143- if (A_jj_iabs > A_jj_rabs) { return false ; } // corner case: real and imag are close to zero, and imag is dominant
150+ if ( A_jj_r <= T (0 ) ) { return false ; } // real should be positive
151+ if (arma_isfinite (A_jj_r) == false ) { return false ; }
152+
153+ if (A_jj_iabs > tol ) { return false ; } // imag should be approx zero
154+ if (A_jj_iabs > A_jj_rabs) { return false ; } // corner case: real and imag are close to zero, and imag is dominant
155+
156+ if (A_jj_r >= tol) { diag_below_tol = false ; }
144157
145158 max_diag = (A_jj_r > max_diag) ? A_jj_r : max_diag;
146159
147160 A_col += N;
148161 }
149162
163+ if (diag_below_tol) { return false ; } // assume matrix is suspect if all diagonal elements are close to zero
164+
150165 const T square_max_diag = max_diag * max_diag;
151166
152167 if (arma_isfinite (square_max_diag) == false ) { return false ; }
@@ -264,15 +279,21 @@ is_approx_sym_worker(const Mat<eT>& A)
264279 const eT* A_mem = A.memptr ();
265280 const eT* A_col = A_mem;
266281
282+ bool diag_below_tol = true ;
283+
267284 for (uword j=0 ; j < N; ++j)
268285 {
269- const eT& A_jj = A_col[j];
286+ const eT A_jj = A_col[j];
270287
271288 if (arma_isfinite (A_jj) == false ) { return false ; }
272289
290+ if (std::abs (A_jj) >= tol) { diag_below_tol = false ; }
291+
273292 A_col += N;
274293 }
275294
295+ if (diag_below_tol) { return false ; } // assume matrix is suspect if all diagonal elements are close to zero
296+
276297 A_col = A_mem;
277298
278299 const uword Nm1 = N-1 ;
@@ -324,6 +345,8 @@ is_approx_sym_worker(const Mat<eT>& A)
324345 const eT* A_mem = A.memptr ();
325346 const eT* A_col = A_mem;
326347
348+ bool diag_below_tol = true ;
349+
327350 // ensure diagonal has approx real-only elements
328351 for (uword j=0 ; j < N; ++j)
329352 {
@@ -338,9 +361,13 @@ is_approx_sym_worker(const Mat<eT>& A)
338361
339362 if (arma_isfinite (A_jj_r) == false ) { return false ; }
340363
364+ if (A_jj_rabs >= tol) { diag_below_tol = false ; }
365+
341366 A_col += N;
342367 }
343368
369+ if (diag_below_tol) { return false ; } // assume matrix is suspect if all diagonal elements are close to zero
370+
344371 A_col = A_mem;
345372
346373 const uword Nm1 = N-1 ;
0 commit comments