bugGNU Octave - Bugs: bug #59488, eigs error in dsaupd for zero...

 
 

bug #59488: eigs error in dsaupd for zero matrix

Submitted by:  JD <jdbancal>
Submitted on:  Thu 19 Nov 2020 09:36:17 PM UTC  
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Regression
Status:  Works For Me Assigned to:  None
Originator Name:  Open/Closed:  Open
Release:  6.0.92 Operating System:  GNU/Linux

Add a New Comment (Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

( Jump to the original submission)

Fri 20 Nov 2020 05:17:15 PM UTC, comment #12: 

Thanks for testing.  4.4.1 and 5.2.0 are the same thing, just with a nicer decoding of the message.  I checked the code in libinterp and liboctave and it doesn't appear that Octave changes the vector in any way.  This means that the error message from 4.4.1 ("error -9 in dsaupd") is correct, but attributing it to an empty starting vector is incorrect in this case.

There is a dependence on the version of the external ARPACK library.  In my case, I am using 3.5.0+real-2 and it just works.

It seems like ARPACK may be getting more brittle in the types of input matrices it can accept.  That is annoying, but means we at Octave probably have to program defensively around it.  I would suggest the nnz() test in _eigs_.cc.

I would probably put this change on the development branch rather than the about-to-be-released stable branch.  It is a corner case, and there is an easy workaround.  If Octave throws this error then the programmer already knows that all eigenvalues are zero and there is no need for any computation.

Rik <rik5>
Project Administrator
Fri 20 Nov 2020 04:51:41 PM UTC, comment #11: 

With 4.2.2:

octave:1> opts.v0 = [1;0;1];
octave:2> eigs (zeros(3), 1, 'lm', opts)
ans = 0

With 4.4.1:

octave:1> opts.v0 = [1;0;1];
octave:2> eigs (zeros(3), 1, 'lm', opts)
error: __eigs__: eigs: error -9 in dsaupd
error: called from
    eigs at line 285 column 18

With 5.2.0:

octave:1> opts.v0 = [1;0;1];
octave:2> eigs (zeros(3), 1, 'lm', opts)
error: __eigs__: eigs: error in dsaupd: Starting vector is zero
This should not happen.  Please, see https://www.gnu.org/software/octave/bugs.html, and file a bug report
error: called from
    eigs at line 285 column 20

All on the same system with the same version of ARPACK installed from the Debian package:

devnull:810> dpkg -l '*arpack*'
Desired=Unknown/Install/Remove/Purge/Hold
| Status=Not/Inst/Conf-files/Unpacked/halF-conf/Half-inst/trig-aWait/Trig-pend
|/ Err?=(none)/Reinst-required (Status,Err: uppercase=bad)
||/ Name                 Version      Architecture Description
+++-====================-============-============-============================================================================
ii  libarpack2:amd64     3.7.0-3      amd64        Fortran77 subroutines to solve large scale eigenvalue problems
ii  libarpack2-dev:amd64 3.7.0-3      amd64        Fortran77 subroutines to solve large scale eigenvalue problems (development)

Also, for the one that succeeds, eig is called instead of _eigs_, so ARPACK isn't even used in that case.

John W. Eaton <jwe>
Project Administrator
Fri 20 Nov 2020 04:48:55 PM UTC, comment #10: 

Here is what I get for the code you mentioned:

octave:1> opts.v0 = [1;0;1];
octave:2> eigs (zeros(3), 1, 'lm', opts)
error: eigs: error in dsaupd: Starting vector is zero
This should not happen.  Please, see https://www.gnu.org/software/octave/bugs.html, and file a bug report
error: called from
    eigs at line 298 column 20

JD <jdbancal>
Fri 20 Nov 2020 04:31:17 PM UTC, comment #9: 

Just to try and disambiguate whether there is a problem with Octave or a problem upstream in the way ARPACK now generates its starting vector could someone try

opts.v0 = [1;0;1];
eigs (zeros(3), 1, 'lm', opts)

Rik <rik5>
Project Administrator
Fri 20 Nov 2020 04:21:44 PM UTC, comment #8: 

The issue of interacting with Fortran and BLAS/LAPACK libraries when the matrices contain non-finite values is a large one.  If we want to tackle that, we should develop a strategy first.  Right now we have patched things up in various ways as problems have arisen (such as in inverse routines).  The danger here is that the code base grows substantially (more maintenance burden) to contain a lot of if/elseif checking for specific conditions.  And we may forget to port each case over to a newly discovered example of a bad C++/Fortran interaction.

I'm not even particularly in love with the "solution" we have for the inverse of a matrix.  Our current approach is to calculate the condition number and check it for Inf/NaN before calling the main inversion routines.  This is a performance hit if all the user wanted was the inverse.  If try/catch blocks are cheap in C++ maybe we should just be wrapping the external library calls in them and using the catch block to translate the undecipherable message "parameter number 4 is invalid" to something nice like "non-finite value encountered in matrix".

For empty matrices, it probably makes sense to have a shortcut path in the code.  Besides generating the correct result, it will be a performance win as well.  eigs() is designed to be called on sparse matrices and the nnz() function is cheap on sparse matrices.  So

if (nnz () == 0)
...
endif

Rik <rik5>
Project Administrator
Fri 20 Nov 2020 03:13:43 PM UTC, comment #7: 

In addition to matrices of zeros, do we also need to have special cases for matrices that contain Inf or NaN?  For example, I see the following:

octave> x = sprand (20, 20, 0.1); x(10) = NaN; eigs (x, 1)
DLASCL: parameter number 4 is invalid
error: Fortran procedure terminated by call to XERBLA
error: called from
    eigs at line 285 column 20

octave> x = sprand (20, 20, 0.1); x(10) = Inf; eigs (x, 1)
error: __eigs__: eigs: error in dneupd: DNAUPD did not find any eigenvalues to sufficient accuracy.
This should not happen.  Please, see https://www.gnu.org/software/octave/bugs.html, and file a bug report
error: called from
    eigs at line 285 column 20

John W. Eaton <jwe>
Project Administrator
Fri 20 Nov 2020 03:10:11 PM UTC, comment #6: 

It works for me with 4.2.2 but fails with 4.4.1, so the bug is not new.

Bisecting for old versions of Octave on a current system doesn't work very well because of changes in tools and things, so I tried to just run the debugger on my 4.2.2 and 4.4.1 builds and I found that dsaupd is not called for this problem with version 4.2.2.  So I stepped through eigs.m in both versions and found that in 4.2.2, _eigs_ was not called.

Next, I extracted all versions of eigs.m from the repo since the time around version 4.2.2 and tried the example with each one until I found the first that failed:

## list of revisions with changes to eigs.m
hg log scripts/*/eigs.m --no-merges --template "{rev}\n" > /tmp/revs

## extract eigs.m from each of those revisions
mkdir /tmp/eigs
for rev in $(cat /tmp/revs) ; do hg update -r $rev ; cp scripts/sparse/eigs.m /tmp/eigs/$rev-eigs.m ; done

## test each until we find a failure
## in this case, *-eigs.m is sorts correctly
## in chronological order so the first that fails
## should be the first bad revision
cd /tmp/eigs
for f in *-eigs.m ; do echo $f ; cp $f eigs.m ; octave --eval "eigs (zeros (3), 1)"; rm -f eigs.m ; done

The first revision that showed the error for me was

# HG changeset patch
# User Marco Caliari <marco.caliari@univr.it>
# Date 1524560299 -7200
#      Tue Apr 24 10:58:19 2018 +0200
# Branch stable
# Node ID 1a632692a58e3204d8597404cdee8d7200deac88
# Parent  4c98a9e5ce2557e64256a8f02d857157187a80a3
Use eig in eigs when p is equal to matrix dimension n (bug #53719)

However, even the previous version that doesn't fail for zeros(3) will fail for matrices of zeros that are larger than the limit eigs may also use to determine when to call eig instead of _eigs_.  So we may just need to have a special case for matrices that are all zeros.

John W. Eaton <jwe>
Project Administrator
Fri 20 Nov 2020 12:48:53 PM UTC, comment #5: 

My version of libarpack2 and libarpack2-dev (on ubuntu 20.04.1) is 3.7.0-3

JD <jdbancal>
Fri 20 Nov 2020 12:09:24 PM UTC, comment #4: 

Maybe a different version of ARPACK?

Markus M├╝tzel <mmuetzel>
Project Member
Fri 20 Nov 2020 11:38:28 AM UTC, comment #3: 

I get the same error with dev Octave on Mageia 7.

Philip Nienhuis <philipnienhuis>
Project Member
Fri 20 Nov 2020 10:46:47 AM UTC, comment #2: 

Thanks for checking!
I just tried on the latest hg, namely

changeset:   29079:c43e748de388
bookmark:    @
tag:         tip
parent:      29077:b85f72ec844a
parent:      29078:90a8e6fd85c3
user:        John W. Eaton <jwe@octave.org>
date:        Thu Nov 19 18:59:10 2020 -0500
summary:     maint: merge stable to default.

and I obtain the same error. The only configuration options I used are "--without-qt --without-opengl". Any idea why we see obtain results?

JD <jdbancal>
Thu 19 Nov 2020 10:49:19 PM UTC, comment #1: 

Could you try again with a recent tip from the stable branch of Octave?  I am at cset

parent: 29075:dc3ee9616267 tip
eigs.m: Fix error message when second input is empty matrix (bug #59486)
branch: stable

and this code works fine.

Rik <rik5>
Project Administrator
Thu 19 Nov 2020 09:36:17 PM UTC, original submission:  

octave:1> eigs(zeros(3), 1)
error: eigs: error in dsaupd: Starting vector is zero
This should not happen.  Please, see https://www.gnu.org/software/octave/bugs.html, and file a bug report
error: called from
    eigs at line 298 column 20

Used to work on Octave 4.2

JD <jdbancal>

 

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

Attach Files:
   
   
Comment:
   

No files currently attached

 

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 mmuetzel (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by jdbancal (Submitted the item)
  •  

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

     

     

    Follows 1 latest change.

    Date Changed by Updated Field Previous Value => Replaced by
    2020-11-19 rik5 StatusNone => Works For Me

    Back to the top


    Powered by Savane 3.5