[Bug-gsl] [bug #28767] gsl_linalg_SV_decomp gives NaN when tb, tab, dt = 0

[Brian Gough]

URL:
  <http://savannah.gnu.org/bugs/?28767>

                 Summary: gsl_linalg_SV_decomp gives NaN when tb, tab, dt = 0
                 Project: GNU Scientific Library
            Submitted by: bjg
            Submitted on: Fri 29 Jan 2010 10:04:47 PM GMT
                Category: Runtime error
                Severity: 3 - Normal
        Operating System: 
                  Status: Confirmed
             Assigned to: None
             Open/Closed: Open
                 Release: 1.13
         Discussion Lock: Any

    _______________________________________________________

Details:

--------------------------------------------------------------
GSL version number: 1.9
hardware and operating system, OS details:
  Linux phit1 2.6.25.20-0.5-default #1 SMP 2009-08-14 01:48:11 +0200 x86_64
  x86_64 x86_64 GNU/Linux
gcc version:
  gcc version 4.3.1 20080507 (prerelease) [gcc-4_3-branch revision 135036]
  (SUSE Linux)
genre: linalg, SVD, gsl_linalg_SV_decomp
description of the bug behavior: NaN results in gsl_linalg_SV_decomp
short program which exercises the bug: attached
--------------------------------------------------------------
        
Dear GSL devolopers,

your library is a great help for me and I use it very much. Thank you very
much for your great work!

Recently, I encountered a problem in the SVD routine of GSL version 1.9:
Certain matrices A with very many zeros lead to an SVD with NaN ("not a
number") results in the vector of singular values S as well as in the left
and
right singular vectors.

I tracked the error down to the function trailing_eigenvalue in svdstep.c:
If tb, tab and dt are all zero, the calculated eigenvalue is undefined:

mu = tb - (tab * tab) / (dt + hypot (dt, tab)) = 0 - 0 / 0;
 
This NaN result is then propagated in the function qrstep via create_givens
to
all the other results.

I attach a test program reproducing the problem when applied to a certain
matrix of 462 rows and 421 columns i.e.  194502 elements of which 189694 are
zeros and 4808 are ones, saved in a file of 389 kbyte. The distribution of
ones in the matrix is not random. I can send this matrix file to whoever is
interested (This mailing list refuses attachments of this size).

Here's the output of this program on my machine with the following
architecture:

Linux phit1 2.6.25.20-0.5-default #1 SMP 2009-08-14 01:48:11 +0200 x86_64
x86_64
x86_64 GNU/Linux

--------------------------------------------------------------

(phit1 ~/pj)  gcc `gsl-config --cflags` svdbug.c `gsl-config --libs` -o
svdbug &&
./svdbug
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000000000000000000000000000000000000000000000100010001000000000000000
010000000000000000000000000000000000000000000000000000000000000000000000000000000
000000001000000111000000000000000000000100000000000000100000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000000000000000000000000000000000000000000000000001000000000000000000
000000000000000000000000000000000000000000000000100000000000000001000000000000000
000000000000011000000000000000000000000000000000000000000000000000000000000000000
000000100000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000100000000000000000000000000000000000000000000000000000001000100000000000000
000000000000100000000000000000000000000000000000000000000000000001000000000000000
000000000000011000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000000000000000000000000000000000000000001000000001000000000000000000
000000000000000000000000000000000000000000000000000000000000000001000000000000000
000000000000011000000000000000000000000000000000000000000000000000000000000000000
001000010000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000100000000000000000000000000000000000000000000100000000000000000000
000000000000000000000000000000000000000000000000100000000000000000000000000000000
000000000011100000000000000000000000000000000000000010000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000100000000000000000000000000000001000000000000100000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
010000000011100000100000000000000000000000000000000100000000000000000000000000000
000000000000010000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000100000000000000000000000000000001000000000000100000000000000000100
000000000000000000000000000000000000000000000000000000000000000000000000000000000
010000000011100000100000000000000000000000000000000100000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000000000000000000000000000000000000000000011000010001000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
001000001000000000100000000000000000110100000000000101000100010000000000000000100
000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000
000000000000000000000000000000000000010000000001000000000000100000000000000000000
001000000000000000000000000001000000000000000000000000000000000000000000000000000
010000000111100000001000000001000000000000000000001000010000000000000000000000000
000000000000010000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000010000000000000000000000000000000000000000000000010
0000000000000000
m = 462, n = 421, 189694 zeros, 4808 ones, 194502 elements
first 10 singular values:
0, 0, nan, nan, nan, nan, nan, nan, nan, nan
(phit1 ~/pj)
--------------------------------------------------------------

Kind regards,
  Bruno Daniel






    _______________________________________________________

File Attachments:


-------------------------------------------------------
Date: Fri 29 Jan 2010 10:04:47 PM GMT  Name: svdbug.c  Size: 1kB   By: bjg

<http://savannah.gnu.org/bugs/download.php?file_id=19564>

    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?28767>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/