Ticket #3794: eigenfunctions.patch

a/sage/matrix/matrix2.pyx

@@ -2424 +2424 @@
-                break
-        return f
-            
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-             (e, V)
-
-
+
+
+
+
+
+
+
+
+
+
-
-        The eigenspaces are returned sorted by the corresponding characteristic
-        polynomials, where polynomials are sorted in dictionary order starting
-        with constant terms. 
-
-        NOTE: 
-            Calling this method of a matrix over an inexact base ring will
-            raise a NotImplementedError.  If you want to force this anyway,
-            pass the option even_if_inexact=True.
-        
-
-        INPUT:
-            var -- variable name used to represent elements of
@@ -2455 +2469 @@
-        algorithm that is in dual_eigenvector in sage/modular/hecke/module.py.
-
-        EXAMPLES:
-
-
-
-
-
-
+left 
+
+
+
+
+_left
-            [
-            (0, Vector space of degree 3 and dimension 1 over Rational Field
-            User basis matrix:
@@ -2469 +2483 @@
-            User basis matrix:
-            [            1 1/15*a1 + 2/5 2/15*a1 - 1/5])
-            ]
-
+
+
+
+
+
+
+
+
+
+
-            sage: delta = e*v - v*A
-            sage: abs(abs(delta)) < 1e-10
-            True
@@ -2479 +2502 @@
-            [0 1 2]
-            [3 4 5]
-            [6 7 8]
-
+_left
-            [
-            (0, Vector space of degree 3 and dimension 1 over Rational Field
-            User basis matrix:
@@ -2489 +2512 @@
-            [            1 1/15*a1 + 2/5 2/15*a1 - 1/5])
-            ]
-
-
+left 
-        $T_2$ on level 43 modular symbols. 
-            sage: A = ModularSymbols(43).T(2).matrix(); A
-            [ 3  0  0  0  0  0 -1]
@@ -2505 +2528 @@
-            x^7 + x^6 - 12*x^5 - 16*x^4 + 36*x^3 + 52*x^2 - 32*x - 48
-            sage: factor(f)
-            (x - 3) * (x + 2)^2 * (x^2 - 2)^2
-(
+_left(algebraic_multiplicity=True
-            [
-            (3, Vector space of degree 7 and dimension 1 over Rational Field
-            User basis matrix:
-
+, 1
-            (-2, Vector space of degree 7 and dimension 2 over Rational Field
-            User basis matrix:
-            [ 0  1  0  1 -1  1 -1]
-
+, 2
-            (a2, Vector space of degree 7 and dimension 2 over Number Field in a2 with defining polynomial x^2 - 2
-            User basis matrix:
-            [      0       1       0      -1 -a2 - 1       1      -1]
-
-
-
-
+
+
+
+
+
-        of order 11:
-            sage: A = ModularSymbols(43, base_ring=GF(11), sign=1).T(2).matrix(); A
-            [ 3  9  0  0]
@@ -2531 +2555 @@
-            Finite Field of size 11
-            sage: A.charpoly()
-            x^4 + 10*x^3 + 3*x^2 + 2*x + 1
-
+_left
-            [
-            (9, Vector space of degree 4 and dimension 1 over Finite Field of size 11
-            User basis matrix:
@@ -2545 +2569 @@
-            ]
-    
-        TESTS:
+        Warnings are issued if the generic algorithm is used over inexact fields.  Garbage may result in these cases because of numerical precision issues.
-           sage: R=RealField(30)
-           sage: M=matrix(R,2,[2,1,1,1])
-
-
-
-
-
+
+
+
+
+
+
+
+
+
-           sage: R=ComplexField(30)
-           sage: N=matrix(R,2,[2,1,1,1])
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-
-        # minpoly is rarely implemented and is unreliable (leading to hangs) via linbox when implemented
-        # as of 2007-03-25.
@@ -2610 +2627 @@
-                A = self.change_ring(F) - alpha
-            W = A.kernel()
-            i = i + 1
-
+, e
-        V = Sequence(V, cr=True)
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-        return V
+
+
+
+    def eigenvectors_left(self):
+        r"""        
+        Compute the left eigenvectors of a matrix.
+
+        For each distinct eigenvalue, returns a list of the form
+        (e,V,n)
+        where e is the eigenvalue, V is a list of eigenvectors forming
+        a basis for the corresponding left eigenspace, and n is the
+        algebraic multiplicity of the eigenvalue.
+
+        EXAMPLES:
+        We compute the left eigenvectors of a $3\times 3$ rational matrix.
+            sage: A = matrix(QQ,3,3,range(9)); A
+            [0 1 2]
+            [3 4 5]
+            [6 7 8]
+            sage: es = A.eigenvectors_left(); es
+            [(0, [
+            (1, -2, 1)
+            ], 1),
+            (-1.348469228349535?, [(1, 0.3101020514433644?, -0.3797958971132713?)], 1),
+            (13.348469228349534?, [(1, 1.289897948556636?, 1.5797958971132712?)], 1)]
+            sage: eval, [evec], mult = es[0]
+            sage: delta = eval*evec - evec*A
+            sage: abs(abs(delta)) < 1e-10
+            True
+        """
+        x = self.fetch('eigenvectors_left')
+        if not x is None:
+            return x
+
+        if not self.base_ring().is_exact():
+            from warnings import warn
+            warn("Using generic algorithm for an inexact ring, which may result in garbarge from numerical precision issues.")
+
+        V = []
+        from sage.rings.qqbar import QQbar
+        from sage.categories.homset import hom
+        eigenspaces = self.eigenspaces_left(algebraic_multiplicity=True)
+        evec_list=[]
+        n = self._nrows
+        evec_eval_list = []
+        for ev in eigenspaces:
+            eigval = ev[0]
+            eigbasis = ev[1].basis()
+            eigmult = ev[2]
+            if hasattr(eigval, 'galois_conjugates'):
+                eigval_conj = eigval.galois_conjugates(QQbar)
+                for e in eigval_conj:
+                    m = hom(eigval.parent(), e.parent(), e)
+                    space = (e.parent())**n
+                    evec_list = [(space)([m(i) for i in v]) for v in eigbasis]
+                    evec_eval_list.append( (e, evec_list, eigmult))
+            else:
+                evec_eval_list.append((eigval, eigbasis, eigmult))
+        return evec_eval_list
+
+    left_eigenvectors = eigenvectors_left
+
+    def eigenvectors_right(self):
+        r"""        
+        Compute the right eigenvectors of a matrix.
+
+        For each distinct eigenvalue, returns a list of the form
+        (e,V,n)
+        where e is the eigenvalue, V is a list of eigenvectors forming
+        a basis for the corresponding right eigenspace, and n is the
+        algebraic multiplicity of the eigenvalue.
+
+        EXAMPLES:
+        We compute the right eigenvectors of a $3\times 3$ rational matrix.
+            sage: A = matrix(QQ,3,3,range(9)); A
+            [0 1 2]
+            [3 4 5]
+            [6 7 8]
+            sage: es = A.eigenvectors_right(); es
+            [(0, [
+            (1, -2, 1)
+            ], 1),
+            (-1.348469228349535?, [(1, 0.1303061543300932?, -0.7393876913398137?)], 1),
+            (13.348469228349534?, [(1, 3.069693845669907?, 5.139387691339814?)], 1)]
+            sage: eval, [evec], mult = es[0]
+            sage: delta = eval*evec - A*evec
+            sage: abs(abs(delta)) < 1e-10
+            True
+        """
+        return self.transpose().eigenvectors_left()
+
+    right_eigenvectors = eigenvectors_right
+
+    def eigenmatrix_left(self):
+        r"""
+        Return matrices D and P, where D is a diagonal matrix of
+        eigenvalues and P is the corresponding matrix where the rows
+        are corresponding eigenvectors so that P*self = D*P.
+
+        EXAMPLES:
+            sage: A = matrix(QQ,3,3,range(9)); A
+            [0 1 2]
+            [3 4 5]
+            [6 7 8]
+            sage: D, P = A.eigenmatrix_left()
+            sage: D
+            [                  0                   0                   0]
+            [                  0 -1.348469228349535?                   0]
+            [                  0                   0 13.348469228349534?]
+            sage: P
+            [                   1                   -2                    1]
+            [                   1  0.3101020514433644? -0.3797958971132713?]
+            [                   1   1.289897948556636?  1.5797958971132712?]
+            sage: P*A == D*P
+            True
+        """
+        from sage.misc.flatten import flatten
+        evecs = self.eigenvectors_left()
+        D = sage.matrix.constructor.diagonal_matrix(flatten([[e[0]]*e[2] for e in evecs]))
+        rows = []
+        for e in evecs:
+            rows.extend(e[1]+[e[1][0].parent().zero_vector()]*(e[2]-len(e[1])))
+        P = sage.matrix.constructor.matrix(rows)
+        return D,P
+
+    left_eigenmatrix = eigenmatrix_left
+
+    def eigenmatrix_right(self):
+        r"""
+        Return matrices D and P, where D is a diagonal matrix of
+        eigenvalues and P is the corresponding matrix where the columns
+        are corresponding eigenvectors so that self*P = P*D.
+
+        EXAMPLES:
+            sage: A = matrix(QQ,3,3,range(9)); A
+            [0 1 2]
+            [3 4 5]
+            [6 7 8]
+            sage: D, P = A.eigenmatrix_right()
+            sage: D
+            [                  0                   0                   0]
+            [                  0 -1.348469228349535?                   0]
+            [                  0                   0 13.348469228349534?]
+            sage: P
+            [                   1                    1                    1]
+            [                  -2  0.1303061543300932?   3.069693845669907?]
+            [                   1 -0.7393876913398137?   5.139387691339814?]
+            sage: A*P == P*D
+            True
+        """
+        D,P=self.transpose().eigenmatrix_left()
+        return D,P.transpose()
+
+    right_eigenmatrix = eigenmatrix_right
+
-
-
-    #####################################################################################