[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty s
From: |
Rik |
Subject: |
[Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix. |
Date: |
Sat, 18 May 2019 17:45:42 -0400 (EDT) |
User-agent: |
Mozilla/5.0 (Windows NT 10.0; WOW64; Trident/7.0; rv:11.0) like Gecko |
Update of bug #56232 (project octave):
Status: Patch Reviewed => Fixed
_______________________________________________________
Follow-up Comment #62:
Thanks to Marco for taking the lead in resolving this problem. I made some
modifications for syntax or for performance, and checked it in here
(https://hg.savannah.gnu.org/hgweb/octave/rev/6e18f0ce268c).
Documenting my additions, I changed the BIST tests which were using some side
effects of our implementation to work, to use more ordinary syntax.
Instead of
+%!test
+%!warning <matrix singular> A = inv (sparse ([1, 2;0 ,0]));
+%! assert (A, sparse ([Inf, Inf; 0, 0]));
I wrote
+%!test
+%! fail ("A = inv (sparse ([1, 2;0 ,0]))", "warning", "matrix singular");
+%! assert (A, sparse ([Inf, Inf; 0, 0]));
For the diagonal matrices, the implementation relied on creating a local
buffer that could hold up to nnz() elements. The algorithm also needed to
cycle over the entire number of diagonal elements before arriving at a
decision about whether the matrix was singular. The initial implementation is
below.
info = 0;
OCTAVE_LOCAL_BUFFER (octave_idx_type, nzelem, len);
octave_idx_type count_nz = 0;
for (octave_idx_type i = 0; i < len; i++)
{
if (elem (i, i) == 0.0)
retval.elem (i, i) = octave::numeric_limits<double>::Inf ();
else
{
retval.elem (i, i) = 1.0 / elem (i, i);
nzelem[count_nz] = i;
count_nz++;
}
}
if (count_nz == 0)
{
(*current_liboctave_error_handler)
("inverse of the null matrix not defined");
}
else if (count_nz < len)
{
info = -1;
for (octave_idx_type i = 0; i < count_nz; i++)
retval.elem (nzelem[i], nzelem[i]) = octave::numeric_limits<double>::Inf
();
}
Instead, I found that it was possible with two count variables to eliminate
the OCTAVE_LOCAL_BUFFER and to short-circuit out of the loop as soon as a zero
is detected. For a very small speed-up, I also used xelem() to access
elements where I could be certain that I was not going to go outside the
bounds of the array.
info = 0;
octave_idx_type len = r; // alias for readability
octave_idx_type z_count = 0; // zeros
octave_idx_type nz_count = 0; // non-zeros
for (octave_idx_type i = 0; i < len; i++)
{
if (xelem (i, i) == 0.0)
{
z_count++;
if (nz_count > 0)
break;
}
else
{
nz_count++;
if (z_count > 0)
break;
retval.elem (i, i) = 1.0 / xelem (i, i);
}
}
if (nz_count == 0)
{
(*current_liboctave_error_handler)
("inverse of the null matrix not defined");
}
else if (z_count > 0)
{
info = -1;
element_type *data = retval.fortran_vec ();
std::fill (data, data + len, octave::numeric_limits<double>::Inf ());
}
For Sparse and CompleSparse arrays, the original additional code was
+ if (rcond == 0.)
+ {
+ for (octave_idx_type i = 0; i < cols () + 1; i++)
+ ret.cidx (i) = cidx (i);
+
+ for (octave_idx_type i = 0; i < nnz (); i++)
+ {
+ ret.data (i) = octave::numeric_limits<double>::Inf ();
+ ret.ridx (i) = ridx (i);
+ }
+
+ return ret;
+ }
Instead of using hand-rolled for loops, I decided on functions from the C++
standard library.
+ if (rcond == 0.0)
+ {
+ // Return all Inf matrix with sparsity pattern of input.
+ octave_idx_type nz = nnz ();
+ ret = SparseComplexMatrix (rows (), cols (), nz);
+ std::fill (ret.xdata (), ret.xdata () + nz,
+ octave::numeric_limits<double>::Inf ());
+ std::copy_n (ridx (), nz, ret.xridx ());
+ std::copy_n (cidx (), cols () + 1, ret.xcidx ());
+
+ return ret;
+ }
_______________________________________________________
Reply to this item at:
<https://savannah.gnu.org/bugs/?56232>
_______________________________________________
Message sent via Savannah
https://savannah.gnu.org/
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., (continued)
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Kai Torben Ohlhus, 2019/05/10
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Rik, 2019/05/13
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Rik, 2019/05/13
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Marco Caliari, 2019/05/14
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Marco Caliari, 2019/05/14
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Marco Caliari, 2019/05/15
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Rik, 2019/05/15
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Kai Torben Ohlhus, 2019/05/15
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Marco Caliari, 2019/05/16
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix., Kai Torben Ohlhus, 2019/05/17
- [Octave-bug-tracker] [bug #56232] Octave crash when inverting an empty sparse matrix.,
Rik <=