From bug-octave-request at bevo dot che dot wisc dot edu Wed Dec 18 08:01:34 2002
Subject: Re: 0*NaN within matrix multiplication
From: Paul Kienzle <pkienzle at jazz dot ncnr dot nist dot gov>
To: "Lippert, Ross A." <Ross dot Lippert at celera dot com>
Cc: "John W. Eaton" <jwe at bevo dot che dot wisc dot edu>,    Schloegl Alois <alois.schloegl@tugraz.at>, bug-octave@bevo.che.wisc.edu
Date: Wed, 18 Dec 2002 09:01:32 -0500

If we are going to deal with it I think we should deal with it in
the BLAS code itself.  You could start by writing a small test
program which demostrates the problem and sending it to 
lapack at cs dot utk dot edu dot 

Paul Kienzle
pkienzle at users dot sf dot net

On Wed, Dec 18, 2002 at 08:52:44AM -0500, Lippert, Ross A. wrote:
> 
> Does anyone think it is a bad enough bug to want to have
> some sort of coding around this?
> 
> -r
> 
> -----Original Message-----
> From: John W. Eaton [mailto:jwe at bevo dot che dot wisc dot edu]
> Sent: Tuesday, December 17, 2002 7:40 PM
> To: Lippert, Ross A.
> Cc: Paul Kienzle; Schloegl Alois; bug-octave at bevo dot che dot wisc dot edu
> Subject: RE: 0*NaN within matrix multiplication
> 
> 
> On 10-Dec-2002, Lippert, Ross A. <Ross dot Lippert at celera dot com> wrote:
> 
> | Ahh.  In that case the bug is smaller (and cuter):
> | octave:17> [nan]*[0]
> | ans = NaN
> | octave:18> [nan nan]*[0;0]
> | ans = 0
> | 
> | Yes, I agree that the above is a bad thing(tm).
> 
> This is what you get with the reference implementation of DGEMM (it
> has some logic to avoid the multiplication if one of the operands is
> zero), but not with the ATLAS implementation (it may also avoid
> multiplication by zero, but apparently not if the other operand is NaN
> or Inf):
> 
>   octave:1> [nan nan]*[0;0]
>   ans = NaN
>   octave:2> [inf inf]*[0;0]
>   ans = NaN
> 
>   $ ldd `type -p octave`
>     liboctinterp.so => /usr/lib/octave-2.1.40/liboctinterp.so (0x40014000)
>     liboctave.so => /usr/lib/octave-2.1.40/liboctave.so (0x40431000)
>     libcruft.so => /usr/lib/octave-2.1.40/libcruft.so (0x40640000)
>     liblapack.so.2 => /usr/lib/atlas/liblapack.so.2 (0x406a0000)
>     libblas.so.2 => /usr/lib/atlas/libblas.so.2 (0x40bcd000)
>     libfftw.so.2 => /usr/lib/libfftw.so.2 (0x40e6b000)
>     libreadline.so.4 => /lib/libreadline.so.4 (0x40e9f000)
>     libncurses.so.5 => /lib/libncurses.so.5 (0x40ec9000)
>     libdl.so.2 => /lib/libdl.so.2 (0x40f07000)
>     libm.so.6 => /lib/libm.so.6 (0x40f0a000)
>     libstdc++-libc6.2-2.so.3 => /usr/lib/libstdc++-libc6.2-2.so.3 (0x40f2b000)
>     libc.so.6 => /lib/libc.so.6 (0x40f74000)
>     /lib/ld-linux.so.2 => /lib/ld-linux.so.2 (0x40000000)
> 
> I don't know what the various vendor versions of the blas do.
> 
> jwe



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------