]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_algebra.py
eja: partially fix the product_on_basis() in Cartesian products.
[sage.d.git] / mjo / eja / eja_algebra.py
index 631d695aa75208a65cb457a4ca11280f44823765..cbe5c08a344a95248e0b7c15e94de80b0c662816 100644 (file)
@@ -84,25 +84,51 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
     INPUT:
 
-      - basis -- a tuple of basis elements in "matrix form," which
-        must be the same form as the arguments to ``jordan_product``
-        and ``inner_product``. In reality, "matrix form" can be either
-        vectors, matrices, or a Cartesian product (ordered tuple)
-        of vectors or matrices. All of these would ideally be vector
-        spaces in sage with no special-casing needed; but in reality
-        we turn vectors into column-matrices and Cartesian products
-        `(a,b)` into column matrices `(a,b)^{T}` after converting
-        `a` and `b` themselves.
-
-      - jordan_product -- function of two ``basis`` elements (in
-        matrix form) that returns their jordan product, also in matrix
-        form; this will be applied to ``basis`` to compute a
-        multiplication table for the algebra.
-
-      - inner_product -- function of two ``basis`` elements (in matrix
-        form) that returns their inner product. This will be applied
-        to ``basis`` to compute an inner-product table (basically a
-        matrix) for this algebra.
+      - ``basis`` -- a tuple; a tuple of basis elements in "matrix
+        form," which must be the same form as the arguments to
+        ``jordan_product`` and ``inner_product``. In reality, "matrix
+        form" can be either vectors, matrices, or a Cartesian product
+        (ordered tuple) of vectors or matrices. All of these would
+        ideally be vector spaces in sage with no special-casing
+        needed; but in reality we turn vectors into column-matrices
+        and Cartesian products `(a,b)` into column matrices
+        `(a,b)^{T}` after converting `a` and `b` themselves.
+
+      - ``jordan_product`` -- a function; afunction of two ``basis``
+        elements (in matrix form) that returns their jordan product,
+        also in matrix form; this will be applied to ``basis`` to
+        compute a multiplication table for the algebra.
+
+      - ``inner_product`` -- a function; a function of two ``basis``
+        elements (in matrix form) that returns their inner
+        product. This will be applied to ``basis`` to compute an
+        inner-product table (basically a matrix) for this algebra.
+
+      - ``field`` -- a subfield of the reals (default: ``AA``); the scalar
+        field for the algebra.
+
+      - ``orthonormalize`` -- boolean (default: ``True``); whether or
+        not to orthonormalize the basis. Doing so is expensive and
+        generally rules out using the rationals as your ``field``, but
+        is required for spectral decompositions.
+
+    SETUP::
+
+        sage: from mjo.eja.eja_algebra import random_eja
+
+    TESTS:
+
+    We should compute that an element subalgebra is associative even
+    if we circumvent the element method::
+
+        sage: set_random_seed()
+        sage: J = random_eja(field=QQ,orthonormalize=False)
+        sage: x = J.random_element()
+        sage: A = x.subalgebra_generated_by(orthonormalize=False)
+        sage: basis = tuple(b.superalgebra_element() for b in A.basis())
+        sage: J.subalgebra(basis, orthonormalize=False).is_associative()
+        True
+
     """
     Element = FiniteDimensionalEJAElement
 
@@ -112,7 +138,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
                  inner_product,
                  field=AA,
                  orthonormalize=True,
-                 associative=False,
+                 associative=None,
                  cartesian_product=False,
                  check_field=True,
                  check_axioms=True,
@@ -169,7 +195,20 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
 
         category = MagmaticAlgebras(field).FiniteDimensional()
-        category = category.WithBasis().Unital()
+        category = category.WithBasis().Unital().Commutative()
+
+        if associative is None:
+            # We should figure it out. As with check_axioms, we have to do
+            # this without the help of the _jordan_product_is_associative()
+            # method because we need to know the category before we
+            # initialize the algebra.
+            associative = all( jordan_product(jordan_product(bi,bj),bk)
+                               ==
+                               jordan_product(bi,jordan_product(bj,bk))
+                               for bi in basis
+                               for bj in basis
+                               for bk in basis)
+
         if associative:
             # Element subalgebras can take advantage of this.
             category = category.Associative()
@@ -413,6 +452,16 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         """
         return "Associative" in self.category().axioms()
 
+    def _is_commutative(self):
+        r"""
+        Whether or not this algebra's multiplication table is commutative.
+
+        This method should of course always return ``True``, unless
+        this algebra was constructed with ``check_axioms=False`` and
+        passed an invalid multiplication table.
+        """
+        return all( x*y == y*x for x in self.gens() for y in self.gens() )
+
     def _is_jordanian(self):
         r"""
         Whether or not this algebra's multiplication table respects the
@@ -420,7 +469,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         We only check one arrangement of `x` and `y`, so for a
         ``True`` result to be truly true, you should also check
-        :meth:`is_commutative`. This method should of course always
+        :meth:`_is_commutative`. This method should of course always
         return ``True``, unless this algebra was constructed with
         ``check_axioms=False`` and passed an invalid multiplication table.
         """
@@ -430,6 +479,92 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
                     for i in range(self.dimension())
                     for j in range(self.dimension()) )
 
+    def _jordan_product_is_associative(self):
+        r"""
+        Return whether or not this algebra's Jordan product is
+        associative; that is, whether or not `x*(y*z) = (x*y)*z`
+        for all `x,y,x`.
+
+        This method should agree with :meth:`is_associative` unless
+        you lied about the value of the ``associative`` parameter
+        when you constructed the algebra.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (random_eja,
+            ....:                                  RealSymmetricEJA,
+            ....:                                  ComplexHermitianEJA,
+            ....:                                  QuaternionHermitianEJA)
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(4, orthonormalize=False)
+            sage: J._jordan_product_is_associative()
+            False
+            sage: x = sum(J.gens())
+            sage: A = x.subalgebra_generated_by()
+            sage: A._jordan_product_is_associative()
+            True
+
+        ::
+
+            sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
+            sage: J._jordan_product_is_associative()
+            False
+            sage: x = sum(J.gens())
+            sage: A = x.subalgebra_generated_by(orthonormalize=False)
+            sage: A._jordan_product_is_associative()
+            True
+
+        ::
+
+            sage: J = QuaternionHermitianEJA(2)
+            sage: J._jordan_product_is_associative()
+            False
+            sage: x = sum(J.gens())
+            sage: A = x.subalgebra_generated_by()
+            sage: A._jordan_product_is_associative()
+            True
+
+        TESTS:
+
+        The values we've presupplied to the constructors agree with
+        the computation::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: J.is_associative() == J._jordan_product_is_associative()
+            True
+
+        """
+        R = self.base_ring()
+
+        # Used to check whether or not something is zero.
+        epsilon = R.zero()
+        if not R.is_exact():
+            # I don't know of any examples that make this magnitude
+            # necessary because I don't know how to make an
+            # associative algebra when the element subalgebra
+            # construction is unreliable (as it is over RDF; we can't
+            # find the degree of an element because we can't compute
+            # the rank of a matrix). But even multiplication of floats
+            # is non-associative, so *some* epsilon is needed... let's
+            # just take the one from _inner_product_is_associative?
+            epsilon = 1e-15
+
+        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]
+                    diff = (x*y)*z - x*(y*z)
+
+                    if diff.norm() > epsilon:
+                        return False
+
+        return True
+
     def _inner_product_is_associative(self):
         r"""
         Return whether or not this algebra's inner product `B` is
@@ -439,11 +574,14 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
         this algebra was constructed with ``check_axioms=False`` and
         passed an invalid Jordan or inner-product.
         """
+        R = self.base_ring()
 
-        # Used to check whether or not something is zero in an inexact
-        # ring. This number is sufficient to allow the construction of
-        # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
-        epsilon = 1e-16
+        # Used to check whether or not something is zero.
+        epsilon = R.zero()
+        if not R.is_exact():
+            # This choice is sufficient to allow the construction of
+            # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
+            epsilon = 1e-15
 
         for i in range(self.dimension()):
             for j in range(self.dimension()):
@@ -453,12 +591,8 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
                     z = self.gens()[k]
                     diff = (x*y).inner_product(z) - x.inner_product(y*z)
 
-                    if self.base_ring().is_exact():
-                        if diff != 0:
-                            return False
-                    else:
-                        if diff.abs() > epsilon:
-                            return False
+                    if diff.abs() > epsilon:
+                        return False
 
         return True
 
@@ -472,7 +606,8 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         SETUP::
 
-            sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+            sage: from mjo.eja.eja_algebra import (random_eja,
+            ....:                                  JordanSpinEJA,
             ....:                                  HadamardEJA,
             ....:                                  RealSymmetricEJA)
 
@@ -505,18 +640,13 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         TESTS:
 
-        Ensure that we can convert any element of the two non-matrix
-        simple algebras (whose matrix representations are columns)
-        back and forth faithfully::
+        Ensure that we can convert any element back and forth
+        faithfully between its matrix and algebra representations::
 
             sage: set_random_seed()
-            sage: J = HadamardEJA.random_instance()
-            sage: x = J.random_element()
-            sage: J(x.to_vector().column()) == x
-            True
-            sage: J = JordanSpinEJA.random_instance()
+            sage: J = random_eja()
             sage: x = J.random_element()
-            sage: J(x.to_vector().column()) == x
+            sage: J(x.to_matrix()) == x
             True
 
         We cannot coerce elements between algebras just because their
@@ -532,7 +662,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
             Traceback (most recent call last):
             ...
             ValueError: not an element of this algebra
-
         """
         msg = "not an element of this algebra"
         if elt in self.base_ring():
@@ -800,7 +929,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
 
         # And to each subsequent row, prepend an entry that belongs to
         # the left-side "header column."
-        M += [ [self.gens()[i]] + [ self.product_on_basis(i,j)
+        M += [ [self.gens()[i]] + [ self.gens()[i]*self.gens()[j]
                                     for j in range(n) ]
                for i in range(n) ]
 
@@ -1471,6 +1600,13 @@ class RationalBasisEJA(FiniteDimensionalEJA):
             if not all( all(b_i in QQ for b_i in b.list()) for b in basis ):
                 raise TypeError("basis not rational")
 
+        super().__init__(basis,
+                         jordan_product,
+                         inner_product,
+                         field=field,
+                         check_field=check_field,
+                         **kwargs)
+
         self._rational_algebra = None
         if field is not QQ:
             # There's no point in constructing the extra algebra if this
@@ -1484,17 +1620,11 @@ class RationalBasisEJA(FiniteDimensionalEJA):
                                        jordan_product,
                                        inner_product,
                                        field=QQ,
+                                       associative=self.is_associative(),
                                        orthonormalize=False,
                                        check_field=False,
                                        check_axioms=False)
 
-        super().__init__(basis,
-                         jordan_product,
-                         inner_product,
-                         field=field,
-                         check_field=check_field,
-                         **kwargs)
-
     @cached_method
     def _charpoly_coefficients(self):
         r"""
@@ -1767,9 +1897,9 @@ class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA):
 
     In theory, our "field" can be any subfield of the reals::
 
-        sage: RealSymmetricEJA(2, field=RDF)
+        sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
         Euclidean Jordan algebra of dimension 3 over Real Double Field
-        sage: RealSymmetricEJA(2, field=RR)
+        sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
         Euclidean Jordan algebra of dimension 3 over Real Field with
         53 bits of precision
 
@@ -1858,10 +1988,15 @@ class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA):
         # if the user passes check_axioms=True.
         if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
 
-        super(RealSymmetricEJA, self).__init__(self._denormalized_basis(n),
-                                               self.jordan_product,
-                                               self.trace_inner_product,
-                                               **kwargs)
+        associative = False
+        if n <= 1:
+            associative = True
+
+        super().__init__(self._denormalized_basis(n),
+                         self.jordan_product,
+                         self.trace_inner_product,
+                         associative=associative,
+                         **kwargs)
 
         # TODO: this could be factored out somehow, but is left here
         # because the MatrixEJA is not presently a subclass of the
@@ -1951,7 +2086,7 @@ class ComplexMatrixEJA(MatrixEJA):
             True
 
         """
-        super(ComplexMatrixEJA,cls).real_embed(M)
+        super().real_embed(M)
         n = M.nrows()
 
         # We don't need any adjoined elements...
@@ -1998,7 +2133,7 @@ class ComplexMatrixEJA(MatrixEJA):
             True
 
         """
-        super(ComplexMatrixEJA,cls).real_unembed(M)
+        super().real_unembed(M)
         n = ZZ(M.nrows())
         d = cls.dimension_over_reals()
         F = cls.complex_extension(M.base_ring())
@@ -2035,9 +2170,9 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA):
 
     In theory, our "field" can be any subfield of the reals::
 
-        sage: ComplexHermitianEJA(2, field=RDF)
+        sage: ComplexHermitianEJA(2, field=RDF, check_axioms=True)
         Euclidean Jordan algebra of dimension 4 over Real Double Field
-        sage: ComplexHermitianEJA(2, field=RR)
+        sage: ComplexHermitianEJA(2, field=RR, check_axioms=True)
         Euclidean Jordan algebra of dimension 4 over Real Field with
         53 bits of precision
 
@@ -2146,10 +2281,15 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA):
         # if the user passes check_axioms=True.
         if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
 
-        super(ComplexHermitianEJA, self).__init__(self._denormalized_basis(n),
-                                                  self.jordan_product,
-                                                  self.trace_inner_product,
-                                                  **kwargs)
+        associative = False
+        if n <= 1:
+            associative = True
+
+        super().__init__(self._denormalized_basis(n),
+                         self.jordan_product,
+                         self.trace_inner_product,
+                         associative=associative,
+                         **kwargs)
         # TODO: this could be factored out somehow, but is left here
         # because the MatrixEJA is not presently a subclass of the
         # FDEJA class that defines rank() and one().
@@ -2232,7 +2372,7 @@ class QuaternionMatrixEJA(MatrixEJA):
             True
 
         """
-        super(QuaternionMatrixEJA,cls).real_embed(M)
+        super().real_embed(M)
         quaternions = M.base_ring()
         n = M.nrows()
 
@@ -2287,7 +2427,7 @@ class QuaternionMatrixEJA(MatrixEJA):
             True
 
         """
-        super(QuaternionMatrixEJA,cls).real_unembed(M)
+        super().real_unembed(M)
         n = ZZ(M.nrows())
         d = cls.dimension_over_reals()
 
@@ -2332,9 +2472,9 @@ class QuaternionHermitianEJA(ConcreteEJA, QuaternionMatrixEJA):
 
     In theory, our "field" can be any subfield of the reals::
 
-        sage: QuaternionHermitianEJA(2, field=RDF)
+        sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
         Euclidean Jordan algebra of dimension 6 over Real Double Field
-        sage: QuaternionHermitianEJA(2, field=RR)
+        sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
         Euclidean Jordan algebra of dimension 6 over Real Field with
         53 bits of precision
 
@@ -2452,10 +2592,16 @@ class QuaternionHermitianEJA(ConcreteEJA, QuaternionMatrixEJA):
         # if the user passes check_axioms=True.
         if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
 
-        super(QuaternionHermitianEJA, self).__init__(self._denormalized_basis(n),
-                                                     self.jordan_product,
-                                                     self.trace_inner_product,
-                                                     **kwargs)
+        associative = False
+        if n <= 1:
+            associative = True
+
+        super().__init__(self._denormalized_basis(n),
+                         self.jordan_product,
+                         self.trace_inner_product,
+                         associative=associative,
+                         **kwargs)
+
         # TODO: this could be factored out somehow, but is left here
         # because the MatrixEJA is not presently a subclass of the
         # FDEJA class that defines rank() and one().
@@ -2681,10 +2827,17 @@ class BilinearFormEJA(ConcreteEJA):
 
         n = B.nrows()
         column_basis = tuple( b.column() for b in FreeModule(ZZ, n).basis() )
-        super(BilinearFormEJA, self).__init__(column_basis,
-                                              jordan_product,
-                                              inner_product,
-                                              **kwargs)
+
+        # TODO: I haven't actually checked this, but it seems legit.
+        associative = False
+        if n <= 2:
+            associative = True
+
+        super().__init__(column_basis,
+                         jordan_product,
+                         inner_product,
+                         associative=associative,
+                         **kwargs)
 
         # The rank of this algebra is two, unless we're in a
         # one-dimensional ambient space (because the rank is bounded
@@ -2789,7 +2942,7 @@ class JordanSpinEJA(BilinearFormEJA):
 
         # But also don't pass check_field=False here, because the user
         # can pass in a field!
-        super(JordanSpinEJA, self).__init__(B, **kwargs)
+        super().__init__(B, **kwargs)
 
     @staticmethod
     def _max_random_instance_size():
@@ -2847,10 +3000,12 @@ class TrivialEJA(ConcreteEJA):
         if "orthonormalize" not in kwargs: kwargs["orthonormalize"] = False
         if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
 
-        super(TrivialEJA, self).__init__(basis,
-                                         jordan_product,
-                                         inner_product,
-                                         **kwargs)
+        super().__init__(basis,
+                         jordan_product,
+                         inner_product,
+                         associative=True,
+                         **kwargs)
+
         # The rank is zero using my definition, namely the dimension of the
         # largest subalgebra generated by any element.
         self.rank.set_cache(0)
@@ -3041,6 +3196,65 @@ class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct,
         self.one.set_cache(self._cartesian_product_of_elements(ones))
         self.rank.set_cache(sum(J.rank() for J in algebras))
 
+    def _monomial_to_generator(self, mon):
+        r"""
+        Convert a monomial index into a generator index.
+
+        SETUP::
+
+            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.
+        factor = mon[0]
+        idx_in_factor = 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):
+        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 (QuaternionHermitianEJA,
+            ....:                                  RealSymmetricEJA)
+
+        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
+
+        """
+        l = self._monomial_to_generator(i)
+        m = self._monomial_to_generator(j)
+        return FiniteDimensionalEJA.product_on_basis(self, l, m)
+
     def matrix_space(self):
         r"""
         Return the space that our matrix basis lives in as a Cartesian