bugGNU Octave - Bugs: bug #45387, Matrix product treats NaN times...

 
 

bug #45387: Matrix product treats NaN times zero as zero

Submitter:  None
Submitted:  Tue 23 Jun 2015 08:44:35 PM UTC
   
 
Category:  Documentation Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  None
Originator Name:  Pete Originator Email:  -email is unavailable-
Open/Closed:  * Closed Release:  * 4.0.0
Operating System:  * Any Fixed Release:  None
Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Sat 25 Jul 2015 04:13:09 AM UTC, comment #14: 

Closing as fixed in the documentation since the underlying bug is out of Octave's control.

Mike Miller <mtmiller>
Group Member
Wed 08 Jul 2015 02:11:13 AM UTC, comment #13: 

I added a warning about this issue to the documentation (http://hg.savannah.gnu.org/hgweb/octave/rev/3300cb76cccd).  That's probably as much as we want to do for this.

Rik <rik5>
Group administrator
Tue 07 Jul 2015 03:56:43 PM UTC, comment #12: 

Again, not many people will probably read it, but we could document this is the Known Bugs section in the manual.

Rik <rik5>
Group administrator
Mon 06 Jul 2015 11:11:29 PM UTC, comment #11: 

It would be possible to perform a check for this problem when Octave starts and issue a warning, but it would probably just be more alarming to new users than ignoring the problem.

John W. Eaton <jwe>
Group administrator
Mon 06 Jul 2015 03:37:34 PM UTC, comment #10: 

Well, to clarify further, Octave packages in Debian are installed with OpenBLAS or ATLAS by default for users, but when they are built they are compiled against netlib BLAS. IMHO nothing should change at the configure/make stage. I don't think there's anything wrong with building against netlib BLAS, and even using it, works fine (yeah, slow) for most purposes. This is a very narrowly scoped bug/feature, doesn't render the entire netlib library useless. I think documentation is the right fix, if this bug report isn't documentation enough I am not sure where but it could be documented in the manual as Rik suggested.

Mike Miller <mtmiller>
Group Member
Mon 06 Jul 2015 07:30:18 AM UTC, comment #9: 

Then, it remains only people compiling from sources. It would be possible to check if netlib BLAS is used (for instance, by testing this bug with NaN propagation) and then stop the configuration and make a warning. If one really wants to use netlib BLAS, there should be a configure option, something like --use-netlib-blas.

Marco

P.S. By the way, it is not possible to compile from sources with ATLAS (umfpack is not found). But this is another bug.

Marco Caliari <caliari>
Group Member
Sat 04 Jul 2015 11:27:41 PM UTC, comment #8: 

I can comment on Debian and Ubuntu. The ATLAS and OpenBLAS libraries are both listed as "Recommends" dependencies of Octave, so one of them will be installed by default for all users of Octave, and I believe OpenBLAS will be preferred for all new installs. If an expert user disables recommended dependencies or intentionally uninstalls ATLAS and OpenBLAS, then they might wind up with netlib BLAS instead.

IIRC the Fedora packagers have chosen a hard dependency on ATLAS so there is no chance of using netlib BLAS even by accident.

Mike Miller <mtmiller>
Group Member
Sat 04 Jul 2015 12:13:17 PM UTC, comment #7: 

What should we say about reference BLAS and where should we say it?  I don't think NEWS is the right place for this because that file is documenting changes.  And anything we put in NEWS today would get aged out into the file NEWS.4 when Octave changes to version 5.0.

In Appendix G, Installing Octave, we list dependencies.  This is the wording for BLAS:


BLAS
 Basic Linear Algebra Subroutine library (http://www.netlib.org/blas). Ac-
celerated blas libraries such as ATLAS (http://math-atlas.sourceforge.net)
are recommeded for better performance.


Maybe we should be downgrading reference BLAS to the lowest recommended alternative.  However, most people don't start by reading Appendix G.  In fact, most of the time I don't think anyone starts with the manual.  They just install the software from their package manager and they get whatever library the Linux distribution is supplying.  Maybe we need to start by asking some of the distro packagers, like Debian, to make OpenBLAS or ATLAS a dependency for Octave.  I'm not sure if that can work, but Mike Miller can comment on that.



Rik <rik5>
Group administrator
Fri 03 Jul 2015 07:17:48 AM UTC, comment #6: 

It is reference BLAS. If you look at line 281 here

http://www.netlib.org/lapack/explore-html/dc/da8/dgemv_8f_source.html

you see that if the entry of the vector is zero, the multiplication is skipped. I tried the attached fortran file and I confirm the correct behavior for openblas and ATLAS and the wrong for reference.

I think that something should be written, otherwise a Linux user may have the impression that 1) octave cannot properly handle NaN and 2) (much more important IMHO) octave is much slower than Matlab.

Marco

(file #34369)

Marco Caliari <caliari>
Group Member
Thu 02 Jul 2015 03:27:16 PM UTC, comment #5: 

Do we have confirmation that it is definitely the BLAS library and not anything Octave is doing in how it performs the library call?

Is there a bug report against the reference BLAS library for this particular behavior that we can cite?

Should we be warning users in the NEWS file, or somewhere else, about this issue?

During the Windows install, the default is OpenBLAS so you have to work harder to get the reference BLAS library.  But on Linux distros, as this bug was reported on, the default BLAS library often is the reference version.

I can also confirm that switching to OpenBLAS fixes the issue.

Rik <rik5>
Group administrator
Thu 02 Jul 2015 01:39:35 PM UTC, comment #4: 

Confirmed here also, with netlib blas I get


octave:1> [nan, 1] * [0; 1]
ans = NaN
octave:2> [nan, 1; 0, 0] * [0; 1]
ans =

   1
   0


and with openBLAS I get


octave:1> [nan, 1] * [0; 1]
ans = NaN
octave:2> [nan, 1; 0, 0] * [0; 1]
ans =

   NaN
     0


Mike Miller <mtmiller>
Group Member
Thu 02 Jul 2015 09:09:00 AM UTC, comment #3: 

It is a problem of (standard) libblas.so.3. In fact, if I select openblas or atlas, I see the correct result [NaN;1].

Marco

P.S. If I'm right, why are you using the slow standard BLAS?

Marco Caliari <caliari>
Group Member
Wed 24 Jun 2015 05:47:03 PM UTC, comment #2: 

More information, this behavior goes way back.  I have tested through stable versions back to 3.2.4 and I always get the wrong result.  This rules out automatic broadcasting which I thought might be related, but was added later than 3.2.4.

Rik <rik5>
Group administrator
Wed 24 Jun 2015 04:51:45 AM UTC, comment #1: 

Confirmed.

There is some sort of inconsistency because a scalar by scalar multiplication


NaN * 0
ans = NaN


and vector by vector multiplication works as you have demonstrated.  Also, matrix by matrix multiplication yields NaNs in the correct places too.


x2 * x2
ans =

   NaN   NaN
   NaN     3


It seems that only matrix by vector multiplication is inconsistent.  And that has something to do with the element being exactly 0 as you have noticed.  Shifting the 0 to eps gives the correct answer


y2 = [eps; 1];
x2 * y2
ans =

      NaN
   1.0000


In fact, its even stranger than that as vector by matrix works.


y' * x2
ans =

   NaN     1



Rik <rik5>
Group administrator
Tue 23 Jun 2015 08:44:35 PM UTC, original submission:  

Start with this


x1 = [NaN, 2]
y = [0;1]
r1 = x1*y


Result is correct
r1 = NaN

But, if we instead use

x2 = [NaN, 2;1,1]
y = [0;1]
r2 = x2*y


The result is
    r2 = [2;1]
which seems wrong (or at least inconsistent with result r1). I would have expected
    r2 = [NaN;1]
which is what Matlab produces.

To summarize,
    x1*y = NaN times 0 plus 2 times 1 = NaN
but
    x2*y = [NaN times 0 + 2 times 1; 1 times 0 + 1 times 1]
            = [2;1]
instead of [NaN;1].


Anonymous

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #34369:  main_dgemv.f added by caliari (523B - text/x-fortran)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by jwe (Posted a comment)
  • -email is unavailable- added by mtmiller (Posted a comment)
  • -email is unavailable- added by caliari (Posted a comment)
  • -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 group members can vote.

     

    Follow 10 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2015-08-14 rik5 Open/ClosedOpen Closed
    2015-07-25 mtmiller CategoryInterpreter Documentation
        StatusWont Fix Fixed
        Operating SystemGNU/Linux Any
    2015-07-23 rik5 StatusConfirmed Wont Fix
    2015-07-07 rik5 Severity4 - Important 3 - Normal
    2015-07-03 caliari Attached File- Added main_dgemv.f, #34369
    2015-06-24 rik5 Severity3 - Normal 4 - Important
    2015-06-24 rik5 StatusNone Confirmed
        SummaryMatrix Product treats NaN times zero as zero Matrix product treats NaN times zero as zero

    Back to the top

    Powered by Savane 3.13-3230.
    Corresponding source code