# bug #36562: Sparse matrix operations with NaN do not follow IEEE standard

 Submitter: None Submitted: Wed 30 May 2012 03:12:51 PM UTC Category: Libraries Severity: 3 - Normal Priority: 5 - Normal Item Group: Incorrect Result Status: Confirmed Assigned to: None Originator Name: Stuart Archibald Originator Email: -email is unavailable- Open/Closed:  * Open Release:  * dev Operating System:  * Any
* Mandatory Fields

Comment Type & Canned Response:

## Discussion

Mon 21 Nov 2016 07:42:22 PM UTC, comment #7:

At least we understand what behavior needs to be implemented.  At a minimum, Octave needs to check whether each operand contains a NaN and it needs to specifically guarantee that those locations result in NaN at the output.  The subsequent decision, about whether to convert to full if it would save space, should probably be left up to Octave, or to the sparse_auto_mutate code.

Pseudo-code which might work for Z = X op Y:

list1 = isnan (X);
list2 = isnan (Y);
tmp = X op Y;
tmp(list1) = NaN;
tmp(list2) = NaN;
Z = tmp;

Of course, all of that would need to be in C++.

Another choice would be to look carefully at the code which performs the operation 'op'.  For sparse matrices, I believe the code is equivalent to

nonzero1 = nonzeros (X);
nonzero2 = nonzeros (Y);
elements_to_operate_on = intersection (nonzero1, nonzero2);
Z = X(elements_to_operate_on) op Y(elements_to_operate_on);

The "intersection" function for sparse matrices in C++ could be changed to include elements where only one of the matrices had a NaN.

Rik <rik5>
Mon 21 Nov 2016 05:49:15 PM UTC, comment #6:

Matlab R2016a:

>> xs = sparse ([1 0; 0 2]);
ys = sparse ([1 NaN; 0 0]);
xs .* ys

ans =

(1,1)        1
(1,2)      NaN

Rik, I understand your suggestion may not be worse than the current status, but it seems the memory issue does happen in Matlab.

Carlo de Falco <cdf>
Mon 21 Nov 2016 04:44:23 PM UTC, comment #5:

@Carlo: The memory problem already exists.  After one of the operations involving NaN, every single element will become NaN and require storage.  In fact, the situation will be worse than a full matrix because there will be the extra overhead of the sparse matrix implementation.

In general, even without sparse_auto_mutate set to true, Octave will convert to full for operations that are likely to result in a full matrix.  From the documentation,

22.1.4.2 Return Types of Operators and Functions

The two basic reasons to use sparse matrices are to reduce the memory usage and to not have
to do calculations on zero elements. The two are closely related in that the computation
time on a sparse matrix operator or function is roughly linear with the number of nonzero
elements.

Therefore, there is a certain density of nonzero elements of a matrix where it no longer
makes sense to store it as a sparse matrix, but rather as a full matrix. For this reason
operators and functions that have a high probability of returning a full matrix will always
return one. For example adding a scalar constant to a sparse matrix will almost always
make it a full matrix, and so the example,

speye (3) + 0

1 0 0
0 1 0
0 0 1
returns a full matrix as can be seen.

For addition and subtraction, the sparse matrices are already converted to full and the NaN is added/subtracted and the operation follows IEEE guidelines.  For Sparse Matrix OPERATION Scalar NaN, the code would be easy to modify to check for NaN and then return

NaN (size (Xsparse))

For Matrix OPERATION Matrix operations it becomes a bit annoying because one would have to check for isnan() on every element.  This may be what Matlab does.  Can you test the following?

xs = sparse ([1 0; 0 2]);
ys = sparse ([1 NaN; 0 0]);
xs .* ys

Rik <rik5>
Mon 21 Nov 2016 10:48:04 AM UTC, comment #4:

Rik :

>> Msparse = sparse ([1 0; 0 3]);
NaN .* Msparse

ans =

(1,1)      NaN
(2,1)      NaN
(1,2)      NaN
(2,2)      NaN

so your solution does not seem to be compatible.

Anyway, I wouldn't like that solution even if it were compatible
as any error in filling the elements of a large sparse matrix would result in out-of-memory.

Carlo de Falco <cdf>
Sun 20 Nov 2016 03:25:31 PM UTC, comment #3:

Any obvious solution would be for NaN to force all matrices into a full implementation.

What does Matlab do for

Msparse = sparse ([1 0; 0 3]);
NaN .* Msparse

Rik <rik5>
Sat 19 Nov 2016 09:45:39 PM UTC, comment #2:

This issue still exists with Octave 4.2.0.

Hartmut <hardy>
Thu 31 May 2012 08:43:26 PM UTC, comment #1:

Confirmed on a recent (5/31/12) development tip.

The issue is more general than just divide so I have changed the bug Summary.

At present, I can verify that multiply and divide are affected while addition and subtraction behaves correctly.

The issue appears to be a side effect of the way that sparse matrices improve calculation performance.  In general, math operations like multiply and divide only need to be performed on the non-zero elements of a sparse matrix.  However, a NaN is different and should affect all elements of the matrix.

Below is some simple code to illustrate the problem.

Mfull = [1 0; 0 3];
Msparse = sparse (Mfull);

NaN .* Mfull
ans =

NaN   NaN
NaN   NaN

NaN .* Msparse
ans =

Compressed Column Sparse (rows = 2, cols = 2, nnz = 2 [50%])

(1, 1) -> NaN
(2, 2) -> NaN

NaN ./ Mfull
ans =

NaN   NaN
NaN   NaN

NaN ./ Msparse

ans =

NaN  -Inf
-Inf   NaN

Rik <rik5>
Wed 30 May 2012 03:12:51 PM UTC, original submission:

rdivide with NaN as first argument does not behave as expected with sparse array as second argument. NaN/0 = NaN, which implies NaN rdivided by sparse should give NaN's in the positions of the the zero pattern as well as NaN's in the non-zero pattern.

Example

octave:1> rdivide(NaN,sparse([1,3],[1,3],[1,3]))
ans =

NaN  -Inf  -Inf
-Inf  -Inf  -Inf
-Inf  -Inf   NaN

I would expect the answer to be:

ans =
NaN   NaN   NaN
NaN   NaN   NaN
NaN   NaN   NaN

Anonymous

## Attached Files

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:

Comment:

Attached Files
file #25956:  sparseNaNbug.m added by None (68B - text/x-objcsrc - Script demonstrates the problem)

## Dependencies

Depends on the following items: None found

Items that depend on this one: None found

Carbon-Copy List
• -email is unavailable- added by cdf (Posted a comment)
• -email is unavailable- added by hardy (Posted a comment)
• -email is unavailable- added by mtmiller (Updated the item)
• -email is unavailable- added by rik5 (Posted a comment)
• -email is unavailable- added by None (Submitted the item)
•

There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

Only project members can vote.