Opened 12 years ago

Closed 11 years ago

#10795 closed defect (fixed)

Fix and upgrade double dense matrix QR decomposition

Reported by: Rob Beezer Owned by: jason, was
Priority: minor Milestone: sage-5.2
Component: linear algebra Keywords: sd40.5
Cc: Jason Grout, Dan Drake Merged in: sage-5.2.beta1
Authors: Rob Beezer Reviewers: Martin Raum, Dan Drake
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Status badges

Description (last modified by Rob Beezer)

The Q matrix of a QR decomposition should be unitary, hence invertible. For zero-column trivial cases, this is broken.

sage: A = zero_matrix(CDF, 5, 0)
sage: Q, R = A.QR()
sage: Q
[0 0 0 0 0]
[0 0 0 0 0]
[0 0 0 0 0]
[0 0 0 0 0]
[0 0 0 0 0]

Besides a bugfix this patch will upgrade the documentation to make it clear how this routine works over the complex numbers. In particular, SciPy routines are using a Hermitian inner product - documentation upgrade will reflect that.

Apply:

  1. trac_10795-QR-decomposition-double-dense-v2.patch
  2. trac_10795-QR-decomposition-formatting.patch

Attachments (2)

trac_10795-QR-decomposition-double-dense-v2.patch (8.9 KB) - added by Rob Beezer 11 years ago.
trac_10795-QR-decomposition-formatting.patch (3.9 KB) - added by Rob Beezer 11 years ago.

Download all attachments as: .zip

Change History (32)

comment:1 Changed 12 years ago by Rob Beezer

Authors: Rob Beezer
Status: newneeds_review

Just one change to the code, to return an identity matrix in a trivial case. Everything else is documentation. Unfortunately, this file is not included in the documentation (yet). You can check that it builds without warnings, then look at it in the notebook to verify the contents and appearance.

comment:2 Changed 11 years ago by Martin Raum

Reviewers: Martin Raum
Status: needs_reviewpositive_review

Everything is OK with this patch.

comment:3 Changed 11 years ago by Jeroen Demeyer

Milestone: sage-4.7.1sage-4.7.2

comment:4 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-4.7.2.alpha2
Resolution: fixed
Status: positive_reviewclosed

comment:5 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-4.7.2.alpha2
Resolution: fixed
Status: closednew

On various systems, there are doctest failures:

hawk (OpenSolaris? 06/2009 i86pc Xeon W3580):

**********************************************************************
File "/export/home/buildbot/build/sage/hawk-1/hawk_full/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 1766:
    sage: Q
Expected:
    [ -0.359210604054   0.569326179705   0.368048420509   0.641385845805]
    [  0.179605302027  -0.144590775798   0.925041158846  -0.301884576418]
    [  0.179605302027  -0.704880032016  0.0774617736597   0.681825307224]
    [  0.898026510134   0.397624633445 -0.0532812182975   0.180566192161]
Got:
    [-0.359210604054  0.569326179705  0.409076682956  0.616028985113]
    [ 0.179605302027 -0.144590775798 -0.683704756325   0.69237507842]
    [ 0.179605302027 -0.704880032016  0.575201135981  0.374205463771]
    [ 0.898026510134  0.397624633445  0.185331397251 0.0330954856071]
**********************************************************************
File "/export/home/buildbot/build/sage/hawk-1/hawk_full/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 1771:
    sage: R
Expected:
    [    -5.56776436283       2.6940795304      -2.6940795304]
    [                 0     -3.56958477752      3.56958477752]
    [                 0                  0 -9.93013661299e-16]
    [                 0                  0                  0]
Got:
    [   -5.56776436283      2.6940795304     -2.6940795304]
    [                0    -3.56958477752     3.56958477752]
    [                0                 0 -4.4408920985e-16]
    [                0                 0                 0]
**********************************************************************

bsd (OS X 10.8.0 x86_64) 64-bit:

**********************************************************************
File "/Users/buildbot/build/sage/bsd-1/bsd_64_full/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 1766:
    sage: Q
Expected:
    [ -0.359210604054   0.569326179705   0.368048420509   0.641385845805]
    [  0.179605302027  -0.144590775798   0.925041158846  -0.301884576418]
    [  0.179605302027  -0.704880032016  0.0774617736597   0.681825307224]
    [  0.898026510134   0.397624633445 -0.0532812182975   0.180566192161]
Got:
    [ -0.359210604054   0.569326179705   0.146337376237   0.724859169325]
    [  0.179605302027  -0.144590775798   0.973035382602 0.00613084354877]
    [  0.179605302027  -0.704880032016  -0.142125402809   0.671331844787]
    [  0.898026510134   0.397624633445  -0.107647045464   0.154451130063]
**********************************************************************
File "/Users/buildbot/build/sage/bsd-1/bsd_64_full/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 1771:
    sage: R
Expected:
    [    -5.56776436283       2.6940795304      -2.6940795304]
    [                 0     -3.56958477752      3.56958477752]
    [                 0                  0 -9.93013661299e-16]
    [                 0                  0                  0]
Got:
    [   -5.56776436283      2.6940795304     -2.6940795304]
    [                0    -3.56958477752     3.56958477752]
    [                0                 0 6.28036983474e-16]
    [                0                 0                 0]
**********************************************************************

cleo (RHEL 5.3 ia64 Itanium 2):

**********************************************************************
File "/home/buildbot/build/sage/cleo-1/cleo_full/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 1766:
    sage: Q
Expected:
    [ -0.359210604054   0.569326179705   0.368048420509   0.641385845805]
    [  0.179605302027  -0.144590775798   0.925041158846  -0.301884576418]
    [  0.179605302027  -0.704880032016  0.0774617736597   0.681825307224]
    [  0.898026510134   0.397624633445 -0.0532812182975   0.180566192161]
Got:
    [ -0.359210604054   0.569326179705   0.146859067236   0.724753652911]
    [  0.179605302027  -0.144590775798   0.973039543327 0.00543048428303]
    [  0.179605302027  -0.704880032016  -0.141642164231   0.671433967908]
    [  0.898026510134   0.397624633445  -0.107535848925   0.154528570726]
**********************************************************************
File "/home/buildbot/build/sage/cleo-1/cleo_full/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 1771:
    sage: R
Expected:
    [    -5.56776436283       2.6940795304      -2.6940795304]
    [                 0     -3.56958477752      3.56958477752]
    [                 0                  0 -9.93013661299e-16]
    [                 0                  0                  0]
Got:
    [   -5.56776436283      2.6940795304     -2.6940795304]
    [                0    -3.56958477752     3.56958477752]
    [                0                 0 9.81879219451e-16]
    [                0                 0                 0]
**********************************************************************

comment:6 Changed 11 years ago by Jeroen Demeyer

Status: newneeds_review

It seems that this patch causes failures in a non-reproducible way. Is there some randomness in the algorithm? For example, I just got this failure on sage.math.washington.edu (a machine on which I remember the test succeeding):

sage -t -long  -force_lib devel/sage/sage/matrix/matrix_double_dense.pyx
**********************************************************************
File
"/mnt/usb1/scratch/buildbot/sage/sage-1/sage_binary/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx",
line 1766:
    sage: Q
Expected:
    [ -0.359210604054   0.569326179705   0.368048420509   0.641385845805]
    [  0.179605302027  -0.144590775798   0.925041158846  -0.301884576418]
    [  0.179605302027  -0.704880032016  0.0774617736597   0.681825307224]
    [  0.898026510134   0.397624633445 -0.0532812182975   0.180566192161]
Got:
    [-0.359210604054  0.569326179705  0.616028985113  0.409076682956]
    [ 0.179605302027 -0.144590775798   0.69237507842 -0.683704756325]
    [ 0.179605302027 -0.704880032016  0.374205463771  0.575201135981]
    [ 0.898026510134  0.397624633445 0.0330954856071  0.185331397251]
**********************************************************************
File
"/mnt/usb1/scratch/buildbot/sage/sage-1/sage_binary/build/sage-4.7.2.alpha2/devel/sage-main/sage/matrix/matrix_double_dense.pyx",
line 1771:
    sage: R
Expected:
    [    -5.56776436283       2.6940795304      -2.6940795304]
    [                 0     -3.56958477752      3.56958477752]
    [                 0                  0 -9.93013661299e-16]
    [                 0                  0                  0]
Got:
    [   -5.56776436283      2.6940795304     -2.6940795304]
    [                0    -3.56958477752     3.56958477752]
    [                0                 0 -4.4408920985e-16]
    [                0                 0                 0]
**********************************************************************

comment:7 Changed 11 years ago by Jeroen Demeyer

Status: needs_reviewneeds_work

comment:8 Changed 11 years ago by Rob Beezer

Description: modified (diff)
Keywords: beginner removed
Status: needs_workneeds_review

"numerical" patch applies accumulated techniques for these numerical computations and should address the doctest failures.

comment:9 Changed 11 years ago by Jason Grout

Cc: Jason Grout added

Same question on here as on #10791 -- is it appropriate to set this back to positive review (assuming the tests still pass on Rob's machine)?

comment:10 Changed 11 years ago by Rob Beezer

Applies and passes all tests on 4.8.alpha3.

comment:11 Changed 11 years ago by Jason Grout

Status: needs_reviewpositive_review

Based on jdmeyer's response to my question on #10791: the fixing patch looks reasonable and Rob tested it. I'd test it as well, except I'm currently building python 2.7 and testing it.

comment:12 Changed 11 years ago by Jason Grout

Sorry, I meant jdemeyer's response...

comment:13 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-4.8.alpha5
Resolution: fixed
Status: positive_reviewclosed

comment:14 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-4.8.alpha5
Resolution: fixed
Status: closednew

This sometimes gives a doctest failure:

sage -t -long  -force_lib devel/sage/sage/matrix/matrix_double_dense.pyx
**********************************************************************
File "/mnt/usb1/scratch/buildbot/sage/sage-1/sage_binary/build/sage-4.8.alpha5/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 2243:
    sage: Q.round(6).zero_at(10^-6)
Expected:
    [-0.359211  0.569326  0.368048  0.641386]
    [ 0.179605 -0.144591  0.925041 -0.301885]
    [ 0.179605  -0.70488  0.077462  0.681825]
    [ 0.898027  0.397625 -0.053281  0.180566]
Got:
    [-0.359211  0.569326 -0.631992 -0.383955]
    [ 0.179605 -0.144591   -0.6643  0.711014]
    [ 0.179605  -0.70488  -0.39705 -0.559676]
    [ 0.898027  0.397625 -0.040527 -0.183849]
**********************************************************************

It's not clear why it happens. Because on the same machine I have one Sage build where it works, one where it doesn't with the exact same source, but different gcc version.

Does this code use any external libraries (ATLAS comes to mind)?

comment:15 Changed 11 years ago by Jeroen Demeyer

Status: newneeds_review

On Fedora 15-64 (flavius):

sage -t -long  -force_lib devel/sage/sage/matrix/matrix_double_dense.pyx
**********************************************************************
File "/home/buildbot/build/sage/flavius-1/flavius_full/build/sage-4.8.alpha5/devel/sage-main/sage/matrix/matrix_double_dense.pyx", line 2243:
    sage: Q.round(6).zero_at(10^-6)
Expected:
    [-0.359211  0.569326  0.368048  0.641386]
    [ 0.179605 -0.144591  0.925041 -0.301885]
    [ 0.179605  -0.70488  0.077462  0.681825]
    [ 0.898027  0.397625 -0.053281  0.180566]
Got:
    [-0.359211  0.569326  0.616029  0.409077]
    [ 0.179605 -0.144591  0.692375 -0.683705]
    [ 0.179605  -0.70488  0.374205  0.575201]
    [ 0.898027  0.397625  0.033095  0.185331]
**********************************************************************

comment:16 Changed 11 years ago by Jeroen Demeyer

Status: needs_reviewneeds_work

comment:17 Changed 11 years ago by Rob Beezer

Thanks, Jeroen. Sorry for the delay in getting back to this one.

Yes, this uses SciPy, which in turn will use ATLAS. It seems a 2-dimensional eigenspace is being returned with different basis vectors.

This is mostly a documentation patch, so I think it would be best to make the tests more functional - checking that returned matrices have the desired properties - and less about the actual entries of the matrices. I'll put that on my list.

Rob

comment:18 Changed 11 years ago by Jeroen Demeyer

* bump *

Changed 11 years ago by Rob Beezer

comment:19 Changed 11 years ago by Rob Beezer

Cc: Dan Drake added
Keywords: sd40.5 added
Status: needs_workneeds_review
  • Jeroen - thanks for the bump.
  • Needed a rebase, so did that and also tidied up a few things. Almost entirely documentation, with one code change. Will need a new review since the documentation has changed.
  • Doctest failures. Seems different platforms give different bases for the dimension 2 left kernel of the matrix. No harm in just removing it and explaining the situation, along with definitive evidence that the result behaves as expected.

comment:20 Changed 11 years ago by Rob Beezer

Description: modified (diff)

comment:21 Changed 11 years ago by Martin Raum

This looks good, and I think it is reasonable to use the Q as it is.

Could you please go over the patch and correct the Q and Q to Q (and I think in some places R is missing the quotes also). Then I can give this a positive review.

Changed 11 years ago by Rob Beezer

comment:22 Changed 11 years ago by Rob Beezer

Description: modified (diff)

Dear Martin,

Thanks for participating remotely in Sage Days 40.5. ;-)

"formatting" patch is an add-on, so you can see the changes. I did the "OUTPUT" section entirely with math quotes, other than the result objects. The rest is in code quotes.

Updated the output description, which needed work. Added some left kernel explanation for Dan Drake, who was also looking at this here at SD 40.5.

Thanks again, Rob

comment:23 Changed 11 years ago by Dan Drake

Reviewers: Martin RaumMartin Raum, Dan Drake
Status: needs_reviewpositive_review

This is a big improvement over what we previously had. Positive review.

Just for the patchbot: apply trac_10795-QR-decomposition-double-dense-v2.patch trac_10795-QR-decomposition-formatting.patch

comment:24 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-5.1.beta5
Resolution: fixed
Status: positive_reviewclosed

comment:25 Changed 11 years ago by Jeroen Demeyer

As reported on sage-devel, this gives doctest failures on OS X Lion, see #13140.

comment:26 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-5.1.beta5
Milestone: sage-5.1sage-5.2
Resolution: fixed
Status: closednew

Unmerging due to #13140.

comment:27 Changed 11 years ago by John Palmieri

Should this be marked "positive review" again, since #13140 is ready?

comment:28 in reply to:  27 Changed 11 years ago by Dan Drake

Status: newneeds_review

Replying to jhpalmieri:

Should this be marked "positive review" again, since #13140 is ready?

I think so. If that was the only reason this got unmerged, and that ticket is now fixed, this should again be ready.

comment:29 Changed 11 years ago by Dan Drake

Status: needs_reviewpositive_review

comment:30 Changed 11 years ago by Jeroen Demeyer

Merged in: sage-5.2.beta1
Resolution: fixed
Status: positive_reviewclosed
Note: See TracTickets for help on using tickets.