Mercurial > hg > octave-nkf
diff liboctave/dSparse.cc @ 5876:565d0cd4d9d0
[project @ 2006-07-01 19:42:06 by dbateman]
author | dbateman |
---|---|
date | Sat, 01 Jul 2006 19:42:07 +0000 |
parents | 6b9cec830d72 |
children | 50d43cdbec80 |
line wrap: on
line diff
--- a/liboctave/dSparse.cc +++ b/liboctave/dSparse.cc @@ -827,6 +827,14 @@ octave_idx_type colXp = retval.xcidx(i); octave_idx_type colUp = cidx(j); octave_idx_type rpX, rpU; + + if (cidx(j) == cidx(j+1)) + { + (*current_liboctave_error_handler) + ("division by zero"); + goto inverse_singular; + } + do { OCTAVE_QUIT; @@ -847,11 +855,17 @@ (colXp<cx) && (colUp<nz)); // get A(m,m) - colUp = cidx(j+1) - 1; + if (typ == MatrixType::Upper) + colUp = cidx(j+1) - 1; + else + colUp = cidx(j) - 1; double pivot = data(colUp); - if (pivot == 0.) - (*current_liboctave_error_handler) - ("division by zero"); + if (pivot == 0. || colUp != j) + { + (*current_liboctave_error_handler) + ("division by zero"); + goto inverse_singular; + } if (v != 0.) { @@ -868,10 +882,17 @@ } // get A(m,m) - octave_idx_type colUp = cidx(i+1) - 1; + octave_idx_type colUp; + if (typ == MatrixType::Upper) + colUp = cidx(i+1) - 1; + else + colUp = cidx(i) - 1; double pivot = data(colUp); - if (pivot == 0.) - (*current_liboctave_error_handler) ("division by zero"); + if (pivot == 0. || colUp != i) + { + (*current_liboctave_error_handler) ("division by zero"); + goto inverse_singular; + } if (pivot != 1.0) for (octave_idx_type j = cx_colstart; j < cx; j++) @@ -929,20 +950,35 @@ } // get A(m,m) - double pivot = data(cidx(jidx+1) - 1); + double pivot; + if (typ == MatrixType::Permuted_Upper) + pivot = data(cidx(jidx+1) - 1); + else + pivot = data(cidx(jidx) - 1); if (pivot == 0.) - (*current_liboctave_error_handler) - ("division by zero"); + { + (*current_liboctave_error_handler) + ("division by zero"); + goto inverse_singular; + } work[j] = v / pivot; } // get A(m,m) - octave_idx_type colUp = cidx(perm[iidx]+1) - 1; + octave_idx_type colUp; + if (typ == MatrixType::Permuted_Upper) + colUp = cidx(perm[iidx]+1) - 1; + else + colUp = cidx(perm[iidx]) - 1; + double pivot = data(colUp); - if (pivot == 0.) - (*current_liboctave_error_handler) - ("division by zero"); + if (pivot == 0.) + { + (*current_liboctave_error_handler) + ("division by zero"); + goto inverse_singular; + } octave_idx_type new_cx = cx; for (octave_idx_type j = iidx; j < nr; j++) @@ -993,6 +1029,9 @@ } return retval; + + inverse_singular: + return SparseMatrix(); } SparseMatrix