bugGNU Octave - Bugs: bug #49523, eigs (R2013b and later) uses NaNs,...

 
 

bug #49523: eigs (R2013b and later) uses NaNs, not zeros, for unconverged Ritz values

Submitted by:  None
Submitted on:  Thu 03 Nov 2016 01:17:22 PM UTC  
 
Category: Octave FunctionSeverity: 3 - Normal
Priority: 5 - NormalItem Group: Matlab Compatibility
Status: FixedAssigned to: None
Originator Name: Tim MitchellOriginator Email: -unavailable-
Open/Closed: ClosedRelease: dev
Operating System: Any

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 Jump to the original submission)

Fri 01 Dec 2017 07:10:06 PM UTC, comment #54:

I pushed the fix of using isreal() to determine the expected value for the specific failing BIST test. See http://hg.savannah.gnu.org/hgweb/octave/rev/1b17a5d82e3f.

Hopefully the builbots will be happy now.

Rik <rik5>
Project Administrator
Fri 01 Dec 2017 07:22:56 AM UTC, comment #53:

I agree with Rik's proposal in comment #52: even if we can find a magic combination of v0, tol, maxit working with the existing combinations of blas/lapack/arpack/whatever, it is not sure it will work for any future combination. Here the important thing is to get either NaN or Nan+1iNaN (and not NaN+1i0).

To comment #37: -14 in XYeupd is "did not find any eigenvalues to sufficient accuracy". It is just non-convergence with the given parameters and it could be expected in our cases since the maximum number of iterations was strongly decreased. We could add a string to the error message like: "try to increase the tolerance, the maximum number of iterations or the number of requested eigenvalues". The third hint seems strange, but arpack uses internally Krylov spaces (spaces in which to look for eigenpairs) of default dimension opts.p = 2 * number_of_requested_eigenvalues.

Marco Caliari <caliari>
Project Member
Fri 01 Dec 2017 05:36:15 AM UTC, comment #52:

Maybe Carlo has an idea on how to relax this. I found that by increasing the maximum number of iterations, opts.maxit, ARPACK would converge and find all the eigenvalues. That is not what is needed here.

I restored maxit to 1, and then started changing the initialization vector v0. For some vectors there was rapid convergence and I got the three eigenvalues listed in comment #50 (although they were all real, no imaginary parts).

Maybe there is a specific initial vector which is always far enough from a solution that the answer from comment #46 is obtained.

Otherwise, we could work around it by coding the BIST test like this

Rik <rik5>
Project Administrator
Fri 01 Dec 2017 05:00:41 AM UTC, comment #51:

The last representation in comment #50 is the best. The imaginary component is very small, but non-zero, so Octave can't discard it.

Rik <rik5>
Project Administrator
Fri 01 Dec 2017 01:51:17 AM UTC, comment #50:

This is probably the best (format long e):

d =

1.69714256958507e-02 + 2.37387114336600e-17i
1.69714256958507e-02 - 2.37387114336600e-17i
6.79928617498172e-02 + 0.00000000000000e+00i
NaN + NaNi

Dmitri A. Sergatskov <dasergatskov>
Fri 01 Dec 2017 01:49:02 AM UTC, comment #49:

This is with "format bit"

Dmitri A. Sergatskov <dasergatskov>
Fri 01 Dec 2017 01:45:05 AM UTC, comment #48:

octave:8> d = eigs (A, 4, "sm", opts)
warning: eigs: Only 3 of the 4 requested eigenvalues converged
warning: called from
eigs at line 265 column 18
d =

0.0169714256958507 + 0.0000000000000000i
0.0169714256958507 - 0.0000000000000000i
0.0679928617498172 + 0.0000000000000000i
NaN + NaNi

Dmitri A. Sergatskov <dasergatskov>
Fri 01 Dec 2017 01:43:47 AM UTC, comment #47:

Well that is the problem. Could you redo the test with openblas (comment #45) and use 'format long' before displaying the variable d? It looks like Octave couldn't narrow the variable from complex to real because there was a very small imaginary value.

Rik <rik5>
Project Administrator
Fri 01 Dec 2017 01:18:45 AM UTC, comment #46:

With atlas it is the same as on i7 (with any blas):

LD_PRELOAD=/usr/lib64/atlas/libtatlas.so ./run-octave -f -q --no-gui
octave:1> test eigs
PASSES 179 out of 179 tests
octave:2> A = toeplitz ([0, 1, zeros(1, 8)], [0, -1, zeros(1, 8)]);
octave:3> A(1, 1) = 1;
octave:4> A = kron (A, eye (10)) + kron (eye (10), A);
octave:5> opts.v0 = (1:100)';
octave:6> opts.maxit = 1;
octave:7> d = eigs (A, 4, "sm", opts)
warning: eigs: Only 2 of the 4 requested eigenvalues converged
warning: called from
eigs at line 265 column 18
d =

0.016971
0.016971
NaN
NaN

Dmitri A. Sergatskov <dasergatskov>
Fri 01 Dec 2017 01:13:07 AM UTC, comment #45:

With openblas:

octave:2> A = toeplitz ([0, 1, zeros(1, 8)], [0, -1, zeros(1, 8)]);
octave:3> A(1, 1) = 1;
octave:4> A = kron (A, eye (10)) + kron (eye (10), A);
octave:5> opts.v0 = (1:100)';
octave:6> opts.maxit = 1;
octave:7> d = eigs (A, 4, "sm", opts)
warning: eigs: Only 3 of the 4 requested eigenvalues converged
warning: called from
eigs at line 265 column 18
d =

0.016971 + 0.000000i
0.016971 - 0.000000i
0.067993 + 0.000000i
NaN + NaNi

Dmitri A. Sergatskov <dasergatskov>
Fri 01 Dec 2017 01:09:42 AM UTC, comment #44:

If you could run interactively and get the value for d(4) that would help.

In the message from assert, the 'O' just stands for observed and the 'E' stands for expected. The message 'complex != real' explains what is going on. The two asserts are

and it is failing on the first one with

The expected is just NaN, a real number. The observed was a complex value. You can see that the second test is designed to ensure that d(4) is a real value by checking that the imaginary part was equal to 0. Perhaps a better test would be to check directly with isreal.

Rik <rik5>
Project Administrator
Fri 01 Dec 2017 12:41:03 AM UTC, comment #43:

I meant to say openblas is 0.2.20 (the latest stable release).

Dmitri.

Dmitri A. Sergatskov <dasergatskov>
Fri 01 Dec 2017 12:38:17 AM UTC, comment #42:

This erros happens with either openblas-0.2.0 (custom compiled for this machine) or reference blas ("standard" Fedora package).
No problem with atlas.
The build bot is AMD FX-8350; my workstation is i7-2600 and does
not show this error (openblas is different , but reference
blas is exactly the same).

I am not sure what to do about it. I can set buildbot to
use atlas.

On related note -- the error seems strange to me. Should not
observed and expected values be numbers rather than letters?

Dmitri.
--

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 11:28:53 PM UTC, comment #41:

A-a-a-nd I see that fedora buildbot still fails with:
http://buildbot.octave.org:8010/builders/clang-fedora/builds/1015/steps/test/logs/stdio

          • testif HAVE_ARPACK

A = toeplitz ([0, 1, zeros(1, 8)], [0, -1, zeros(1, 8)]);
A(1, 1) = 1;
A = kron (A, eye (10)) + kron (eye (10), A);
opts.v0 = (1:100)';
opts.maxit = 1;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(4), NaN);
assert (imag (d(4)), 0);
!!!!! test failed
ASSERT errors for: assert (d (4),NaN)

Location | Observed | Expected | Reason
() O E complex != real

it could be related to openblas library it uses. I will investigate.

Dmitri.
--

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 10:54:36 PM UTC, comment #40:

I pushed Marco's relaxed BIST test changeset here: http://hg.savannah.gnu.org/hgweb/octave/rev/30bbef7bf919.

Closing report again.

Rik <rik5>
Project Administrator
Thu 30 Nov 2017 08:59:08 PM UTC, comment #39:

I fixed the incorrect error reporting here http://hg.savannah.gnu.org/hgweb/octave/rev/30240b146061.

I will review eigs3.diff and push that shortly as well.

Rik <rik5>
Project Administrator
Thu 30 Nov 2017 08:46:43 PM UTC, comment #38:

Yes, with eigs3.diff eigs passes all tests.

Dmitri.
--

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 08:45:52 PM UTC, comment #37:

@Marco: I was just making a guess. Since eigs is returning NaN for sum of the eigenvalues I assume that some of the eigenvectors must also have NaN in them and that maybe this was upsetting zneupd. I suppose the real way to figure out the problem is to go look at the Fortran code and figure out what the error code of -14 means.

I will take care of changing the error messages to match the Fortran call involved.

Rik <rik5>
Project Administrator
Thu 30 Nov 2017 08:33:06 PM UTC, comment #36:

@Rik: not sure I understood your comments. znaupd quits with a value of ido in {-1,1,2,99}. 99 means "done" and info is checked. Of course the error message is wrong. Which vector could have norm NaN before zneupd?

Marco Caliari <caliari>
Project Member
Thu 30 Nov 2017 07:24:11 PM UTC, comment #35:

@Dmitri: Can you test Marco's patch on your failing system? I don't want to commit it to Mercurial until it is clear that it solves the problem.

@Marco: Are we sure that there is no way for the code in eigs-base.cc to exit from znaupd with an error code that is not checked before calling zneupd? Brief documentation for zneupd is at http://www.caam.rice.edu/software/ARPACK/UG/node44.html#2332.

According to that, "On the final return from znaupd (indicated by ido = 99), the error flag info must be checked. If info = 0, then no fatal errors have occurred and it is time to post-process using zneupd to get eigenvalues of the original problem and the corresponding eigenvectors if desired."

The error reporting looks wrong too.

This looks like it should reference znaupd, not dsaupd.

Also, what happens if the norm of one of the vectors is NaN and it attempts to normalize the lengths in zneupd? That might throw an error.

Rik <rik5>
Project Administrator
Thu 30 Nov 2017 07:07:30 PM UTC, comment #34:

Please remove eigs2.diff and apply eigs3.diff.

@Rik: in fact, I relaxed the tests.

(file #42536)

Marco Caliari <caliari>
Project Member
Thu 30 Nov 2017 06:40:12 PM UTC, comment #33:

OK, so with the patch is only one test is failing, although the
errors are different depending on blas version.

Ref. blas:

LD_PRELOAD=/usr/lib64/libblas.so ./run-octave -f -q --no-gui
octave:1> test ("eigs", "normal", "results_blas_1.txt")

>>>>> /d2/home/dima/src/octave/gcc_def/../scripts/sparse/eigs.m

PASSES 178 out of 179 tests
octave:2> test ("eigs", "normal", "results_blas_2.txt")

>>>>> /d2/home/dima/src/octave/gcc_def/../scripts/sparse/eigs.m

PASSES 178 out of 179 tests
octave:3> test eigs

          • testif HAVE_ARPACK

A = 1i * magic (100) / 100 + eye (100);
opts.v0 = (1:100)';
opts.maxit = 1;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(3:4), [NaN+1iNaN; NaN+1iNaN]);
!!!!! test failed
eigs: error -14 in zneupd

Same with openblas:

LD_PRELOAD=/usr/lib64/libopenblaso.so ./run-octave -f -q --no-gui
octave:1> test ("eigs", "normal", "results_openblas_1.txt")

>>>>> /d2/home/dima/src/octave/gcc_def/../scripts/sparse/eigs.m

PASSES 178 out of 179 tests
octave:2> test ("eigs", "normal", "results_openblas_2.txt")

>>>>> /d2/home/dima/src/octave/gcc_def/../scripts/sparse/eigs.m

PASSES 178 out of 179 tests
octave:3> test eigs

          • testif HAVE_ARPACK

A = 1i * magic (100) / 100 + eye (100);
opts.v0 = (1:100)';
opts.maxit = 1;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(3:4), [NaN+1iNaN; NaN+1iNaN]);
!!!!! test failed
eigs: error -14 in zneupd

Same test, different error with atlas:

LD_PRELOAD=/usr/lib64/atlas/libtatlas.so ./run-octave -f -q --no-gui
octave:1> test eigs

          • testif HAVE_ARPACK

A = magic (100) / 100 + eye (100);
opts.v0 = (1:100)';
opts.maxit = 10;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(4), NaN+1i*NaN);
!!!!! test failed
ASSERT errors for: assert (d (4),NaN + 1i * NaN)

Location | Observed | Expected | Reason
() 1-1.219e-14i NaN+NaNi 'NaN' mismatch
octave:2> test ("eigs", "normal", "results_atlas_1.txt")

>>>>> /d2/home/dima/src/octave/gcc_def/../scripts/sparse/eigs.m

PASSES 178 out of 179 tests

Dmitri.
--

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 05:07:23 PM UTC, comment #32:

@Marco: Is it possible to pick the runtime parameters such as the initialization vector v0 and the maximum number of iterations is such that at least some of the returned values are NaN? I know you have initialized both of those items, but I wonder whether we should attempt to optimize those, or maybe we should relax the tests. For example, this is a failing test

What if we relaxed the conditions for the assert so that it was

This would still verify that unconverged Ritz values are returned as NaN, but wouldn't require that exactly 3 Ritz values be unconverged.

Rik <rik5>
Project Administrator
Thu 30 Nov 2017 04:54:24 PM UTC, comment #31:

@Dmitri: You should be able to run all of the tests and dump the log into a file with

This will stick the logfile into the filename "results.txt".

Rik <rik5>
Project Administrator
Thu 30 Nov 2017 04:14:15 PM UTC, comment #30:

Before your patch 6 eigs test were failing on Fedora and 5 on debian (jwe buildbot). You can see logs at

http://buildbot.octave.org:8010/waterfall

There are still multiple failures (some appear vary depending on
blas version used), but I could not figure out how to get
"test" to report all the failures in interactive mode.

Dmitri.
--

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 04:02:31 PM UTC, comment #29:

@Dmitri: Any other test not working?

Marco Caliari <caliari>
Project Member
Thu 30 Nov 2017 03:59:33 PM UTC, comment #28:

octave:1> A = 1i * magic (100) + eye (100);
octave:2> opts.v0 = (1:100)';
octave:3> opts.maxit = 7;
octave:4> d = eigs (A, 10, "sm", opts)
warning: eigs: Only 3 of the 10 requested eigenvalues converged
warning: called from
eigs at line 265 column 18
d =

1.00000 - 0.00000i
1.00000 + 0.00000i
1.00000 - 0.00000i
NaN + NaNi
NaN + NaNi
NaN + NaNi
NaN + NaNi
NaN + NaNi
NaN + NaNi
NaN + NaNi

(Same with either ref blas, atlas, openblas)

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 03:52:26 PM UTC, comment #27:

@Dmitri: can you tell me the result of

But is it worth/possible to try to construct a set of examples working (that is, with eigs failing) with any possible combination of blas/lapack/arpack libraries?

Marco Caliari <caliari>
Project Member
Thu 30 Nov 2017 02:57:02 PM UTC, comment #26:

The patch changed some errors, e.g. with the reference blas
I am getting this now:

test eigs

          • testif HAVE_ARPACK

A = 1i * magic (100) / 100 + eye (100);
opts.v0 = (1:100)';
opts.maxit = 1;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(3:4), [NaN+1iNaN; NaN+1iNaN]);
!!!!! test failed
eigs: error -14 in zneupd

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 08:06:58 AM UTC, comment #25:

@Rik, Dmitri: are the results really intermittent or 100% reproducible?
In the first case, no idea, since in the tests I initialized the initial vector which should be the only otherwise random choice.
I tried standard blas, atlas and openblas on a Ubuntu systems, I was able to reproduce some failures. Convergence clearly depends on the blas/lapack libraries (and the very stringent tolerance, which is eps by default). I relaxed the tests and now they succeed with all my libraries. Please try the patch.

(file #42528)

Marco Caliari <caliari>
Project Member
Thu 30 Nov 2017 05:53:33 AM UTC, comment #24:

I tried arpack 3.4 (from Fedora 26) arpack 3.6 (master source from https://github.com/opencollab/arpack-ng).

Same results.

Dmitri.
--

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 05:33:43 AM UTC, comment #23:

I also get different test failures depending on the blas version used:

LD_PRELOAD=/usr/lib64/libopenblas.so.0 ./run-octave -f -q --no-gui
octave:1> test eigs

          • testif HAVE_ARPACK

A = toeplitz ([-2, 1, zeros(1, 8)]);
A = kron (A, eye (10)) + kron (eye (10), A);
opts.v0 = (1:100)';
opts.maxit = 3;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "lm", opts);
assert (d(2:4), [NaN; NaN; NaN]);
!!!!! test failed
ASSERT errors for: assert (d (2:4),[NaN; NaN; NaN])

Location | Observed | Expected | Reason
(1) -7.365 NaN 'NaN' mismatch

LD_PRELOAD=/usr/lib64/atlas/libsatlas.so ./run-octave -f -q --no-gui
octave:1> test eigs

          • testif HAVE_ARPACK

A = magic (100) / 100 + eye (100);
opts.v0 = (1:100)';
opts.maxit = 1;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(3:4), [NaN+1iNaN; NaN+1iNaN]);
!!!!! test failed
ASSERT errors for: assert (d (3:4),[NaN + 1i * NaN; NaN + 1i * NaN])

Location | Observed | Expected | Reason
(1) 1+1.783e-14i NaN+NaNi 'NaN' mismatch

LD_PRELOAD=/usr/lib64/libblas.so.3 ./run-octave -f -q --no-gui
octave:1> test eigs

          • testif HAVE_ARPACK

A = toeplitz ([0, 1, zeros(1, 8)], [0, -1, zeros(1, 8)]);
A(1, 1) = 1;
A = kron (A, eye (10)) + kron (eye (10), A);
opts.v0 = (1:100)';
opts.maxit = 1;
warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
d = eigs (A, 4, "sm", opts);
assert (d(3:4), [NaN; NaN]);
assert (imag (d(3:4)), [0; 0]);
!!!!! test failed
ASSERT errors for: assert (d (3:4),[NaN; NaN])

Location | Observed | Expected | Reason
(1) 0.067993 NaN 'NaN' mismatch

Those are 100% reproducible.
I use
arpack-devel-3.5.0-4.fc27.x86_64
arpack-3.5.0-4.fc27.x86_64
(same as fedora buildbots)

Dmitri.

Dmitri A. Sergatskov <dasergatskov>
Thu 30 Nov 2017 05:22:15 AM UTC, comment #22:

Have to re-open the bug temporarily. The BIST tests are creating intermittent failures on the automated buildbot testing system. Dmitri first reported this. See log files at http://buildbot.octave.org:8010/waterfall.

Some sample errors are:

This is stressing the ARPACK library. And maybe it indicates the ARPACK library used by the buildbots is bad.

Rik <rik5>
Project Administrator
Wed 29 Nov 2017 09:12:18 PM UTC, comment #21:

I checked in Marco's changeset on the development branch here (http://hg.savannah.gnu.org/hgweb/octave/rev/0c1417a8d9d0).

@Marco: You were right that it isn't a good idea to disable all warnings during the BIST tests. Instead, I created a new warning ID, "Octave:eigs:UnconvergedEigenvalues", and I temporarily disable just this one warning during the BIST tests. See cset http://hg.savannah.gnu.org/hgweb/octave/rev/d3d38a443df8.

Marking report as fixed and closing it.

Rik <rik5>
Project Administrator
Mon 27 Nov 2017 02:35:13 PM UTC, comment #20:

Great, thank you so much for addressing the issue I reported!

Anonymous
Mon 27 Nov 2017 02:04:10 PM UTC, comment #19:

I changed my mind and implemented the NaN/NaN+1i*NaN version. Unconverged eigenvalues and eigenvectors are always in the last positions. I used warning("off","all") in the tests, but maybe there is a better way to suppress the specific warning abount unconverged values.

(file #42497)

Marco Caliari <caliari>
Project Member
Fri 24 Nov 2017 05:34:09 PM UTC, comment #18:

@Rik: I'm going to implement the NaN + i*0 version for both eigenvalues and eigenvectors. Now that eigs in octave is going to be in a good shape, I discovered that matlab switched to "A Krylov Schur algorithm for large eigenproblems". I found only a simple implementation around (https://github.com/dingxiong/KrylovSchur), so I think that matlab implemented it from scratch.

Marco Caliari <caliari>
Project Member
Thu 23 Nov 2017 06:30:30 PM UTC, comment #17:

@Marco: I had originally proposed setting just the real part to NaN. But, I agree with you that it looks nicer if NaN + NaNi is returned for complex unconverged eigenvalues.

I don't think it matters from a performance standpoint if we post-process the list of eigenvalues because the number is going to be small, say 50 eigenvalues, relative to the size of the sparse matrix itself which could have thousands to millions of rows or columns.

I'll accept whatever you feel is best.

Rik <rik5>
Project Administrator
Thu 23 Nov 2017 03:57:45 PM UTC, comment #16:

@Rik: Matlab now gives NaN+0.0000i, if there are other complex eigenvalues. Is it what you proposed in comment #12?

Marco Caliari <caliari>
Project Member
Thu 23 Nov 2017 03:42:04 PM UTC, comment #15:

@Marco: Any progress? This bug seems easy enough to knock off now that the behavior (return NaN, or Nan+Nani, for unconverged Ritz values and full NaN vectors for corresponding eigenvectors) has been agreed upon.

Rik <rik5>
Project Administrator
Tue 14 Nov 2017 05:53:18 AM UTC, comment #14:

@Marco: I just checked in 3 of your csets improving the behavior of eigs(). This bug report seems like a natural one to finish off as well. I think you can code Octave to return Nan+NaNi for unconverged values.

Rik <rik5>
Project Administrator
Mon 17 Apr 2017 03:50:10 AM UTC, comment #13:

@Marco: I'm going to change my advice in comment #12. It seems to me that eigs is generally only used when you want a few eigenvalues from a large sparse matrix. For a low number, say N < 50, it really isn't much of a performance hit to iterate over the entire return list to see if there are any true complex eigenvalues or if there are only complex eigenvalues where the imaginary part is zero. Thus, I think you should go ahead and return NaN if all the values are real or Nan+NaNi if there are complex values.

Rik <rik5>
Project Administrator
Wed 12 Apr 2017 04:31:01 PM UTC, comment #12:

The unconverged Ritz values represent an exceptional case. It is important to warn about them, but the exact form of the warning is not that critical. My suggestion was to set only the real portion to NaN. After that, I would rely on Octave's automatic narrowing to change a complex matrix, with all imaginary parts equal to zero, to a real matrix. If the matrix can't be narrowed because there are non-zero complex values then you do get

but that is okay for me since the NaN stands out as the first thing in the entry, and isnan returns true if either the real or imaginary part is NaN.

Rik <rik5>
Project Administrator
Wed 12 Apr 2017 07:37:11 AM UTC, comment #11:

@Rik: If I set them to Complex (NaN, 0.), the result will be somethinkg like (there are other complex eigenvalues)

or (all the remaining eigenvalues are real)

If I set them to Complex (NaN, NaN), the result will be somethinkg like

or

If we want

or

we have to check the imaginary part of the converged eigenvalues. I mean, it is not possible to set "without the complex [imaginary] part".

Marco Caliari <caliari>
Project Member
Tue 11 Apr 2017 03:15:01 PM UTC, comment #10:

I would just set the unconverged values to NaN without the complex part. The real part being NaN is enough to indicate that there is something wrong with this particular value.

Rik <rik5>
Project Administrator
Tue 11 Apr 2017 07:52:57 AM UTC, comment #9:

In the nonsymmetric case, the nonconverged Ritz values are set to NaN+NaNi in Matlab. You can try

The converged values are complex with imaginary part zero. If you omit opts, you see that all the converged eigenvalues are double. If we follow Matlab and comment #5, in case of nonconverged Ritz values we should set them and the corresponding eigenvectors to NaN+NaNi. An alternative would be: if the converged Ritz values are double, set the nonconverged and the corresponding eigenvectors to (double) NaN. This would save some memory.

Marco Caliari <caliari>
Project Member
Mon 10 Apr 2017 05:28:53 PM UTC, comment #8:

@Rik: Interesting. On R2016b, on my hardware, I consistently see that allocating NaNs is about 10 times slower than just zeros. Allocating ones is also about 10 times slower than just zeros on my hardware.

In any case, I don't feel strongly either way on how to handle the values of the "non-converged" eigenvectors and I do generally agree with trying to prevent garbage in, garbage out. I was just pointing out that there was a potential downside to explicitly setting all the entries in the eigenvectors.

Anonymous
Mon 10 Apr 2017 04:46:22 PM UTC, comment #7:

I don't see any real difference in creating matrices with all zeros or all NaNs. See below

I don't think performance is an issue here. And I have always tried to code Octave to avoid Garbage In / Garbage Out phenomena. I believe many Octave users are professional scientists first, and coders second. I don't think they will necessarily realize that they need to code defensively and check their eigenvalues before going on to use the eigenvectors. For that reason, I would prefer to return NaNs for the eigenvectors to immediately inform the coder that something is amiss.

Rik <rik5>
Project Administrator
Mon 10 Apr 2017 04:22:19 PM UTC, comment #6:

Original submitter here.

@Marco: The MathWorks has their own private fork of the original ARPACK source with various changes of their own, which probably explains why you are seeing different behaviors regarding the eigenvectors of non-converged Ritz values in Matlab compared to ARPACK and arpack-ng.

@Marco, @Rik: an argument against using nans for eigenvecotors whose Ritz values did not converge is that returning a large number of nans will take longer than just returning zeros (or using whatever values are already set for these eigenvectors). On my laptop, setting up an array of 1 million nans takes about 10 times longer than 1 million zeros. Perhaps it's a negligible cost but I personally don't see that it's necessary to use nans. The user already has to check which, if any, eigenvalues are nans. At that point, it should be clear to the user which eigenvectors are good and which are to be ignored. That said, I could see a user using an eigenvector without checking if its Ritz value converged!

In my eigsPlus package, I avoided this issue altogether by only returning the converged eigenvalues (and eigenvectors). ;-) I believe I also suggested this to The MathWorks long ago but they didn't want to change that.

Anonymous
Mon 10 Apr 2017 04:01:24 PM UTC, comment #5:

@Marco: I agree with you. If one of the eigenvalues is NaN, it seems like the eigenvector for that eigenvalue should also be NaN.

Rik <rik5>
Project Administrator
Mon 10 Apr 2017 02:59:48 PM UTC, comment #4:

I checked that ARPACK gives the zero eigenvector for not converged Ritz values. I cannot understand the result in Matlab. For non-converged Ritz values I would return a NaN eigenvector. And I will propose this change to arpack-ng.

Marco Caliari <caliari>
Project Member
Fri 07 Apr 2017 12:29:47 PM UTC, comment #3:

For the same example below, if you call

in Matlab, you get D with the entry NaN, but V with four columns of numbers. Should we copy such a weird result?

Marco Caliari <caliari>
Project Member
Fri 07 Apr 2017 08:55:17 AM UTC, comment #2:

I created a reproducible example for this issue.

Marco Caliari <caliari>
Project Member
Thu 03 Nov 2016 03:38:05 PM UTC, comment #1:

This seems like a good idea. The Octave devs are working really hard right now to get the 4.2.0 release out, but after that this could be taken up on the development branch and made a part of the 4.4 release.

Rik <rik5>
Project Administrator
Thu 03 Nov 2016 01:17:22 PM UTC, original submission:

Since R2013b, Matlab's version of eigs uses NaNs as placeholders for the Ritz values, not zeros, which was the older convention and is what Octave appears to still be using (I tested 4.0.0).

The reason for the change to NaNs is under the older scheme using zeros, if not all the requested eigenvalues converged, it could be difficult to sort out the converged eigenvalues from the placeholder zeros, especially if 0 should be an eigenvalue of the matrix in question. Furthermore, if a user is interested in computing the spectral abscissa of a stable matrix, they may try calling
but this would misleading return 0 in the case when not all requested eigenvalues converge.

In fact, I am the person who pointed these issues out to The MathWorks in the first place and at the time, I made the suggestion for them to use NaNs instead, after I had created my own modified version of eigs doing just that. It would be great if Octave could adopt this convention as well.

Anonymous

 

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

Attach File(s):
   
   
Comment:
   

Attached Files
file #42536:  eigs3.diff added by caliari (3KiB - application/x-tex)
file #42528:  eigs2.diff added by caliari (2KiB - application/x-tex)
file #42497:  eigs.diff added by caliari (20KiB - application/x-tex)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -unavailable- added by dasergatskov (Posted a comment)
  • -unavailable- added by caliari (Posted a comment)
  • -unavailable- added by rik5 (Posted a comment)
  • -unavailable- added by None (Submitted the item)
  •  

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

    Only project members can vote.

     

    Please enter the title of George Orwell's famous dystopian book (it's a date):

     

     

    Follow 12 latest changes.

    Date Changed By Updated Field Previous Value => Replaced By
    Thu 30 Nov 2017 10:54:36 PM UTCrik5StatusIn Progress=>Fixed
      Open/ClosedOpen=>Closed
    Thu 30 Nov 2017 07:07:30 PM UTCcaliariAttached File-=>Added eigs3.diff, #42536
    Thu 30 Nov 2017 08:06:58 AM UTCcaliariAttached File-=>Added eigs2.diff, #42528
    Thu 30 Nov 2017 05:34:37 AM UTCrik5Open/ClosedClosed=>Open
    Thu 30 Nov 2017 05:22:15 AM UTCrik5StatusFixed=>In Progress
    Wed 29 Nov 2017 09:12:18 PM UTCrik5StatusPatch Submitted=>Fixed
      Open/ClosedOpen=>Closed
    Mon 27 Nov 2017 02:09:33 PM UTCcaliariStatusConfirmed=>Patch Submitted
    Mon 27 Nov 2017 02:04:10 PM UTCcaliariAttached File-=>Added eigs.diff, #42497
    Mon 10 Apr 2017 04:01:24 PM UTCrik5Release4.0.0=>dev
    Thu 03 Nov 2016 03:38:05 PM UTCrik5StatusNone=>Confirmed

    Back to the top


    Powered by Savane 3.1-cleanup1