]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_element.py
eja: update todo, and rename "natural" to "matrix".
[sage.d.git] / mjo / eja / eja_element.py
index 926f2bf57a612c71eab53b693daa1e3f559ad6bc..6547668965a690b883a5ac59757ef1e625016604 100644 (file)
@@ -1,14 +1,7 @@
-# -*- coding: utf-8 -*-
-
 from sage.matrix.constructor import matrix
 from sage.modules.free_module import VectorSpace
 from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
 
-# TODO: make this unnecessary somehow.
-from sage.misc.lazy_import import lazy_import
-lazy_import('mjo.eja.eja_algebra', 'FiniteDimensionalEuclideanJordanAlgebra')
-lazy_import('mjo.eja.eja_element_subalgebra',
-            'FiniteDimensionalEuclideanJordanElementSubalgebra')
 from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
 from mjo.eja.eja_utils import _mat2vec
 
@@ -183,7 +176,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             True
 
         """
-        p = self.parent().characteristic_polynomial()
+        p = self.parent().characteristic_polynomial_of()
         return p(*self.to_vector())
 
 
@@ -346,6 +339,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         SETUP::
 
             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+            ....:                                  TrivialEJA,
             ....:                                  random_eja)
 
         EXAMPLES::
@@ -364,6 +358,17 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: x.det()
             -1
 
+        The determinant of the sole element in the rank-zero trivial
+        algebra is ``1``, by three paths of reasoning. First, its
+        characteristic polynomial is a constant ``1``, so the constant
+        term in that polynomial is ``1``. Second, the characteristic
+        polynomial evaluated at zero is again ``1``. And finally, the
+        (empty) product of its eigenvalues is likewise just unity::
+
+            sage: J = TrivialEJA()
+            sage: J.zero().det()
+            1
+
         TESTS:
 
         An element is invertible if and only if its determinant is
@@ -382,10 +387,16 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: x,y = J.random_elements(2)
             sage: (x*y).det() == x.det()*y.det()
             True
-
         """
         P = self.parent()
         r = P.rank()
+
+        if r == 0:
+            # Special case, since we don't get the a0=1
+            # coefficient when the rank of the algebra
+            # is zero.
+            return P.base_ring().one()
+
         p = P._charpoly_coefficients()[0]
         # The _charpoly_coeff function already adds the factor of -1
         # to ensure that _charpoly_coefficients()[0] is really what
@@ -420,11 +431,11 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: while not x.is_invertible():
             ....:     x = J.random_element()
             sage: x_vec = x.to_vector()
-            sage: x0 = x_vec[0]
+            sage: x0 = x_vec[:1]
             sage: x_bar = x_vec[1:]
-            sage: coeff = ~(x0^2 - x_bar.inner_product(x_bar))
-            sage: inv_vec = x_vec.parent()([x0] + (-x_bar).list())
-            sage: x_inverse = coeff*inv_vec
+            sage: coeff = x0.inner_product(x0) - x_bar.inner_product(x_bar)
+            sage: x_inverse = x_vec.parent()(x0.list() + (-x_bar).list())
+            sage: if not coeff.is_zero(): x_inverse = x_inverse/coeff
             sage: x.inverse() == J.from_vector(x_inverse)
             True
 
@@ -491,6 +502,14 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         if not self.is_invertible():
             raise ValueError("element is not invertible")
 
+        if self.parent()._charpoly_coefficients.is_in_cache():
+            # We can invert using our charpoly if it will be fast to
+            # compute. If the coefficients are cached, our rank had
+            # better be too!
+            r = self.parent().rank()
+            a = self.characteristic_polynomial().coefficients(sparse=False)
+            return (-1)**(r+1)*sum(a[i+1]*self**i for i in range(r))/self.det()
+
         return (~self.quadratic_representation())(self)
 
 
@@ -507,6 +526,11 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         whether or not the paren't algebra's zero element is a root
         of this element's minimal polynomial.
 
+        That is... unless the coefficients of our algebra's
+        "characteristic polynomial of" function are already cached!
+        In that case, we just use the determinant (which will be fast
+        as a result).
+
         Beware that we can't use the superclass method, because it
         relies on the algebra being associative.
 
@@ -537,6 +561,11 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             else:
                 return False
 
+        if self.parent()._charpoly_coefficients.is_in_cache():
+            # The determinant will be quicker than computing the minimal
+            # polynomial from scratch, most likely.
+            return (not self.det().is_zero())
+
         # In fact, we only need to know if the constant term is non-zero,
         # so we can pass in the field's zero element instead.
         zero = self.base_ring().zero()
@@ -778,8 +807,9 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
             sage: set_random_seed()
             sage: J = JordanSpinEJA.random_instance()
+            sage: n = J.dimension()
             sage: x = J.random_element()
-            sage: x == x.coefficient(0)*J.one() or x.degree() == 2
+            sage: x.degree() == min(n,2) or (x == x.coefficient(0)*J.one())
             True
 
         TESTS:
@@ -857,13 +887,18 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         TESTS:
 
         The minimal polynomial of the identity and zero elements are
-        always the same::
+        always the same, except in trivial algebras where the minimal
+        polynomial of the unit/zero element is ``1``::
 
             sage: set_random_seed()
-            sage: J = random_eja(nontrivial=True)
-            sage: J.one().minimal_polynomial()
+            sage: J = random_eja()
+            sage: mu = J.one().minimal_polynomial()
+            sage: t = mu.parent().gen()
+            sage: mu + int(J.is_trivial())*(t-2)
             t - 1
-            sage: J.zero().minimal_polynomial()
+            sage: mu = J.zero().minimal_polynomial()
+            sage: t = mu.parent().gen()
+            sage: mu + int(J.is_trivial())*(t-1)
             t
 
         The degree of an element is (by one definition) the degree
@@ -881,7 +916,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         two here so that said elements actually exist::
 
             sage: set_random_seed()
-            sage: n_max = max(2, JordanSpinEJA._max_test_case_size())
+            sage: n_max = max(2, JordanSpinEJA._max_random_instance_size())
             sage: n = ZZ.random_element(2, n_max)
             sage: J = JordanSpinEJA(n)
             sage: y = J.random_element()
@@ -907,7 +942,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         and in particular, a re-scaling of the basis::
 
             sage: set_random_seed()
-            sage: n_max = RealSymmetricEJA._max_test_case_size()
+            sage: n_max = RealSymmetricEJA._max_random_instance_size()
             sage: n = ZZ.random_element(1, n_max)
             sage: J1 = RealSymmetricEJA(n)
             sage: J2 = RealSymmetricEJA(n,normalize_basis=False)
@@ -931,20 +966,21 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
                 # in the "normal" case without us having to think about it.
                 return self.operator().minimal_polynomial()
 
-        A = self.subalgebra_generated_by()
+        A = self.subalgebra_generated_by(orthonormalize_basis=False)
         return A(self).operator().minimal_polynomial()
 
 
 
-    def natural_representation(self):
+    def to_matrix(self):
         """
-        Return a more-natural representation of this element.
+        Return an (often more natural) representation of this element as a
+        matrix.
 
-        Every finite-dimensional Euclidean Jordan Algebra is a
-        direct sum of five simple algebras, four of which comprise
-        Hermitian matrices. This method returns the original
-        "natural" representation of this element as a Hermitian
-        matrix, if it has one. If not, you get the usual representation.
+        Every finite-dimensional Euclidean Jordan Algebra is a direct
+        sum of five simple algebras, four of which comprise Hermitian
+        matrices. This method returns a "natural" matrix
+        representation of this element as either a Hermitian matrix or
+        column vector.
 
         SETUP::
 
@@ -956,7 +992,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: J = ComplexHermitianEJA(3)
             sage: J.one()
             e0 + e3 + e8
-            sage: J.one().natural_representation()
+            sage: J.one().to_matrix()
             [1 0 0 0 0 0]
             [0 1 0 0 0 0]
             [0 0 1 0 0 0]
@@ -969,7 +1005,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: J = QuaternionHermitianEJA(3)
             sage: J.one()
             e0 + e5 + e14
-            sage: J.one().natural_representation()
+            sage: J.one().to_matrix()
             [1 0 0 0 0 0 0 0 0 0 0 0]
             [0 1 0 0 0 0 0 0 0 0 0 0]
             [0 0 1 0 0 0 0 0 0 0 0 0]
@@ -982,11 +1018,14 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             [0 0 0 0 0 0 0 0 0 1 0 0]
             [0 0 0 0 0 0 0 0 0 0 1 0]
             [0 0 0 0 0 0 0 0 0 0 0 1]
-
         """
-        B = self.parent().natural_basis()
-        W = self.parent().natural_basis_space()
-        return W.linear_combination(zip(B,self.to_vector()))
+        B = self.parent().matrix_basis()
+        W = self.parent().matrix_space()
+
+        # This is just a manual "from_vector()", but of course
+        # matrix spaces aren't vector spaces in sage, so they
+        # don't have a from_vector() method.
+        return W.linear_combination( zip(B, self.to_vector()) )
 
 
     def norm(self):
@@ -1062,16 +1101,18 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: set_random_seed()
             sage: x = JordanSpinEJA.random_instance().random_element()
             sage: x_vec = x.to_vector()
+            sage: Q = matrix.identity(x.base_ring(), 0)
             sage: n = x_vec.degree()
-            sage: x0 = x_vec[0]
-            sage: x_bar = x_vec[1:]
-            sage: A = matrix(AA, 1, [x_vec.inner_product(x_vec)])
-            sage: B = 2*x0*x_bar.row()
-            sage: C = 2*x0*x_bar.column()
-            sage: D = matrix.identity(AA, n-1)
-            sage: D = (x0^2 - x_bar.inner_product(x_bar))*D
-            sage: D = D + 2*x_bar.tensor_product(x_bar)
-            sage: Q = matrix.block(2,2,[A,B,C,D])
+            sage: if n > 0:
+            ....:     x0 = x_vec[0]
+            ....:     x_bar = x_vec[1:]
+            ....:     A = matrix(x.base_ring(), 1, [x_vec.inner_product(x_vec)])
+            ....:     B = 2*x0*x_bar.row()
+            ....:     C = 2*x0*x_bar.column()
+            ....:     D = matrix.identity(x.base_ring(), n-1)
+            ....:     D = (x0^2 - x_bar.inner_product(x_bar))*D
+            ....:     D = D + 2*x_bar.tensor_product(x_bar)
+            ....:     Q = matrix.block(2,2,[A,B,C,D])
             sage: Q == x.quadratic_representation().matrix()
             True
 
@@ -1229,7 +1270,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: (J0, J5, J1) = J.peirce_decomposition(c1)
             sage: (f0, f1, f2) = J1.gens()
             sage: f0.spectral_decomposition()
-            [(0, 1.000000000000000?*f2), (1, 1.000000000000000?*f0)]
+            [(0, f2), (1, f0)]
 
         """
         A = self.subalgebra_generated_by(orthonormalize_basis=True)
@@ -1285,6 +1326,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             True
 
         """
+        from mjo.eja.eja_element_subalgebra import FiniteDimensionalEuclideanJordanElementSubalgebra
         return FiniteDimensionalEuclideanJordanElementSubalgebra(self, orthonormalize_basis)
 
 
@@ -1300,12 +1342,13 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         TESTS:
 
         Ensure that we can find an idempotent in a non-trivial algebra
-        where there are non-nilpotent elements::
+        where there are non-nilpotent elements, or that we get the dumb
+        solution in the trivial algebra::
 
             sage: set_random_seed()
-            sage: J = random_eja(nontrivial=True)
+            sage: J = random_eja()
             sage: x = J.random_element()
-            sage: while x.is_nilpotent():
+            sage: while x.is_nilpotent() and not J.is_trivial():
             ....:     x = J.random_element()
             sage: c = x.subalgebra_idempotent()
             sage: c^2 == c