bugGNU Octave - Bugs: bug #60357, Failing test for `qz` with LAPACK...

 
 

bug #60357: Failing test for `qz` with LAPACK 3.9.1

Submitter:  Markus Mützel <mmuetzel>
Submitted:  Wed 07 Apr 2021 07:41:38 PM UTC
   
 
Category:  Libraries Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  None
Originator Name:  Open/Closed:  * Closed
Release:  * 6.2.0 Operating System:  * Any
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Sun 25 Apr 2021 04:48:08 PM UTC, comment #25: 

Ok. I grafted those changes to stable here:
https://hg.savannah.gnu.org/hgweb/octave/rev/d75aa2bf4915
https://hg.savannah.gnu.org/hgweb/octave/rev/50aa945836c3

Marking as ready for test again.

Markus Mützel <mmuetzel>
Group administrator
Sun 25 Apr 2021 04:19:02 PM UTC, comment #24: 

I would be okay with just backporting the change.  It's not just a matter of Matlab compatibility, it's actually incorrect that Octave returns only 2 eigenvalues when there are 3 eigenvalues.

Rik <rik5>
Group administrator
Sun 25 Apr 2021 03:00:16 PM UTC, comment #23: 

The test still fails on the stable branch with LAPACK 3.9.1.
I'm not sure if we should backport that change to the stable branch.
Maybe we should better re-work the BIST on the stable branch. Or just mark it with this bug number so that it is marked as a known error.

What do you think?

Re-opening for the stable branch.

Markus Mützel <mmuetzel>
Group administrator
Sun 25 Apr 2021 11:48:54 AM UTC, comment #22: 

Seems to be working on Windows, too.

Closing as fixed.

Markus Mützel <mmuetzel>
Group administrator
Fri 23 Apr 2021 05:53:22 PM UTC, comment #21: 

@Markus: Thanks for making that change.  It's what I expected would have to be done.

Rik <rik5>
Group administrator
Fri 23 Apr 2021 05:33:52 PM UTC, comment #20: 

Thanks, the BIST didn't pass for me with LAPACK 3.9.1 or LAPACK 3.9.0.
I re-phrased it slightly to have it pass for me with both versions on Ubuntu:
https://hg.savannah.gnu.org/hgweb/octave/rev/33556123b892

I haven't checked on Windows yet.

Markus Mützel <mmuetzel>
Group administrator
Fri 23 Apr 2021 03:11:35 PM UTC, comment #19: 

I checked in a change (http://hg.savannah.gnu.org/hgweb/octave/rev/06f4bb661d96) so that qz() always returns the correct number of eigenvalues.  When LAPACK doesn't calculate an eigenvalue, because the rank of the matrix is lower than the size of the matrix, then I changed qz() to return Inf.  This matches what the eig() function does for the same situation which at least is now consistenct.  That is, the eigenvalues calculated by eig() or by qz() are the same.

I changed the BIST tests to reflect the new reality.  However, my guess is that the test still fails with LAPACK 3.9.1 because the eigenvalue returned by LAPACK is large (3x10^15), but not actually Inf.  In any case, I'm marking as "Ready for Test" so Markus can check how the new code behaves with his system.

Rik <rik5>
Group administrator
Tue 20 Apr 2021 03:06:41 PM UTC, comment #18: 

Incidentally, the same issue exists for gsvd().  When calling gsvd (a,b) using the matrices from this bug report, only two generalized eigenvalues are returned.  This is also incorrect and should be fixed.  I checked with Matlab and gsvd should also return a 3x1 vector.

Rik <rik5>
Group administrator
Tue 20 Apr 2021 03:04:58 PM UTC, comment #17: 

Attached is a one-line patch that at least means qz returns 3 eigenvalues.  The patch is


          //tmp.resize (cnt);  // Trim vector to number of return values


Octave was trimming down the eigenvalue matrix which was originally created to size 3x1 to the number of eigenvalues returned from LAPACK.  Since the third eigenvalue can take any value, why not make it 0?  That's what this patch effectively does since the matrix is effectively created with all 0 values to begin with.

(file #51298)

Rik <rik5>
Group administrator
Sun 18 Apr 2021 11:35:34 PM UTC, comment #16: 

Some interesting data points on this bug report.

First, Matlab no longer documents that there is a final output from qz() which is the eigenvalues for Octave.  See https://www.mathworks.com/help/matlab/ref/qz.html.  Maybe that means something?  Maybe the eigenvalue output is unstable?  Should Octave drop support for the eigenvalue output?

Second, for the test case in this bug report, all outputs (aa, bb, q, z, v, w) agree with Matlab.  This part of the calculation seems correct.

Third, although undocumented, you can ask for a seventh output from Matlab's qz() function and it appears to be the eigenvalues.  The output is always 3x1, as expected, because there should be 3 eigenvalues for 3x3 input matrices A and B.  The fact that Octave returns only 2 eigenvalues is incorrect.  It should return 3 eigenvalues, as it does with the LAPACK 3.9.1.

Fourth, there seems to be some convergence/roundoff issues in Matlab.  The first 6 outputs are always the same, but the eigenvalues are sometimes all 0 ([0;0;0]) and sometimes something very small (10^-306 x [1; 1; 1]) which is effectively all zeros.  The two Octave eigenvalues are of normal size ([-4.2, 0.2]) which seems wrong.

Fifth, obviously this is a very hard test case for the QZ algorithm.  The rank of the matrix B is 2 which means the last eigenvalue is not really determined and could be anything.



Rik <rik5>
Group administrator
Fri 09 Apr 2021 11:46:04 AM UTC, comment #15: 

Looks like I misspelled some tags. Next attempt:

Re comment #12:
With LAPACK 3.9.1, I see the following:

> a = [ 10 1 2;
1 2 -1;
1 1 2];
b = reshape (1:9,3,3);
[aa, bb, q, z, v, w, lambda] = qz (a, b);
>> a*v

ans =

5.000000 7.659527 -1.273807
-2.000000 3.806974 -1.526360
0.500000 -0.045579 -1.778912

>> b*v*diag(lambda)

ans =

-3.101857 7.659527 -1.273807
-4.652786 3.806974 -1.526360
-6.203714 -0.045579 -1.778912


If I read the documentation of `qz` correctly, the last two matrices should be the same. But the first column differs by a large amount.
IIUC, that column corresponds to the eigenvalue that differs by several orders of magnitude from the other two values.

`b*v` is close to double `eps` for that eigenvalue:

>> b*v
ans =

-8.8818e-16 -1.7998e+00 -6.8665e+00
-1.3323e-15 -8.9456e-01 -8.2279e+00
-1.7764e-15 1.0710e-02 -9.5893e+00


It might be good if someone with a background in linear algebra could have a look.
Should we postprocess the results from LAPACK and discard the "outlier" eigenvalue?

Markus Mützel <mmuetzel>
Group administrator
Fri 09 Apr 2021 11:43:39 AM UTC, comment #14: 

Re comment #12:
With LAPACK 3.9.1, I see the following:

> a = [ 10  1  2;
       1  2 -1;
       1  1  2];
b = reshape (1:9,3,3);
[aa, bb, q, z, v, w, lambda] = qz (a, b);
>> a*v

ans =

   5.000000   7.659527  -1.273807
  -2.000000   3.806974  -1.526360
   0.500000  -0.045579  -1.778912

>> b*v*diag(lambda)

ans =

  -3.101857   7.659527  -1.273807
  -4.652786   3.806974  -1.526360
  -6.203714  -0.045579  -1.778912
+verbatim+

If I read the documentation of `qz` correctly, the last two matrices should be the same. But the first column differs by a large amount.
IIUC, that column corresponds to the eigenvalue that differs by several orders of magnitude from the other two values.

`b*v` is close to double `eps` for that eigenvalue:
+verbatim+
>> b*v
ans =

  -8.8818e-16  -1.7998e+00  -6.8665e+00
  -1.3323e-15  -8.9456e-01  -8.2279e+00
  -1.7764e-15   1.0710e-02  -9.5893e+00


It might be good if someone with a background in linear algebra could have a look.
Should we postprocess the results from LAPACK and discard the "outlier" eigenvalue?

Markus Mützel <mmuetzel>
Group administrator
Thu 08 Apr 2021 07:13:34 PM UTC, comment #13: 

OK, if I understand correctly, we need two test blocks because we expect a failure on some systems but not others.  If we have just one test that is marked as fixed, we will get a regression on the system where the test fails.  But we've also closed the report that is listed in the expected failure test.  So the update-bug-status script will tag that as fixed and we'll still see a regression.

If we accept that the test fails in a specific way on some systems, do we want to see a failure or should we instead fix the test so it "succeeds" when it fails in the expected way?  If yes, then we can have a single test that is tagged as fixed.  Otherwise, I'm OK with separate tests with one expected to fail, but then the associated report should not be closed with a "Fixed" status so that the update-bug-status script doesn't mark it as fixed in the test suite.

In the particular case of tests marked with bug 52627, that bug is about a bunch of different problems that were showing up as failures on Mac systems.  We could leave that closed as Fixed and open another just for this specific issue with mapper functions.  It could either remain open or be closed with a "Won't Fix" status.  Either way will prevent the update-bug-status script from marking these tests as fixed in the test suite.

John W. Eaton <jwe>
Group administrator
Thu 08 Apr 2021 06:42:34 PM UTC, comment #12: 

Interesting, so we may have been covering for a problem in qz all along.  For a 3x3 matrix there should be exactly 3 eigenvalues, always.  As you say, it seems in 3.9.0 it was returning 2.

Rik <rik5>
Group administrator
Thu 08 Apr 2021 06:32:59 PM UTC, comment #11: 

@Rik: The issue is already visible is comment #1. `lambda` has 3 elements with LAPACK 3.9.1 but only 2 elements with LAPACK 3.9.0.
Looking at the documentation of `qz`, I would have assumed that 3 is the number of elements that we should be expecting. So in a sense, the results with 3.9.0 are unexpected.
But the BIST adds an additional `0` to the vector.

Markus Mützel <mmuetzel>
Group administrator
Thu 08 Apr 2021 06:28:13 PM UTC, comment #10: 

@Markus: Can you run the parts of the BIST test up through the call to qz(), but not the asserts, and post the results as you did in comment #1?

The errors from Octave seem to indicate that even the return size of one of the matrices is incorrect.

Rik <rik5>
Group administrator
Thu 08 Apr 2021 06:22:22 PM UTC, comment #9: 

I introduced these constructs a long time back before we had amended test.m to handle bug numbers and a regression indicator.

The point was to have accurate results from the regression test suite.  There are tests in Octave which pass just fine on Linux, but get different answers or require overly large tolerances on Mac systems.  This is a problem with Mac systems and their libraries, not something that needs to be inflicted on other users.

The solution at the time was to have two BIST test blocks.  One that would run only on non-mac platforms and would pass and contribute to the number of successful tests, and a second test for mac platforms which would fail and contribute to the number of unsuccessful tests.

Rik <rik5>
Group administrator
Thu 08 Apr 2021 06:03:17 PM UTC, comment #8: 

Re comment #6, jwe: IIRC, we split those tests to have them marked as a regression on platforms where we know they once passed. But don't have regressions on platforms where they never passed.
The body of the BISTs should be identical.

I'm not sure though why we skipped the original test on platforms where they never passed, but have the other construct for having the duplicated test pass on platforms that are marked as fixed.

I don't know why the tests are marked with different numbers. I assumed they were discussed in different reports.

Markus Mützel <mmuetzel>
Group administrator
Thu 08 Apr 2021 05:57:18 PM UTC, comment #7: 

We are using `dhgeqz` in "qz.cc".

IIUC, that function was changed for LAPACK 3.9.1:
https://github.com/Reference-LAPACK/lapack/pull/505

I don't know if that is related though.

Markus Mützel <mmuetzel>
Group administrator
Thu 08 Apr 2021 05:50:57 PM UTC, comment #6: 

Markus: I started making the same change as you did:


diff --git a/libinterp/corefcn/mappers.cc b/libinterp/corefcn/mappers.cc
--- a/libinterp/corefcn/mappers.cc
+++ b/libinterp/corefcn/mappers.cc
@@ -133,10 +133,9 @@ Compute the inverse cosine in radians fo
 %! v = [0, pi, pi/2, pi/2];
 %! assert (real (acos (x)), v);

-%!test <52627>
-%! ## Same test code as above, but intended only for test statistics on Mac and
-%! ## Windows.  Their trig/hyperbolic functions have huge tolerances.
-%! if (! ismac () && ! ispc ()), return; endif
+## Same test code as above, but intended only for test statistics on Mac and
+## Windows.  Their trig/hyperbolic functions have huge tolerances.
+%!testif ; ismac () || ispc () <*52627>
 %! x = [1, -1, i, -i] .* 1e150;
 %! v = [0, pi, pi/2, pi/2];
 %! assert (real (acos (x)), v);


and then noticed that for each of these cases, the code in the two sets of tests are now identical.  That may not have always been true?  But unless I'm missing something, it suggests we should merge them.

I also noticed that the two sets of tests use different bug numbers.  We don't currently have a way to specify more than one bug number for a test (build-aux/update-bug-script.sh and test.m would need to be modified to handle that).  I suppose we could just add additional bug numbers in a comment.  For bugs that have been fixed, the real point is to be able to say that any failure is a regression, so the particular bug number is less important than the "*".

John W. Eaton <jwe>
Group administrator
Thu 08 Apr 2021 05:35:34 PM UTC, comment #5: 

Re comment #3, @jwe: I'm attaching a patch that does the change you propose. I can't remember correctly. But I have the vague feeling that we used that construct for a specific reason.

(file #51217)

Markus Mützel <mmuetzel>
Group administrator
Thu 08 Apr 2021 05:15:24 PM UTC, comment #4: 

I get the same (wrong) result with LAPACK 3.9.1 on Linux.

Markus Mützel <mmuetzel>
Group administrator
Thu 08 Apr 2021 05:06:20 PM UTC, comment #3: 

These bugs will be marked with "*" again the next time someone runs "make update-bug-status" and the problem will resurface. Instead of embedding the OS test in the test itself, we should probably be using things like


%!testif ; ismac () || ispc () <*52627>
...


John W. Eaton <jwe>
Group administrator
Thu 08 Apr 2021 04:55:51 PM UTC, comment #2: 

This might be LAPACK problem.  Can you build LAPACK 3.9.1 on Linux and check?

Incidentally, LAPACK 3.9.1 is very recent (7 days old).  See this project page http://performance.netlib.org/lapack/lapack-3.9.1.html.

The bugs that were fixed are listed here in case you see something relevant: https://github.com/Reference-LAPACK/lapack/issues?q=is%3Aclosed+milestone%3A%22LAPACK+3.9.1%22


Rik <rik5>
Group administrator
Thu 08 Apr 2021 04:34:43 PM UTC, comment #1: 

With a clean build from hg id d57c1d781093 and MXE Octave hg id 60a500c74ac7, the error with "zip" disappeared.
A bunch of BISTs erroneously were marked as regressions. Should be fixed with this change:
https://hg.savannah.gnu.org/hgweb/octave/rev/526a5ba817ca

The BIST in qz.cc-tst still fails.
Running the beginning of the BIST:

>> a = [ 10  1  2;
        1  2 -1;
        1  1  2];
>> b = reshape (1:9,3,3);
>> [aa, bb, q, z, v, w, lambda] = qz (a, b)
aa =

   4.4159   6.8422  -4.4675
        0  -4.6396   2.0919
        0        0   2.1964

bb =

    0.0000   -0.9650   -4.8929
         0    1.0902   10.8977
         0         0   11.8398

q =

   0.924500  -0.369800   0.092450
  -0.377574  -0.921698   0.088946
   0.052319  -0.117137  -0.991737

z =

   0.4082   0.8731  -0.2667
  -0.8165   0.2186  -0.5344
   0.4082  -0.4359  -0.8021

v =

   0.500000   0.856123   0.056123
  -1.000000   1.000000  -1.000000
   0.500000  -0.950851  -0.417518

w =

   0.500000  -0.424485   0.052755
  -1.000000  -0.811212  -0.118113
   0.500000   1.000000  -1.000000

lambda =

   3.4924e+15
  -4.2557e+00
   1.8551e-01


The same on Ubuntu 20.10 with the same hg id:

>> a = [ 10  1  2;
        1  2 -1;
        1  1  2];
>> b = reshape (1:9,3,3);
>> [aa, bb, q, z, v, w, lambda] = qz (a, b)
aa =

  -5.2673   0.0102  -8.0910
        0   1.9024   0.0690
        0        0   4.4907

bb =

    1.2377    9.8871   -5.7329
         0   10.2551   -6.9050
         0         0         0

q =

  -8.9548e-01  -4.4507e-01   5.3287e-03
   1.7735e-01  -3.6775e-01  -9.1286e-01
   4.0825e-01  -8.1650e-01   4.0825e-01

z =

   0.5272   0.2294   0.8182
   0.6158  -0.7666  -0.1818
  -0.5855  -0.5997   0.5455

v =

   0.856123   0.056123   0.500000
   1.000000  -1.000000  -1.000000
  -0.950851  -0.417518   0.500000

w =

  -0.424485   0.052755   0.500000
  -0.811212  -0.118113  -1.000000
   1.000000  -1.000000   0.500000

lambda =

  -4.2557
   0.1855

>> version -lapack
ans = Linear Algebra PACKage Version 3.9.0



That looks quite different to me.

If I drop in the "liblapack.dll" from the build from 5th of April (LAPACK version 3.9.0), the results on Windows and Linux match.

Is this a bug in LAPACK 3.9.1? Or an error in our build rule in MXE Octave? Or do we use those functions incorrectly?

CC'ing John D.

Markus Mützel <mmuetzel>
Group administrator
Wed 07 Apr 2021 07:41:38 PM UTC, original submission:  

Since roughly one day ago some BISTs are failing on Windows for me:

>>>>> processing D:\SVN\Octave\test\octave-2021-04-07-19-01-w64 non-ASCII\octave-2021-04-07-19-01-w64\mingw64\share\octave\7.0.0\etc\tests\libinterp\corefcn\qz.cc-tst
***** test
 a = [ 10  1  2;
        1  2 -1;
        1  1  2];
 b = reshape (1:9,3,3);
 [aa, bb, q, z, v, w, lambda] = qz (a, b);
 sz = length (lambda);
 observed = (b * v * diag ([lambda;0])) (:, 1:sz);
 assert ((a*v)(:, 1:sz), observed, norm (observed) * 1e-14);
 observed = (diag ([lambda;0]) * w' * b) (1:sz, :);
 assert ((w'*a)(1:sz, :) , observed, norm (observed) * 1e-13);
 assert (q * a * z, aa, norm (aa) * 1e-14);
 assert (q * b * z, bb, norm (bb) * 1e-14);
!!!!! test failed
operator *: nonconformant arguments (op1 is 3x3, op2 is 4x4)
shared variables   scalar structure containing the fields:

    a =

       1   2
       0   3

    b =

       1   0
       0   0

    c =

       0   1
       0   0


>>>>> processing D:\SVN\Octave\test\octave-2021-04-07-19-01-w64 non-ASCII\octave-2021-04-07-19-01-w64\mingw64\share\octave\7.0.0\m\miscellaneous\zip.m
***** test
 ## test zip together with unzip
 unwind_protect
   filename = tempname ();
   tmp_var  = pi;
   save (filename, "tmp_var");
   dirname = tempname ();
   mkdir (dirname);
   zipfile = tempname ();
   [~, basename, ext] = fileparts (filename);
   filelist = zip (zipfile, [basename ext], tempdir);
   filelist = filelist{1};
   if (! strcmp (filelist, [basename ext]))
     error ("zip archive does not contain expected name!");
   endif
   if (! exist ([zipfile ".zip"], "file"))
     error ("zip file cannot be found!");
   endif
   unzip ([zipfile ".zip"], dirname);
   fid = fopen (filename, "rb");
   assert (fid >= 0);
   orig_data = fread (fid);
   fclose (fid);
   fid = fopen ([dirname filesep basename ext], "rb");
   assert (fid >= 0);
   new_data = fread (fid);
   fclose (fid);
   if (orig_data != new_data)
     error ("unzipped file not equal to original file!");
   endif
 unwind_protect_cleanup
   sts = unlink (filename);
   sts = unlink ([dirname, filesep, basename, ext]);
   sts = unlink (zipfile);
   sts = unlink ([zipfile ".zip"]);
   sts = rmdir (dirname);
 end_unwind_protect
!!!!! test failed
zip: zip failed with exit status = 15


The two tests don't look related. But they started failing with the same build. (That might be incidental.)

I was building with local changes that are stripped now. So, I can only give a rough range.
Last good: build from 5th of April (MXE hg id 05703685a4a8, Octave hg id unknown)
First bad: build from 6th of April

I built from the default MXE branch. IIRC, the update of lapack, lua and curl falls into that range.
I don't remember which revisions my local changes in Octave were based on.

I did incremental builds. I'll check if the errors disappear with a clean build.

Markus Mützel <mmuetzel>
Group administrator

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #51298:  60357.patch added by rik5 (505B - text/x-patch)
file #51217:  bug60357_platform-specific-BISTs.patch added by mmuetzel (4KiB - application/octet-stream)

 

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

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 25 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2021-04-25 rik5 StatusReady For Test Fixed
        Open/ClosedOpen Closed
    2021-04-25 mmuetzel StatusIn Progress Ready For Test
    2021-04-25 mmuetzel StatusFixed In Progress
        Open/ClosedClosed Open
        Releasedev 6.2.0
    2021-04-25 mmuetzel StatusReady For Test Fixed
        Open/ClosedOpen Closed
    2021-04-23 rik5 StatusConfirmed Ready For Test
    2021-04-20 rik5 Attached File- Added 60357.patch, #51298
    2021-04-18 rik5 Item GroupNone Incorrect Result
        StatusNone Confirmed
    2021-04-08 jwe Item GroupIncorrect Result None
        StatusConfirmed None
    2021-04-08 rik5 StatusNone Confirmed
    2021-04-08 rik5 Item GroupNone Incorrect Result
        Operating SystemMicrosoft Windows Any
    2021-04-08 rik5 Operating SystemAny Microsoft Windows
    2021-04-08 mmuetzel Operating SystemMicrosoft Windows Any
    2021-04-08 jwe Operating SystemAny Microsoft Windows
    2021-04-08 mmuetzel Attached File- Added bug60357_platform-specific-BISTs.patch, #51217
    2021-04-08 mmuetzel Operating SystemMicrosoft Windows Any
        SummaryFailing test for `qz` on Windows with LAPACK 3.9.1 Failing test for `qz` with LAPACK 3.9.1
    2021-04-08 mmuetzel CategoryInterpreter Libraries
        SummaryFailing tests on Windows Failing test for `qz` on Windows with LAPACK 3.9.1

    Back to the top

    Powered by Savane 3.13-f8d8.
    Corresponding source code