eja: allow Cartesian products to be returned from random_eja().

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

diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py
index 031359fbfdab4ea8dc852f24ad7cfea8e54bfe4d..5bb4cee119e31cde1fd17accc7cc6e0bdee7f041 100644 (file)
--- a/mjo/eja/eja_algebra.py
+++ b/mjo/eja/eja_algebra.py
@@ -144,17 +144,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
                  check_axioms=True,
                  prefix='e'):
 
                  check_axioms=True,
                  prefix='e'):
 

-        # Keep track of whether or not the matrix basis consists of
-        # tuples, since we need special cases for them damned near
-        # everywhere.  This is INDEPENDENT of whether or not the
-        # algebra is a cartesian product, since a subalgebra of a
-        # cartesian product will have a basis of tuples, but will not
-        # in general itself be a cartesian product algebra.
-        self._matrix_basis_is_cartesian = False

         n = len(basis)
         n = len(basis)

-        if n > 0:
-            if hasattr(basis[0], 'cartesian_factors'):
-                self._matrix_basis_is_cartesian = True

 
         if check_field:
             if not field.is_subring(RR):
 
         if check_field:
             if not field.is_subring(RR):


@@ -163,20 +153,10 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
                 # we've specified a real embedding.
                 raise ValueError("scalar field is not real")
 
                 # we've specified a real embedding.
                 raise ValueError("scalar field is not real")
 

+        from mjo.eja.eja_utils import _change_ring

         # If the basis given to us wasn't over the field that it's
         # supposed to be over, fix that. Or, you know, crash.
         # If the basis given to us wasn't over the field that it's
         # supposed to be over, fix that. Or, you know, crash.

-        if not cartesian_product:
-            # The field for a cartesian product algebra comes from one
-            # of its factors and is the same for all factors, so
-            # there's no need to "reapply" it on product algebras.
-            if self._matrix_basis_is_cartesian:
-                # OK since if n == 0, the basis does not consist of tuples.
-                P = basis[0].parent()
-                basis = tuple( P(tuple(b_i.change_ring(field) for b_i in b))
-                               for b in basis )
-            else:
-                basis = tuple( b.change_ring(field) for b in basis )
-
+        basis = tuple( _change_ring(b, field) for b in basis )

 
         if check_axioms:
             # Check commutativity of the Jordan and inner-products.
 
         if check_axioms:
             # Check commutativity of the Jordan and inner-products.


@@ -213,7 +193,10 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
             # Element subalgebras can take advantage of this.
             category = category.Associative()
         if cartesian_product:
             # Element subalgebras can take advantage of this.
             category = category.Associative()
         if cartesian_product:

-            category = category.CartesianProducts()
+            # Use join() here because otherwise we only get the
+            # "Cartesian product of..." and not the things themselves.
+            category = category.join([category,
+                                      category.CartesianProducts()])

 
         # Call the superclass constructor so that we can use its from_vector()
         # method to build our multiplication table.
 
         # Call the superclass constructor so that we can use its from_vector()
         # method to build our multiplication table.


@@ -361,8 +344,8 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
             sage: if n > 0:
             ....:     i = ZZ.random_element(n)
             ....:     j = ZZ.random_element(n)
             sage: if n > 0:
             ....:     i = ZZ.random_element(n)
             ....:     j = ZZ.random_element(n)

-            ....:     ei = J.gens()[i]
-            ....:     ej = J.gens()[j]
+            ....:     ei = J.monomial(i)
+            ....:     ej = J.monomial(j)

             ....:     ei_ej = J.product_on_basis(i,j)
             sage: ei*ej == ei_ej
             True
             ....:     ei_ej = J.product_on_basis(i,j)
             sage: ei*ej == ei_ej
             True


@@ -473,9 +456,9 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         return ``True``, unless this algebra was constructed with
         ``check_axioms=False`` and passed an invalid multiplication table.
         """
         return ``True``, unless this algebra was constructed with
         ``check_axioms=False`` and passed an invalid multiplication table.
         """

-        return all( (self.gens()[i]**2)*(self.gens()[i]*self.gens()[j])
+        return all( (self.monomial(i)**2)*(self.monomial(i)*self.monomial(j))

                     ==
                     ==

-                    (self.gens()[i])*((self.gens()[i]**2)*self.gens()[j])
+                    (self.monomial(i))*((self.monomial(i)**2)*self.monomial(j))

                     for i in range(self.dimension())
                     for j in range(self.dimension()) )
 
                     for i in range(self.dimension())
                     for j in range(self.dimension()) )
 


@@ -555,9 +538,9 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         for i in range(self.dimension()):
             for j in range(self.dimension()):
                 for k in range(self.dimension()):
         for i in range(self.dimension()):
             for j in range(self.dimension()):
                 for k in range(self.dimension()):

-                    x = self.gens()[i]
-                    y = self.gens()[j]
-                    z = self.gens()[k]
+                    x = self.monomial(i)
+                    y = self.monomial(j)
+                    z = self.monomial(k)

                     diff = (x*y)*z - x*(y*z)
 
                     if diff.norm() > epsilon:
                     diff = (x*y)*z - x*(y*z)
 
                     if diff.norm() > epsilon:


@@ -586,9 +569,9 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         for i in range(self.dimension()):
             for j in range(self.dimension()):
                 for k in range(self.dimension()):
         for i in range(self.dimension()):
             for j in range(self.dimension()):
                 for k in range(self.dimension()):

-                    x = self.gens()[i]
-                    y = self.gens()[j]
-                    z = self.gens()[k]
+                    x = self.monomial(i)
+                    y = self.monomial(j)
+                    z = self.monomial(k)

                     diff = (x*y).inner_product(z) - x.inner_product(y*z)
 
                     if diff.abs() > epsilon:
                     diff = (x*y).inner_product(z) - x.inner_product(y*z)
 
                     if diff.abs() > epsilon:


@@ -636,7 +619,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
             sage: J2 = RealSymmetricEJA(2)
             sage: J = cartesian_product([J1,J2])
             sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
             sage: J2 = RealSymmetricEJA(2)
             sage: J = cartesian_product([J1,J2])
             sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )

-            e(0, 1) + e(1, 2)
+            e1 + e5

 
         TESTS:
 
 
         TESTS:
 


@@ -929,7 +912,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         # And to each subsequent row, prepend an entry that belongs to
         # the left-side "header column."
 
         # And to each subsequent row, prepend an entry that belongs to
         # the left-side "header column."

-        M += [ [self.gens()[i]] + [ self.gens()[i]*self.gens()[j]
+        M += [ [self.monomial(i)] + [ self.monomial(i)*self.monomial(j)

                                     for j in range(n) ]
                for i in range(n) ]
 
                                     for j in range(n) ]
                for i in range(n) ]
 


@@ -1431,7 +1414,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         def L_x_i_j(i,j):
             # From a result in my book, these are the entries of the
             # basis representation of L_x.
         def L_x_i_j(i,j):
             # From a result in my book, these are the entries of the
             # basis representation of L_x.

-            return sum( vars[k]*self.gens()[k].operator().matrix()[i,j]
+            return sum( vars[k]*self.monomial(k).operator().matrix()[i,j]

                         for k in range(n) )
 
         L_x = matrix(F, n, n, L_x_i_j)
                         for k in range(n) )
 
         L_x = matrix(F, n, n, L_x_i_j)


@@ -3018,8 +3001,7 @@ class TrivialEJA(ConcreteEJA):
         return cls(**kwargs)
 
 
         return cls(**kwargs)
 
 

-class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
-                          FiniteDimensionalEJA):
+class CartesianProductEJA(FiniteDimensionalEJA):

     r"""
     The external (orthogonal) direct sum of two or more Euclidean
     Jordan algebras. Every Euclidean Jordan algebra decomposes into an
     r"""
     The external (orthogonal) direct sum of two or more Euclidean
     Jordan algebras. Every Euclidean Jordan algebra decomposes into an


@@ -3115,6 +3097,33 @@ class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
         sage: CP2.is_associative()
         False
 
         sage: CP2.is_associative()
         False
 

+    Cartesian products of Cartesian products work::
+
+        sage: J1 = JordanSpinEJA(1)
+        sage: J2 = JordanSpinEJA(1)
+        sage: J3 = JordanSpinEJA(1)
+        sage: J = cartesian_product([J1,cartesian_product([J2,J3])])
+        sage: J.multiplication_table()
+        +----++----+----+----+
+        | *  || e0 | e1 | e2 |
+        +====++====+====+====+
+        | e0 || e0 | 0  | 0  |
+        +----++----+----+----+
+        | e1 || 0  | e1 | 0  |
+        +----++----+----+----+
+        | e2 || 0  | 0  | e2 |
+        +----++----+----+----+
+        sage: HadamardEJA(3).multiplication_table()
+        +----++----+----+----+
+        | *  || e0 | e1 | e2 |
+        +====++====+====+====+
+        | e0 || e0 | 0  | 0  |
+        +----++----+----+----+
+        | e1 || 0  | e1 | 0  |
+        +----++----+----+----+
+        | e2 || 0  | 0  | e2 |
+        +----++----+----+----+
+

     TESTS:
 
     All factors must share the same base field::
     TESTS:
 
     All factors must share the same base field::


@@ -3142,37 +3151,41 @@ class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
     Element = FiniteDimensionalEJAElement
 
 
     Element = FiniteDimensionalEJAElement
 
 

-    def __init__(self, algebras, **kwargs):
-        CombinatorialFreeModule_CartesianProduct.__init__(self,
-                                                          algebras,
-                                                          **kwargs)
-        field = algebras[0].base_ring()
-        if not all( J.base_ring() == field for J in algebras ):
+    def __init__(self, factors, **kwargs):
+        m = len(factors)
+        if m == 0:
+            return TrivialEJA()
+
+        self._sets = factors
+
+        field = factors[0].base_ring()
+        if not all( J.base_ring() == field for J in factors ):

             raise ValueError("all factors must share the same base field")
 
             raise ValueError("all factors must share the same base field")
 

-        associative = all( m.is_associative() for m in algebras )
+        associative = all( f.is_associative() for f in factors )

 
 

-        # The definition of matrix_space() and self.basis() relies
-        # only on the stuff in the CFM_CartesianProduct class, which
-        # we've already initialized.
-        Js = self.cartesian_factors()
-        m = len(Js)

         MS = self.matrix_space()
         MS = self.matrix_space()

-        basis = tuple(
-            MS(tuple( self.cartesian_projection(i)(b).to_matrix()
-                      for i in range(m) ))
-            for b in self.basis()
-        )
+        basis = []
+        zero = MS.zero()
+        for i in range(m):
+            for b in factors[i].matrix_basis():
+                z = list(zero)
+                z[i] = b
+                basis.append(z)
+
+        basis = tuple( MS(b) for b in basis )

 
         # Define jordan/inner products that operate on that matrix_basis.
         def jordan_product(x,y):
             return MS(tuple(
 
         # Define jordan/inner products that operate on that matrix_basis.
         def jordan_product(x,y):
             return MS(tuple(

-                (Js[i](x[i])*Js[i](y[i])).to_matrix() for i in range(m)
+                (factors[i](x[i])*factors[i](y[i])).to_matrix()
+                for i in range(m)

             ))
 
         def inner_product(x, y):
             return sum(
             ))
 
         def inner_product(x, y):
             return sum(

-                Js[i](x[i]).inner_product(Js[i](y[i])) for i in range(m)
+                factors[i](x[i]).inner_product(factors[i](y[i]))
+                for i in range(m)

             )
 
         # There's no need to check the field since it already came
             )
 
         # There's no need to check the field since it already came


@@ -3192,87 +3205,25 @@ class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
                                       check_field=False,
                                       check_axioms=False)
 
                                       check_field=False,
                                       check_axioms=False)
 

-        ones = tuple(J.one() for J in algebras)
-        self.one.set_cache(self._cartesian_product_of_elements(ones))
-        self.rank.set_cache(sum(J.rank() for J in algebras))
+        ones = tuple(J.one().to_matrix() for J in factors)
+        self.one.set_cache(self(ones))
+        self.rank.set_cache(sum(J.rank() for J in factors))

 
 

-    def _monomial_to_generator(self, mon):
-        r"""
-        Convert a monomial index into a generator index.
-
-        SETUP::
+    def cartesian_factors(self):
+        # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
+        return self._sets

 
 

-            sage: from mjo.eja.eja_algebra import random_eja()
-
-        TESTS::
-
-            sage: J1 = random_eja(field=QQ, orthonormalize=False)
-            sage: J2 = random_eja(field=QQ, orthonormalize=False)
-            sage: J = cartesian_product([J1,J2])
-            sage: all( J.monomial(m)
-            ....:      ==
-            ....:      J.gens()[J._monomial_to_generator(m)]
-            ....:      for m in J.basis().keys() )
-
-        """
-        # The superclass method indexes into a matrix, so we have to
-        # turn the tuples i and j into integers. This is easy enough
-        # given that the first coordinate of i and j corresponds to
-        # the factor, and the second coordinate corresponds to the
-        # index of the generator within that factor.
-        try:
-            factor = mon[0]
-        except TypeError: # 'int' object is not subscriptable
-            return mon
-        idx_in_factor = self._monomial_to_generator(mon[1])
-
-        offset = sum( f.dimension()
-                      for f in self.cartesian_factors()[:factor] )
-        return offset + idx_in_factor
-
-    def product_on_basis(self, i, j):
+    def cartesian_factor(self, i):

         r"""
         r"""

-        Return the product of the monomials indexed by ``i`` and ``j``.
-
-        This overrides the superclass method because here, both ``i``
-        and ``j`` will be ordered pairs.
-
-        SETUP::
-
-            sage: from mjo.eja.eja_algebra import (HadamardEJA,
-            ....:                                  JordanSpinEJA,
-            ....:                                  QuaternionHermitianEJA,
-            ....:                                  RealSymmetricEJA,)
-
-        EXAMPLES::
-
-            sage: J1 = JordanSpinEJA(2, field=QQ)
-            sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
-            sage: J3 = HadamardEJA(1, field=QQ)
-            sage: K1 = cartesian_product([J1,J2])
-            sage: K2 = cartesian_product([K1,J3])
-            sage: list(K2.basis())
-            [e(0, (0, 0)), e(0, (0, 1)), e(0, (1, 0)), e(0, (1, 1)),
-            e(0, (1, 2)), e(1, 0)]
-            sage: sage: g = K2.gens()
-            sage: (g[0] + 2*g[3]) * (g[1] - 4*g[2])
-            e(0, (0, 1)) - 4*e(0, (1, 1))
-
-        TESTS::
-
-            sage: J1 = RealSymmetricEJA(1,field=QQ)
-            sage: J2 = QuaternionHermitianEJA(1,field=QQ)
-            sage: J = cartesian_product([J1,J2])
-            sage: x = sum(J.gens())
-            sage: x == J.one()
-            True
-            sage: x*x == x
-            True
-
+        Return the ``i``th factor of this algebra.

         """
         """

-        l = self._monomial_to_generator(i)
-        m = self._monomial_to_generator(j)
-        return FiniteDimensionalEJA.product_on_basis(self, l, m)
+        return self._sets[i]
+
+    def _repr_(self):
+        # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
+        from sage.categories.cartesian_product import cartesian_product
+        return cartesian_product.symbol.join("%s" % factor
+                                             for factor in self._sets)

 
     def matrix_space(self):
         r"""
 
     def matrix_space(self):
         r"""


@@ -3369,9 +3320,12 @@ class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
             True
 
         """
             True
 
         """

-        Ji = self.cartesian_factors()[i]
-        # Requires the fix on Trac 31421/31422 to work!
-        Pi = super().cartesian_projection(i)
+        offset = sum( self.cartesian_factor(k).dimension()
+                      for k in range(i) )
+        Ji = self.cartesian_factor(i)
+        Pi = self._module_morphism(lambda j: Ji.monomial(j - offset),
+                                   codomain=Ji)
+

         return FiniteDimensionalEJAOperator(self,Ji,Pi.matrix())
 
     @cached_method
         return FiniteDimensionalEJAOperator(self,Ji,Pi.matrix())
 
     @cached_method


@@ -3477,9 +3431,11 @@ class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
             True
 
         """
             True
 
         """

-        Ji = self.cartesian_factors()[i]
-        # Requires the fix on Trac 31421/31422 to work!
-        Ei = super().cartesian_embedding(i)
+        offset = sum( self.cartesian_factor(k).dimension()
+                      for k in range(i) )
+        Ji = self.cartesian_factor(i)
+        Ei = Ji._module_morphism(lambda j: self.monomial(j + offset),
+                                 codomain=self)

         return FiniteDimensionalEJAOperator(Ji,self,Ei.matrix())
 
 
         return FiniteDimensionalEJAOperator(Ji,self,Ei.matrix())
 
 


@@ -3524,17 +3480,15 @@ class RationalBasisCartesianProductEJA(CartesianProductEJA,
 
 RationalBasisEJA.CartesianProduct = RationalBasisCartesianProductEJA
 
 
 RationalBasisEJA.CartesianProduct = RationalBasisCartesianProductEJA
 

-random_eja = ConcreteEJA.random_instance
-
-# def random_eja(*args, **kwargs):
-#     J1 = ConcreteEJA.random_instance(*args, **kwargs)
-
-#     # This might make Cartesian products appear roughly as often as
-#     # any other ConcreteEJA.
-#     if ZZ.random_element(len(ConcreteEJA.__subclasses__()) + 1) == 0:
-#         # Use random_eja() again so we can get more than two factors.
-#         J2 = random_eja(*args, **kwargs)
-#         J = cartesian_product([J1,J2])
-#         return J
-#     else:
-#         return J1
+def random_eja(*args, **kwargs):
+    J1 = ConcreteEJA.random_instance(*args, **kwargs)
+
+    # This might make Cartesian products appear roughly as often as
+    # any other ConcreteEJA.
+    if ZZ.random_element(len(ConcreteEJA.__subclasses__()) + 1) == 0:
+        # Use random_eja() again so we can get more than two factors.
+        J2 = random_eja(*args, **kwargs)
+        J = cartesian_product([J1,J2])
+        return J
+    else:
+        return J1