Opened 11 years ago

Last modified 8 years ago

#12514 new defect

default behavior of matrix([...]).eigenvectors_right() should be sane

Reported by: Punarbasu Purkayastha Owned by: jason, was
Priority: major Milestone: sage-6.4
Component: linear algebra Keywords: matrix eigenvectors
Cc: Eviatar Bach Merged in:
Authors: Reviewers:
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Status badges

Description (last modified by Punarbasu Purkayastha)

I think the default behavior of eigenvectors_right() should be sane. Current behavior varies wildly.

  1. For example this fails:
    sage: matrix([ [2, 1], [1, 2.0]]).eigenvectors_right()
    ---------------------------------------------------------------------------
    NotImplementedError                       Traceback (most recent call last)
    
    /home/punarbasu/Installations/sage-5.0.beta2/devel/sagenb-github/<ipython console> in <module>()
    
    /home/punarbasu/Installations/sage-5.0.beta2/local/lib/python2.7/site-packages/sage/matrix/matrix2.so in sage.matrix.matrix2.Matrix.eigenvectors_right (sage/matrix/matrix2.c:27925)()
    
    /home/punarbasu/Installations/sage-5.0.beta2/local/lib/python2.7/site-packages/sage/matrix/matrix2.so in sage.matrix.matrix2.Matrix.eigenvectors_left (sage/matrix/matrix2.c:27260)()
    
    /home/punarbasu/Installations/sage-5.0.beta2/local/lib/python2.7/site-packages/sage/matrix/matrix2.so in sage.matrix.matrix2.Matrix.eigenspaces_left (sage/matrix/matrix2.c:24126)()
    
    NotImplementedError: eigenspaces cannot be computed reliably for inexact rings such as Real Field with 53 bits of precision,
    consult numerical or symbolic matrix classes for other options
    
  1. [Neglect this case. It is a typo. Keeping it here so that the references from the comments below make sense] This also fails but with a different error:
    sage: matrix( [ [1.0, 0], [0.1, 2.0] ] ).right_eigenvector()
    ---------------------------------------------------------------------------
    AttributeError                            Traceback (most recent call last)
    
    /home/punarbasu/Installations/sage-5.0.beta2/devel/sagenb-github/<ipython console> in <module>()
    
    /home/punarbasu/Installations/sage-5.0.beta2/local/lib/python2.7/site-packages/sage/structure/element.so in sage.structure.element.Element.__getattr__ (sage/structure/element.c:2921)()
    
    /home/punarbasu/Installations/sage-5.0.beta2/local/lib/python2.7/site-packages/sage/structure/parent.so in sage.structure.parent.getattr_from_other_class (sage/structure/parent.c:3302)()
    
    AttributeError: 'sage.matrix.matrix_generic_dense.Matrix_generic_dense' object has no attribute 'right_eigenvector'
    
  1. Using RDF gives a sane output:
    sage: matrix(RDF, [ [2, 1], [1, 2]]).eigenvectors_right()
    [(3.0, [(0.707106781187, 0.707106781187)], 1), (1.0, [(-0.707106781187, 0.707106781187)], 1)]
    
  1. But not using RDF in this case gives an output that appears "broken":
    sage: matrix([ [2, 1], [1, 2]]).eigenvectors_right()                           
    [(3, [
    (1, 1)
    ], 1), (1, [
    (1, -1)
    ], 1)]
    

In my opinion, by default all matrices should be made to belong to RDF or CDF if the field is not specified, so that most of the functionality can be availed.

The output in the last example appears "broken" because at some stage this function is called with cr=True, and this inserts a carriage return after every eigenvector:

Sequence(vecs, universe=V, check = False, immutable=True, cr=True)

Change History (15)

comment:1 Changed 11 years ago by Rob Beezer

(2): This just a typo in the method name used in the example.

matrix( [ [1.0, 0], [0.1, 2.0] ] ).right_eigenvectors()

will behave identically to (1). (Note plural.)

(4): I have never liked the use of Sequence with cr=True for the basis vectors of an eigenspace in these calls. IIRC, it would be easy to strip out, but a bunch of doctests would need reformatting. In fact, right now it is *harder* to doctest this sort of output for large instances since the output from Sequence does not include as many whitespace characters that can be used to rearrange the output. This might be a good question to pose to sage-devel with an eye to removing it?

(1) vs (3): For RDF/CDF matrices there is no notion of a geometric multiplicity for an eigenvalue. You will see that every eigenvalue comes back with a single eigenvector, implying that they have geometric multiplicty one. So changing the base ring to RDF/CDF should be a conscious decision, in my mind. The current "problem" with (1) is that absent a specified base ring, the choice is made to use RR, not RDF, and matrices over RR are not as capable as over RDF:

sage: A = matrix([ [2, 1], [1, 2.0]])
sage: A.base_ring() == RR
True

Maybe the error message in (1) could be improved to specifically suggest RDF, CDF or SR?

comment:2 Changed 11 years ago by Punarbasu Purkayastha

Description: modified (diff)

Indeed, (2) was a typo (now I! How quaint of me to have missed that. I will put a message on sage-devel regarding (4). Regarding (1) there are two things that could be done to improve the default behavior:

  1. Simply use RDF/CDF by default.
  1. Enclose the code in a try statement (which first uses whatever is the default) and return the computations from RDF (if the numbers are not complex) if the NotImplementedError exception occurs.

comment:3 in reply to:  2 ; Changed 11 years ago by Rob Beezer

Replying to ppurka:

I'd prefer not to default to RDF, perhaps just because somebody enters a decimal point. I'm "simulating" that situation below. The output is vastly different over RDF versus ZZ/QQ. I think it is important that a user understand the difference, rather than being a victim of behind-the-scenes choices about which base ring to use when.

sage: entries = [0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0]
sage: A = matrix(ZZ, 5, 5, entries)
sage: A.eigenvectors_right()       
[(4, [
(1, 1, 1, 1, 1)
], 1), (-1, [
(1, 0, 0, 0, -1),
(0, 1, 0, 0, -1),
(0, 0, 1, 0, -1),
(0, 0, 0, 1, -1)
], 4)]
sage: entries[0]=0.0  
sage: B = matrix(RDF, 5, 5, entries)
sage: B.eigenvectors_right()
[(-1.0, [(-0.894427191, 0.22360679775, 0.22360679775, 0.22360679775, 0.22360679775)], 1), (4.0, [(0.4472135955, 0.4472135955, 0.4472135955, 0.4472135955, 0.4472135955)], 1), (-1.0, [(-0.498272879122, -0.290659179488, 0.816613885228, -0.013840913309, -0.013840913309)], 1), (-1.0, [(-0.279399378114, 0.768348289812, 0.246802784, -0.36787584785, -0.36787584785)], 1), (-1.0, [(-0.226762156554, -0.132277924656, -0.0602899384488, -0.438226866039, 0.857556885698)], 1)]

comment:4 in reply to:  3 ; Changed 11 years ago by Punarbasu Purkayastha

Replying to rbeezer:

Replying to ppurka:

I'd prefer not to default to RDF, perhaps just because somebody enters a decimal point. I'm "simulating" that situation below. The output is vastly different over RDF versus ZZ/QQ.

Then how about the other option? Use RDF only if the code when run over RR returns with NotImplementedError? This would leave the ZZ or QQ case alone.

comment:5 in reply to:  4 Changed 11 years ago by Rob Beezer

Replying to ppurka:

Then how about the other option? Use RDF only if the code when run over RR returns with NotImplementedError? This would leave the ZZ or QQ case alone.

Converting from RR to RDF could mean a massive loss of precision (it seems possible, I checked). So there are two situations: (1) A matrix is entered with decimals and no base ring. User gets RR matrix. Eigenvalues are poorly implemented, eigenvectors quit.

(2) A matrix has RR entries to a high precision and the user likes it this way. Eigenvalues gives a warning (good), eigenvectors fails (better than a poor implementaion IMHO).

How about the following to satisfy both situations? Implement all the eigenstuff for RR/CC matrices by doing a conversion immediately, and then using output of the RDF/CDF routines. Maybe with a one-time warning like we do now for eigenvalues over RR/CC and good documentation. That way, if anybody ever gets excited about implementing numerical routines properly for RR/CC, the conversion will get replaced.

comment:6 Changed 11 years ago by Punarbasu Purkayastha

Yes. Doing RDF/CDF for only the eigen* stuff and at the same time giving a warning to the user is a good compromise, in my opinion.

comment:7 in reply to:  6 Changed 11 years ago by Rob Beezer

Replying to ppurka:

Yes. Doing RDF/CDF for only the eigen* stuff and at the same time giving a warning to the user is a good compromise, in my opinion.

Sounds good. Maybe the RDF/CDF changes should go on a new ticket, since this seems like a good place to house any dealings with Sequence(cr=True)? In any event, include me as a cc on anything you do and I'll try to look at it.

Rob

comment:8 Changed 10 years ago by Punarbasu Purkayastha

I was looking at this ticket again. It seems notoriously difficult to fix.

comment:9 Changed 9 years ago by Jeroen Demeyer

Milestone: sage-5.11sage-5.12

comment:10 Changed 9 years ago by For batch modifications

Milestone: sage-6.1sage-6.2

comment:11 Changed 9 years ago by Eviatar Bach

Cc: Eviatar Bach added

I don't know if this is related, but why does

sage: matrix(ZZ, [[6, -5], [1, 2]]).eigenvectors_right()
[(4 - 1*I, [(1, 0.4000000000000000? + 0.2000000000000000?*I)], 1),
 (4 + 1*I, [(1, 0.4000000000000000? - 0.2000000000000000?*I)], 1)]

return inexact eigenvectors? They could be computed exactly. Doing the same in SR returns exact results, but shouldn't matrices over ZZ or QQ try to get exact eigenvectors when possible?

comment:12 Changed 9 years ago by Punarbasu Purkayastha

Those are exact.

sage: E = matrix(ZZ, [[6, -5], [1, 2]]).eigenvectors_right()
sage: aa = E[0][1][0][1].real()             
sage: QQ(aa)
2/5

comment:13 Changed 9 years ago by Punarbasu Purkayastha

The output looks similar, though not exactly the same in terms of relative positioning of elements, for all of ZZ, QQ, and RDF. I think this ticket is no longer valid since 1. is being addressed in #13660 while 4. (which was the only thing left to be addressed) is similar for the rings I checked.

Example for RDF now:

sage: matrix(RDF, [ [2, 1], [1, 2]]).eigenvectors_right()
[(3.0, [(0.707106781187, 0.707106781187)], 1),
 (1.0, [(-0.707106781187, 0.707106781187)], 1)]

sage: entries = [0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0]
sage: A = matrix(RDF, 5, 5, entries)
sage: A.eigenvectors_right()
[(-1.0,
  [(-0.894427191, 0.22360679775, 0.22360679775, 0.22360679775, 0.22360679775)],
  1),
 (4.0,
  [(0.4472135955, 0.4472135955, 0.4472135955, 0.4472135955, 0.4472135955)],
  1),
 (-1.0,
  [(-0.19518001459, -0.439155032827, 0.862045064438, -0.113855008511, -0.113855008511)],
  1),
 (-1.0,
  [(-0.154320361522, 0.810181897988, 0.124313624559, -0.390087580513, -0.390087580513)],
  1),
 (-1.0,
  [(-0.093194221102, -0.209686997479, -0.209686997479, -0.365010699316, 0.877578915377)],
  1)]
Last edited 9 years ago by Punarbasu Purkayastha (previous) (diff)

comment:14 Changed 8 years ago by For batch modifications

Milestone: sage-6.2sage-6.3

comment:15 Changed 8 years ago by For batch modifications

Milestone: sage-6.3sage-6.4
Note: See TracTickets for help on using tickets.