arpack-ng
arpack-ng copied to clipboard
For zero input matrix, fails with "starting vector is zero"
Expected behavior
With ARPACK 3.5.0 this works:
$ octave
octave:1> m = zeros(4,4);
octave:2> eigs(m, 1)
ans = 0
Actual behavior
With ARPACK 3.7.0 and 3.6.0:
$ octave
octave:1> m = zeros(4,4);
octave:2> eigs(m,1)
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
Notes, remarks
This is probably from the same changes as #142: forcing the initial vector to the range of a zero matrix produces a zero initial vector, which fails. The *getv0
routines could perhaps check if resid
is a zero vector, and (in the non-generalized cases?), skip the forcing to the range?
Yes, the problem could be related with #79 (force the initial vector to be in the range). In #142 there is a further modification which should manage rank deficient matrices. This case is the most rank deficient matrix. I think it should be enough simply not to check whether the initial vector is zero (around line 324 of dnaup2.f and related). In fact, after the first try, the initial vector is no more forced to be in the range (this was the outcome of #142).
By the way, this issue should appear also with 3.5.0, by selecting the generalized problem with the zero input matrix and the identity matrix B. But you need plain arpack to check, not octave.
I verified that arpack 3.5.0 with a dummy generalized problem with the zero input matrix and the identity matrix B stops with "starting vector is zero".
Any ideas here? This is also breaking some tests when building the R igraph package with arpack 3.7.0.
-- 1. Error: cluster_leading_eigen works (@test_leading.eigenvector.community.R#
At arpack.c:944 : ARPACK error, Starting vector is zero
Backtrace:
1. igraph::cluster_leading_eigen(g)
I think I found a workaround at the time, give me a week or so.
This is my proposal
diff --git a/SRC/dgetv0.f b/SRC/dgetv0.f
index 1d6dc01..9635720 100644
--- a/SRC/dgetv0.f
+++ b/SRC/dgetv0.f
@@ -258,6 +258,8 @@ c %-----------------------------------------%
c | Back from computing OP*(initial-vector) |
c %-----------------------------------------%
c
+ 10 continue
+
if (first) go to 20
c
c %-----------------------------------------------%
@@ -309,6 +311,11 @@ c %---------------------------------------------%
c | Exit if this is the very first Arnoldi step |
c %---------------------------------------------%
c
+ if (rnorm == zero .and. .not.initv) then
+ call dlarnv (idist, iseed, n, workd(n + 1))
+ go to 10
+ end if
+
if (j .eq. 1) go to 50
c
c %----------------------------------------------------------------
that is, after the computation of the initial vector (in the range of the operator), if it is zero and was not provided by the user, generate a random vector as initial vector and continue. If acceptable, the four *getv0.f should be modified and a proper test added.
Thanks for the patch. It does not fix igraph, but I just realized it isn't failing in dgetv0
, but in dsaupd
.
Thanks, I will investigate why the problem is still there for the symmetric case.
This is the error before my patch
Error with _saupd, info = -9
Check documentation in _saupd
and this is the result after
Ritz values and relative residuals
----------------------------------
Col 1 Col 2
Row 1: 0.00000D+00 NaN
Row 2: 0.00000D+00 NaN
Row 3: 0.00000D+00 NaN
Row 4: 0.00000D+00 NaN
_SDRV1
======
Size of the matrix is 100
The number of Ritz values requested is 4
The number of Arnoldi vectors generated (NCV) is 10
What portion of the spectrum: SM
The number of converged Ritz values is 4
The number of Implicit Arnoldi update iterations taken is 1
The number of OP*x is 20
The convergence criterion is 1.1102230246251565E-016
I tested driver dsdrv1 with the zero input matrix. Of corse there are NANs because the driver computes relative errors. As you see, all the computed eigenvalues are zero. So, my patch fixes the problem for me.
Has a patch been merged for this?
You likely need to shift a bit as the eigne value is zero? Checkout out #408