Opened 4 years ago

Closed 2 years ago

Last modified 23 months ago

#20439 closed defect (fixed)

eigenmatrix_right gives the conjugate of what it should

Reported by: tmonteil Owned by:
Priority: major Milestone: sage-8.1
Component: linear algebra Keywords: sd90
Cc: Merged in:
Authors: Alina Bucur, Renata Paramastri Reviewers: Kiran Kedlaya, Rebecca Lauren Miller, Travis Scrimshaw
Report Upstream: N/A Work issues:
Branch: a70e0db (Commits) Commit:
Dependencies: Stopgaps:

Description (last modified by tmonteil)

Here is the example reported at http://ask.sagemath.org/question/33084/bug-in-eigenmatrix-command/

sage: A = matrix(CDF, [[-2.53634347567,  2.04801738686, -0.0, -62.166145304], [ 0.7, -0.6, 0.0, 0.0], [0.547271128842, 0.0, -0.3015, -21.7532081652], [0.0, 0.0, 0.3, -0.4]])
sage: D, P = A.eigenmatrix_right()

According to the documentation, D,P should satisfy A*P == P*D, however:

sage: (A*P - P*D).norm()
7.454841195108627

The conjugate of P seems to be the correct answer though:

sage: P  = P.conjugate()
sage: (A*P - P*D).norm()
6.716506829378007e-15

The same issue holds for eigenmatrix_left

sage: sage: D, P = A.eigenmatrix_left()
sage: (P*A - D*P).norm()
8.009616737465649
sage: P = P.conjugate()
sage: (P*A - D*P).norm()
8.481237941055673e-15

Change History (20)

comment:1 Changed 4 years ago by tmonteil

  • Description modified (diff)

comment:2 Changed 2 years ago by kedlaya

It appears that both eigenmatrix_right and eigenmatrix_left are computed using eigenvectors_left, which for CDF is implemented in sage/matrix/matrix_double_dense.py. That in turn is calling scipy.linalg.eig, whose documentation

https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.eig.html

states that indeed for left eigenvectors, the reported "eigenvectors" are actually the conjugates thereof.

So I would propose changing eigenvectors_left to account for this.

Last edited 2 years ago by kedlaya (previous) (diff)

comment:3 Changed 2 years ago by kedlaya

  • Keywords sd90 added

comment:4 Changed 2 years ago by rparamastri

  • Branch set to u/rparamastri/eigenmatrix_right_gives_the_conjugate_of_what_it_should

comment:5 Changed 2 years ago by rparamastri

  • Authors set to Renata Paramastri
  • Commit set to cc7c2caef70be32f3b1e73830304c5c05b010932
  • Status changed from new to needs_review

Conjugated the eigenvectors returned by scipy.linalg.eig as suggested here.

comment:6 Changed 2 years ago by kedlaya

  • Status changed from needs_review to needs_work

All tests passed on k8s, so it looks like there were no unexpected side effects. (This is a bit surprising to me, but never mind.) I do have some comments on the doctests, though.

It is a Sage coding convention that lines should be no longer than 79 characters; anything longer should be split into multiple lines. In matrix2.pyx, the new doctest in eigenmatrix_left has a line that violates this convention.

The new doctests should reference this ticket. Sample text to illustrate the formatting:

    This test shows that :trac:`20439` has been resolved::
    ...

I would add a new doctest to left_eigenvectors where the substantive change was made, showing the effect of the change and again referencing this ticket. Repurposing the example from the top of the ticket would suffice.

comment:7 Changed 2 years ago by rlmiller

There should also be a blank line before and after examples, but not after the last example. Also I would reformat the examples so the lines are shorter.

comment:8 Changed 2 years ago by tscrim

To elaborate a bit more on the character limit @kedlaya is referring to, the usual way to split long lines is to do:

        sage: A = matrix(CDF, [[-2.53634347567, 2.04801738686, -0.0, -62.166145304],
        ....:                  [0.7, -0.6, 0.0, 0.0],
        ....:                  [0.547271128842, 0.0, -0.3015, -21.7532081652],
        ....:                  [0.0, 0.0, 0.3, -0.4]])

In some cases, you might have to go over 79 char/line for readability, but we try to stick to that limit as much as possible.

comment:9 Changed 2 years ago by kedlaya

  • Branch changed from u/rparamastri/eigenmatrix_right_gives_the_conjugate_of_what_it_should to u/kedlaya/eigenmatrix_right_gives_the_conjugate_of_what_it_should

comment:10 Changed 2 years ago by kedlaya

  • Authors changed from Renata Paramastri to Alina Bucur, Renata Paramastri
  • Commit changed from cc7c2caef70be32f3b1e73830304c5c05b010932 to 74b4876f93f5dbce3fcd7d182d1250bd21c27817
  • Status changed from needs_work to needs_review

Alina was having trouble with trac just now, so I've pushed her changes.


New commits:

74b4876Edited doctests

comment:11 Changed 2 years ago by alina

  • Reviewers set to ​Travis Scrimshaw, Rebecca Miller

Added ​Travis Scrimshaw, Rebecca Miller to reviewers.

comment:12 follow-up: Changed 2 years ago by tscrim

@kedlaya, you should add yourself to the reviewers list too.

I also have a few additional suggestions for the doctests:

-            sage: for i in range(A.dimensions()[0]):
-            ....:     ((A - spectrum[i][0]*identity_matrix(A.dimensions()[0]))*Matrix(spectrum[i][1]).transpose()).norm()<10^(-2)
+            sage: all(((A - spectrum[i][0]) * Matrix(spectrum[i][1]).transpose()).norm() < 10^(-2)
+            ....:     for i in range(A.nrows())

In fact, all of those < should have a space around them by PEP8 (it also is much easier to read quickly).

Two more spaces after the ....: in these:

            sage: for i in range(len(spectrum)):
            ....:   spectrum[i][1][0]=matrix(CDF, spectrum[i][1]).echelon_form()[0]

comment:13 Changed 2 years ago by git

  • Commit changed from 74b4876f93f5dbce3fcd7d182d1250bd21c27817 to 3465c52f1bae1391676f1074bbc570ced47a3bb1

Branch pushed to git repo; I updated commit sha1. New commits:

3465c52Edited doctests again

comment:14 Changed 2 years ago by kedlaya

  • Reviewers changed from ​Travis Scrimshaw, Rebecca Miller to Kiran Kedlaya, Rebecca Miller, ​Travis Scrimshaw

Added myself as reviewer and edited the doctests as proposed.

comment:15 in reply to: ↑ 12 ; follow-up: Changed 2 years ago by tscrim

Thank you. Almost done.

Replying to tscrim:

In fact, all of those < should have a space around them by PEP8 (it also is much easier to read quickly).

You missed a few:

-(P*A - D*P).norm()<10^(-2)
+(P*A - D*P).norm() < 10^(-2)

Once that is done, you can set a positive review on my behalf.

comment:16 Changed 2 years ago by tscrim

  • Milestone changed from sage-7.2 to sage-8.1

comment:17 Changed 2 years ago by git

  • Commit changed from 3465c52f1bae1391676f1074bbc570ced47a3bb1 to a70e0db8a58042d1562ce616b2cb8407a68f84be

Branch pushed to git repo; I updated commit sha1. New commits:

a70e0dbYet more edits to fix doctest formatting

comment:18 in reply to: ↑ 15 Changed 2 years ago by kedlaya

  • Status changed from needs_review to positive_review

Replying to tscrim:

Thank you. Almost done.

Replying to tscrim:

In fact, all of those < should have a space around them by PEP8 (it also is much easier to read quickly).

You missed a few:

-(P*A - D*P).norm()<10^(-2)
+(P*A - D*P).norm() < 10^(-2)

Once that is done, you can set a positive review on my behalf.

Done and done.

comment:19 Changed 2 years ago by vbraun

  • Branch changed from u/kedlaya/eigenmatrix_right_gives_the_conjugate_of_what_it_should to a70e0db8a58042d1562ce616b2cb8407a68f84be
  • Resolution set to fixed
  • Status changed from positive_review to closed

comment:20 Changed 23 months ago by jdemeyer

  • Commit a70e0db8a58042d1562ce616b2cb8407a68f84be deleted
  • Reviewers changed from Kiran Kedlaya, Rebecca Miller, ​Travis Scrimshaw to Kiran Kedlaya, Rebecca Lauren Miller, Travis Scrimshaw
Note: See TracTickets for help on using tickets.