eja: drop a pointless "solve" in the EJA charpoly system.

[sage.d.git] / mjo / eja / eja_algebra.py

diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py
index ae41b8c6e12415d4b940958ca3798749a4e8c7bf..c35bf256453cfa9a3758c034a47a043745eb0c2b 100644 (file)
--- a/mjo/eja/eja_algebra.py
+++ b/mjo/eja/eja_algebra.py
@@ -1101,6 +1101,21 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         r"""
         The `r` polynomial coefficients of the "characteristic polynomial
         of" function.
         r"""
         The `r` polynomial coefficients of the "characteristic polynomial
         of" function.

+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import random_eja
+
+        TESTS:
+
+        The theory shows that these are all homogeneous polynomials of
+        a known degree::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
+            True
+

         """
         n = self.dimension()
         R = self.coordinate_polynomial_ring()
         """
         n = self.dimension()
         R = self.coordinate_polynomial_ring()


@@ -1136,10 +1151,17 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         # The theory says that only the first "r" coefficients are
         # nonzero, and they actually live in the original polynomial
 
         # The theory says that only the first "r" coefficients are
         # nonzero, and they actually live in the original polynomial

-        # ring and not the fraction field. We negate them because
-        # in the actual characteristic polynomial, they get moved
-        # to the other side where x^r lives.
-        return -A_rref.solve_right(E*b).change_ring(R)[:r]
+        # ring and not the fraction field. We negate them because in
+        # the actual characteristic polynomial, they get moved to the
+        # other side where x^r lives. We don't bother to trim A_rref
+        # down to a square matrix and solve the resulting system,
+        # because the upper-left r-by-r portion of A_rref is
+        # guaranteed to be the identity matrix, so e.g.
+        #
+        #   A_rref.solve_right(Y)
+        #
+        # would just be returning Y.
+        return (-E*b)[:r].change_ring(R)

 
     @cached_method
     def rank(self):
 
     @cached_method
     def rank(self):


@@ -1200,7 +1222,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
             sage: set_random_seed()    # long time
             sage: J = random_eja()     # long time
 
             sage: set_random_seed()    # long time
             sage: J = random_eja()     # long time

-            sage: caches = J.rank()    # long time
+            sage: cached = J.rank()    # long time

             sage: J.rank.clear_cache() # long time
             sage: J.rank() == cached   # long time
             True
             sage: J.rank.clear_cache() # long time
             sage: J.rank() == cached   # long time
             True


@@ -2652,7 +2674,8 @@ class TrivialEJA(ConcreteEJA):
         # inappropriate for us.
         return cls(**kwargs)
 
         # inappropriate for us.
         return cls(**kwargs)
 

-class DirectSumEJA(ConcreteEJA):
+
+class DirectSumEJA(FiniteDimensionalEJA):

     r"""
     The external (orthogonal) direct sum of two other Euclidean Jordan
     algebras. Essentially the Cartesian product of its two factors.
     r"""
     The external (orthogonal) direct sum of two other Euclidean Jordan
     algebras. Essentially the Cartesian product of its two factors.


@@ -2676,6 +2699,10 @@ class DirectSumEJA(ConcreteEJA):
         8
         sage: J.rank()
         5
         8
         sage: J.rank()
         5

+        sage: J.matrix_space()
+        The Cartesian product of (Full MatrixSpace of 2 by 1 dense matrices
+        over Algebraic Real Field, Full MatrixSpace of 3 by 3 dense matrices
+        over Algebraic Real Field)

 
     TESTS:
 
 
     TESTS:
 


@@ -2695,21 +2722,50 @@ class DirectSumEJA(ConcreteEJA):
         if J1.base_ring() != J2.base_ring():
             raise ValueError("algebras must share the same base field")
         field = J1.base_ring()
         if J1.base_ring() != J2.base_ring():
             raise ValueError("algebras must share the same base field")
         field = J1.base_ring()

-        self._factors = (J1, J2)
-        basis = tuple( (a,b) for a in J1.basis() for b in J2.basis() )

 
 

-        def jordan_product(x,y):
-            return (x[0]*y[0], x[1]*y[1])
-
-        def inner_product(x,y):
-            return x[0].inner_product(y[0]) + x[1].inner_product(y[1])
+        M = J1.matrix_space().cartesian_product(J2.matrix_space())
+        self._cartprod_algebra = J1.cartesian_product(J2)

 
 

-        super().__init__(basis, jordan_product, inner_product)
+        self._matrix_basis = tuple( [M((a,0)) for a in J1.matrix_basis()] +
+                                    [M((0,b)) for b in J2.matrix_basis()] )

 
 

+        n = len(self._matrix_basis)
+        self._sets = None
+        CombinatorialFreeModule.__init__(
+                         self,
+                         field,
+                         range(n),
+                         category=self._cartprod_algebra.category(),
+                         bracket=False,
+                         **kwargs)

         self.rank.set_cache(J1.rank() + J2.rank())
 
 
         self.rank.set_cache(J1.rank() + J2.rank())
 
 

-    def factors(self):
+
+    def product(self,x,y):
+        r"""
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+            ....:                                  ComplexHermitianEJA,
+            ....:                                  DirectSumEJA)
+
+        TESTS::
+
+            sage: set_random_seed()
+            sage: J1 = JordanSpinEJA(3, field=QQ)
+            sage: J2 = ComplexHermitianEJA(2, field=QQ, orthonormalize=False)
+            sage: J = DirectSumEJA(J1,J2)
+            sage: J.random_element()*J.random_element() in J
+            True
+
+        """
+        xv = self._cartprod_algebra.from_vector(x.to_vector())
+        yv = self._cartprod_algebra.from_vector(y.to_vector())
+        return self.from_vector((xv*yv).to_vector())
+
+
+    def cartesian_factors(self):

         r"""
         Return the pair of this algebra's factors.
 
         r"""
         Return the pair of this algebra's factors.
 


@@ -2724,12 +2780,13 @@ class DirectSumEJA(ConcreteEJA):
             sage: J1 = HadamardEJA(2, field=QQ)
             sage: J2 = JordanSpinEJA(3, field=QQ)
             sage: J = DirectSumEJA(J1,J2)
             sage: J1 = HadamardEJA(2, field=QQ)
             sage: J2 = JordanSpinEJA(3, field=QQ)
             sage: J = DirectSumEJA(J1,J2)

-            sage: J.factors()
+            sage: J.cartesian_factors()

             (Euclidean Jordan algebra of dimension 2 over Rational Field,
              Euclidean Jordan algebra of dimension 3 over Rational Field)
 
         """
             (Euclidean Jordan algebra of dimension 2 over Rational Field,
              Euclidean Jordan algebra of dimension 3 over Rational Field)
 
         """

-        return self._factors
+        return self._cartprod_algebra.cartesian_factors()
+

 
 #     def projections(self):
 #         r"""
 
 #     def projections(self):
 #         r"""