On OS X 10.7, the following tests fail: 631 - gmres_2 (Failed) 632 - minres_1 (Failed) 634 - levenberg_marquardt (Failed) This is the output from running each test manually: -------- ./gmres_2 Initializing random number generator with seed 1387021094 Repeating each test 10 times Test check_sparse_square_solving(gmres_colmajor_diag) failed in /sw/build.build/eigen3-3.2.0-1/eigen-eigen-ffa86ffb5570/unsupported/test/../../test/sparse_solver.h (39) x.isApprox(refX,test_precision<Scalar>()) Stack: - check_sparse_square_solving(gmres_colmajor_diag) - test_gmres_T<std::complex<double> >() - gmres - Seed: 1387021094 Abort trap: 6 --------- ./minres_1 Initializing random number generator with seed 1387021144 Repeating each test 10 times Test check_sparse_square_solving(minres_colmajor_diag) failed in /sw/build.build/eigen3-3.2.0-1/eigen-eigen-ffa86ffb5570/unsupported/test/../../test/sparse_solver.h (39) x.isApprox(refX,test_precision<Scalar>()) Stack: - check_sparse_square_solving(minres_colmajor_diag) - test_minres_T<double>() - minres - Seed: 1387021144 Abort trap: 6 ---------- levenberg_marquardt Initializing random number generator with seed 1387021162 Repeating each test 10 times actual = 546 expected = 545 Test testNistMGH17() failed in /sw/build.build/eigen3-3.2.0-1/eigen-eigen-ffa86ffb5570/unsupported/test/levenberg_marquardt.cpp (954) test_is_equal(lm.njev(), 545) Stack: - testNistMGH17() - levenberg_marquardt - Seed: 1387021162 Abort trap: 6 ------------ I can provide the crash log that OS X generates if necessary. Interestingly, gmres_2 succeeds about 20% of the time, while the other two tests fail all the time.

I also encountered the minres_1 crash when running on Ubuntu with seed 1389235356.

I can confirm that levenberg_marquardt and minres_1 fail constantly on all my platforms (minres_2 actually does no tests at all). Also, gmres_2 fails very often, which might just be a too strict precision requirement. I added all authors of the test-cases to the cc-list, so maybe they can investigate this.

Hi, I contributed the MINRES solver. The test seems to fail for two reasons: 1) the MINRES solver does not check for zero rhs. 2) the "check_sparse_square_solving" is not appropriate for MINRES, because MINRES requires symmetric (possibly indefinite) problems. This said, MINRES still needs support for complex scalars. Is it possible to get write privileges to the Eigen repository so that I can commit the necessary changes?

Also, I'm not sure why two binaries (miners_1, and miners_2) are created. I'm not really sure what case each binary tests. Any help understanding this is appreciated.

(In reply to comment #3) > Is it possible to get write privileges to the Eigen repository so that I can > commit the necessary changes? The easiest solution is to either provide a patch http://eigen.tuxfamily.org/index.php?title=Mercurial#Generating_a_patch Or make a fork of Eigen and make pull-request: https://bitbucket.org/eigen/eigen/fork If you make regular contributions we can think about giving you write access. (In reply to comment #4) > Also, I'm not sure why two binaries (miners_1, and miners_2) are created. Our test-scripts search all unit-test source files for strings such as CALL_SUBTEST_1 and then compile this for every number which occurred. This is also the case, if the CALL_SUBTEST_2 is just part of a comment (not really often the case and no serious problem). Generally, the idea is to split test files into smaller sub-tests. Details are here: http://eigen.tuxfamily.org/index.php?title=Tests#Writing_unit_tests

Regarding levenberg_marquardt, I bisected which commit started the test to fail: https://bitbucket.org/eigen/eigen/commits/774113905339/ This was a bugfix to make stableNorm() more stable. I have not debugged yet what exactly goes wrong then.

Ok, I will send a patch then. I have a quick questions. Is there a test for symmetric indefinite matrices? Or symmetric positive semi-definite matrices? Something like "check_sparse_sind_solving" or "check_sparse_spsd_solving"?

Created attachment 431 [details] patch for bug 715: unit test fail Here is the proposed patch. Note that MINRES still lacks support for std::complex. Also, some additional optimization may still be possible.

(In reply to comment #8) > Created attachment 431 [details] [review] > patch for bug 715: unit test fail > > Here is the proposed patch. > > Note that MINRES still lacks support for std::complex. Also, some additional > optimization may still be possible. Notice that, even after the patch, the MINRES test fails from time to time, but it does not crash. The failure is due to the low tolerance required by the test.

(In reply to comment #8) > Here is the proposed patch. Patch accepted and pushed upstream: https://bitbucket.org/eigen/eigen/commits/08493d4 > Note that MINRES still lacks support for std::complex. Also, some additional > optimization may still be possible. Some simple optimizations are to use v1.swap(v2) when instead of v1=v2, if v2 is overwritten anyways in the next instruction. There is no need for .noalias() in line 122 (only if a matrix-vector product is involved), but it does not hurt. Regarding code style: Please avoid CAPSLOCK. We use it only for preprocessor macros and for compile-time assertions. If you want to mark problems or open questions inside your code, use FIXME or TODO (followed by normal cased description). For complex support, you might want to create a new bug-entry.

There is another perspective on the unit tests that fail: An important parameter for the convergence rate of Krylov methods is the condition number of the matrix. Although Krylov methods theoretically are known to provide the exact solution after a finite number of iterations in exact arithmetic, in finite precision the convergence rate often is very slow or the algorithm may even fail to converge due to loss of (bi-)orthogonality if the condition number is too large. Currently to my knowledge we are using (regular) random matrices in our unit tests, where the condition numbers are not controlled, so convergence may fail even if there is no bug in the implementation of the methods. Therefore I would like to propose to use random matrices with controlled condition numbers. An example for this technique can be found in the FORTRAN 77 implementation of the LSQR algorithm by Paige an Saunders, see <http://www.stanford.edu/group/SOL/software/lsqr/> and the comments in the code. Random matrices with controlled condition numbers can be constructed by setting up their singular value decomposition (SVD), namely we first construct a diagonal matrix with the singular values on the diagonal, so its condition number (and e.g. positive definiteness in the symmetric case) can easily be prescribed. Then the random component is added by multiplying the diagonal matrix by appropriate orthogonal matrices from the left and with the transpose from the right. Note that in the symmetric case (e.g. for CG) the orthogonal matrices used from the left and the right should be the same to preserve symmetry. Random orthogonal matrices can easily be constructed by Householder's method using random Householder vectors.

Just a quick note: the BiCGSTAB and CG solvers work remarkably well on these random matrices, so I think that a GMRES should also be able to handle them well. However, to achieve this stability in BiCGSTAB I had to implement a couple of tricks to deal with numerical degeneracies. Such tricks are never discussed in books. So if the GMRES implementation is a somewhat direct implementation of a textbook algorithm, it is not surprising that it lacks numerical stability.

(In reply to comment #12) > Just a quick note: the BiCGSTAB and CG solvers work remarkably well on these > random matrices, so I think that a GMRES should also be able to handle them > well. However, to achieve this stability in BiCGSTAB I had to implement a > couple of tricks to deal with numerical degeneracies. Such tricks are never > discussed in books. So if the GMRES implementation is a somewhat direct > implementation of a textbook algorithm, it is not surprising that it lacks > numerical stability. Hi Gael, could you please specify what "tricks" you are referring to? Maybe these can be used in GMRES and MINRES too.

Here is how I proceeded: 1 - pick a matrix that fails (and check that there is indeed nothing really wrong about it). 2 - then track the algorithm to spot the line where bad things is happening (a quantity vanish while it should not, NaN, inf, etc.) 3 - fix it! The last step is of course the most tricky part. IIRC, for BiCGSTAB I had to deal with some divisions by zero, and one had to be fixed by returning 0, and the other by implementing a restart mechanism. Anyway, understanding clearly what's going wrong is the first thing to do.

The problem with gmres_2 seems to arise only when GMRES is applied to a complex case but works well in the case of real data. Probably there is a complex conjugation missing somewhere. Also gmres_1 previously failed quite often, see bug 841 for a fix. Indeed, Gael, you are right that GMRES now usually handles the random matrices well, just as BICGSTAB does. So I do not see that the GMRES implementation is numerically unstable, the reason was a different one, see bug 841. Let me know if you find numerical instabilities.

Created attachment 481 [details] second patch for bug 715: unit test fixed

Here is another fix for the GMRES implementation: * Fix error in calculation of residual at restart. * Use relative residual as stopping criterion. * Improve documentation. Now, as far as I can see, the unit tests no longer fail. Moreover, Gael, BICGSTAB and GMRES work well with random matrices as soon as a preconditioner is used, see also test/bicgstab.cpp and unsupported/test/gmres.cpp.

(In reply to Kolja Brix from comment #16) > Created attachment 481 [details] > second patch for bug 715: unit test fixed Pushed, thanks for the patch! https://bitbucket.org/eigen/eigen/commits/47eb59a2ea6 levenberg_marquardt still fails on my machines, in comment 6 I already bisected the commit which made the tests fail (though I think, it's more a problem how the tests are implemented). Any volunteer to look into that?

hm, the problem in LM is that if the vector contains NaN, then the optimizations in stableNorm skip it, and it will return 0 instead of NaN. I'll fix that, though it feels weird that LM call stableNorm on vectors containing NaNs...

For the NaN in stableNorm, see: https://bitbucket.org/eigen/eigen/commits/c3b3502d75b0

-- GitLab Migration Automatic Message -- This bug has been migrated to gitlab.com's GitLab instance and has been closed from further activity. You can subscribe and participate further through the new bug through this link to our GitLab instance: https://gitlab.com/libeigen/eigen/issues/715.