]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_algebra.py
eja: introduce an intermediate class for "concrete" algebras.
[sage.d.git] / mjo / eja / eja_algebra.py
index 53f61511ed480a242df760d6a1808d91ce2f4607..3f2127d7df9b811a599ba1da01b7b100ab0a5573 100644 (file)
@@ -5,7 +5,7 @@ are used in optimization, and have some additional nice methods beyond
 what can be supported in a general Jordan Algebra.
 """
 
-from itertools import repeat
+from itertools import izip, repeat
 
 from sage.algebras.quatalg.quaternion_algebra import QuaternionAlgebra
 from sage.categories.magmatic_algebras import MagmaticAlgebras
@@ -153,26 +153,6 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule):
         return self.from_vector(coords)
 
 
-    @staticmethod
-    def _max_test_case_size():
-        """
-        Return an integer "size" that is an upper bound on the size of
-        this algebra when it is used in a random test
-        case. Unfortunately, the term "size" is quite vague -- when
-        dealing with `R^n` under either the Hadamard or Jordan spin
-        product, the "size" refers to the dimension `n`. When dealing
-        with a matrix algebra (real symmetric or complex/quaternion
-        Hermitian), it refers to the size of the matrix, which is
-        far less than the dimension of the underlying vector space.
-
-        We default to five in this class, which is safe in `R^n`. The
-        matrix algebra subclasses (or any class where the "size" is
-        interpreted to be far less than the dimension) should override
-        with a smaller number.
-        """
-        return 5
-
-
     def _repr_(self):
         """
         Return a string representation of ``self``.
@@ -660,27 +640,6 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule):
         """
         return  tuple( self.random_element() for idx in xrange(count) )
 
-    @classmethod
-    def random_instance(cls, field=QQ, **kwargs):
-        """
-        Return a random instance of this type of algebra.
-
-        In subclasses for algebras that we know how to construct, this
-        is a shortcut for constructing test cases and examples.
-        """
-        if cls is FiniteDimensionalEuclideanJordanAlgebra:
-            # Red flag! But in theory we could do this I guess. The
-            # only finite-dimensional exceptional EJA is the
-            # octononions. So, we could just create an EJA from an
-            # associative matrix algebra (generated by a subset of
-            # elements) with the symmetric product. Or, we could punt
-            # to random_eja() here, override it in our subclasses, and
-            # not worry about it.
-            raise NotImplementedError
-
-        n = ZZ.random_element(cls._max_test_case_size()) + 1
-        return cls(n, field, **kwargs)
-
 
     def rank(self):
         """
@@ -757,7 +716,57 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule):
     Element = FiniteDimensionalEuclideanJordanAlgebraElement
 
 
-class RealCartesianProductEJA(FiniteDimensionalEuclideanJordanAlgebra):
+class KnownRankEJA(object):
+    """
+    A class for algebras that we actually know we can construct.  The
+    main issue is that, for most of our methods to make sense, we need
+    to know the rank of our algebra. Thus we can't simply generate a
+    "random" algebra, or even check that a given basis and product
+    satisfy the axioms; because even if everything looks OK, we wouldn't
+    know the rank we need to actuallty build the thing.
+
+    Not really a subclass of FDEJA because doing that causes method
+    resolution errors, e.g.
+
+      TypeError: Error when calling the metaclass bases
+      Cannot create a consistent method resolution
+      order (MRO) for bases FiniteDimensionalEuclideanJordanAlgebra,
+      KnownRankEJA
+
+    """
+    @staticmethod
+    def _max_test_case_size():
+        """
+        Return an integer "size" that is an upper bound on the size of
+        this algebra when it is used in a random test
+        case. Unfortunately, the term "size" is quite vague -- when
+        dealing with `R^n` under either the Hadamard or Jordan spin
+        product, the "size" refers to the dimension `n`. When dealing
+        with a matrix algebra (real symmetric or complex/quaternion
+        Hermitian), it refers to the size of the matrix, which is
+        far less than the dimension of the underlying vector space.
+
+        We default to five in this class, which is safe in `R^n`. The
+        matrix algebra subclasses (or any class where the "size" is
+        interpreted to be far less than the dimension) should override
+        with a smaller number.
+        """
+        return 5
+
+    @classmethod
+    def random_instance(cls, field=QQ, **kwargs):
+        """
+        Return a random instance of this type of algebra.
+
+        Beware, this will crash for "most instances" because the
+        constructor below looks wrong.
+        """
+        n = ZZ.random_element(cls._max_test_case_size()) + 1
+        return cls(n, field, **kwargs)
+
+
+class RealCartesianProductEJA(FiniteDimensionalEuclideanJordanAlgebra,
+                              KnownRankEJA):
     """
     Return the Euclidean Jordan Algebra corresponding to the set
     `R^n` under the Hadamard product.
@@ -866,11 +875,7 @@ def random_eja():
         Euclidean Jordan algebra of dimension...
 
     """
-    classname = choice([RealCartesianProductEJA,
-                        JordanSpinEJA,
-                        RealSymmetricEJA,
-                        ComplexHermitianEJA,
-                        QuaternionHermitianEJA])
+    classname = choice(KnownRankEJA.__subclasses__())
     return classname.random_instance()
 
 
@@ -885,17 +890,20 @@ class MatrixEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra):
         # field can have dimension 4 (quaternions) too.
         return 2
 
-    @classmethod
-    def _denormalized_basis(cls, n, field):
-        raise NotImplementedError
-
-    def __init__(self, n, field=QQ, normalize_basis=True, **kwargs):
-        S = self._denormalized_basis(n, field)
-
+    def __init__(self, field, basis, rank, normalize_basis=True, **kwargs):
+        """
+        Compared to the superclass constructor, we take a basis instead of
+        a multiplication table because the latter can be computed in terms
+        of the former when the product is known (like it is here).
+        """
         # Used in this class's fast _charpoly_coeff() override.
         self._basis_normalizers = None
 
-        if n > 1 and normalize_basis:
+        # We're going to loop through this a few times, so now's a good
+        # time to ensure that it isn't a generator expression.
+        basis = tuple(basis)
+
+        if rank > 1 and normalize_basis:
             # We'll need sqrt(2) to normalize the basis, and this
             # winds up in the multiplication table, so the whole
             # algebra needs to be over the field extension.
@@ -904,18 +912,18 @@ class MatrixEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra):
             p = z**2 - 2
             if p.is_irreducible():
                 field = NumberField(p, 'sqrt2', embedding=RLF(2).sqrt())
-                S = [ s.change_ring(field) for s in S ]
+                basis = tuple( s.change_ring(field) for s in basis )
             self._basis_normalizers = tuple(
-                ~(self.natural_inner_product(s,s).sqrt()) for s in S )
-            S = tuple( s*c for (s,c) in zip(S,self._basis_normalizers) )
+                ~(self.natural_inner_product(s,s).sqrt()) for s in basis )
+            basis = tuple(s*c for (s,c) in izip(basis,self._basis_normalizers))
 
-        Qs = self.multiplication_table_from_matrix_basis(S)
+        Qs = self.multiplication_table_from_matrix_basis(basis)
 
         fdeja = super(MatrixEuclideanJordanAlgebra, self)
         return fdeja.__init__(field,
                               Qs,
-                              rank=n,
-                              natural_basis=S,
+                              rank=rank,
+                              natural_basis=basis,
                               **kwargs)
 
 
@@ -930,15 +938,17 @@ class MatrixEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra):
             # with had entries in a nice field.
             return super(MatrixEuclideanJordanAlgebra, self)._charpoly_coeff(i)
         else:
-            # If we didn't unembed first, this number would be wrong
-            # by a power-of-two factor for complex/quaternion matrices.
-            n = self.real_unembed(self.natural_basis_space().zero()).nrows()
-            field = self.base_ring().base_ring() # yeeeeaaaahhh
-            J = self.__class__(n, field, False)
+            basis = ( (b/n) for (b,n) in izip(self.natural_basis(),
+                                              self._basis_normalizers) )
+            field = self.base_ring().base_ring() # yeeeaahhhhhhh
+            J = MatrixEuclideanJordanAlgebra(field,
+                                             basis,
+                                             self.rank(),
+                                             normalize_basis=False)
             (_,x,_,_) = J._charpoly_matrix_system()
             p = J._charpoly_coeff(i)
             # p might be missing some vars, have to substitute "optionally"
-            pairs = zip(x.base_ring().gens(), self._basis_normalizers)
+            pairs = izip(x.base_ring().gens(), self._basis_normalizers)
             substitutions = { v: v*c for (v,c) in pairs }
             return p.subs(substitutions)
 
@@ -1022,29 +1032,21 @@ class RealMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra):
     @staticmethod
     def real_embed(M):
         """
-        Embed the matrix ``M`` into a space of real matrices.
-
-        The matrix ``M`` can have entries in any field at the moment:
-        the real numbers, complex numbers, or quaternions. And although
-        they are not a field, we can probably support octonions at some
-        point, too. This function returns a real matrix that "acts like"
-        the original with respect to matrix multiplication; i.e.
-
-          real_embed(M*N) = real_embed(M)*real_embed(N)
-
+        The identity function, for embedding real matrices into real
+        matrices.
         """
         return M
 
-
     @staticmethod
     def real_unembed(M):
         """
-        The inverse of :meth:`real_embed`.
+        The identity function, for unembedding real matrices from real
+        matrices.
         """
         return M
 
 
-class RealSymmetricEJA(RealMatrixEuclideanJordanAlgebra):
+class RealSymmetricEJA(RealMatrixEuclideanJordanAlgebra, KnownRankEJA):
     """
     The rank-n simple EJA consisting of real symmetric n-by-n
     matrices, the usual symmetric Jordan product, and the trace inner
@@ -1144,7 +1146,7 @@ class RealSymmetricEJA(RealMatrixEuclideanJordanAlgebra):
                 else:
                     Sij = Eij + Eij.transpose()
                 S.append(Sij)
-        return tuple(S)
+        return S
 
 
     @staticmethod
@@ -1152,6 +1154,10 @@ class RealSymmetricEJA(RealMatrixEuclideanJordanAlgebra):
         return 4 # Dimension 10
 
 
+    def __init__(self, n, field=QQ, **kwargs):
+        basis = self._denormalized_basis(n, field)
+        super(RealSymmetricEJA, self).__init__(field, basis, n, **kwargs)
+
 
 class ComplexMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra):
     @staticmethod
@@ -1303,7 +1309,7 @@ class ComplexMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra):
         return RealMatrixEuclideanJordanAlgebra.natural_inner_product(X,Y)/2
 
 
-class ComplexHermitianEJA(ComplexMatrixEuclideanJordanAlgebra):
+class ComplexHermitianEJA(ComplexMatrixEuclideanJordanAlgebra, KnownRankEJA):
     """
     The rank-n simple EJA consisting of complex Hermitian n-by-n
     matrices over the real numbers, the usual symmetric Jordan product,
@@ -1364,6 +1370,7 @@ class ComplexHermitianEJA(ComplexMatrixEuclideanJordanAlgebra):
         True
 
     """
+
     @classmethod
     def _denormalized_basis(cls, n, field):
         """
@@ -1415,8 +1422,12 @@ class ComplexHermitianEJA(ComplexMatrixEuclideanJordanAlgebra):
 
         # Since we embedded these, we can drop back to the "field" that we
         # started with instead of the complex extension "F".
-        return tuple( s.change_ring(field) for s in S )
+        return ( s.change_ring(field) for s in S )
+
 
+    def __init__(self, n, field=QQ, **kwargs):
+        basis = self._denormalized_basis(n,field)
+        super(ComplexHermitianEJA,self).__init__(field, basis, n, **kwargs)
 
 
 class QuaternionMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra):
@@ -1582,7 +1593,8 @@ class QuaternionMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra):
         return RealMatrixEuclideanJordanAlgebra.natural_inner_product(X,Y)/4
 
 
-class QuaternionHermitianEJA(QuaternionMatrixEuclideanJordanAlgebra):
+class QuaternionHermitianEJA(QuaternionMatrixEuclideanJordanAlgebra,
+                             KnownRankEJA):
     """
     The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
     matrices, the usual symmetric Jordan product, and the
@@ -1693,11 +1705,15 @@ class QuaternionHermitianEJA(QuaternionMatrixEuclideanJordanAlgebra):
                     S.append(Sij_J)
                     Sij_K = cls.real_embed(K*Eij - K*Eij.transpose())
                     S.append(Sij_K)
-        return tuple(S)
+        return S
 
 
+    def __init__(self, n, field=QQ, **kwargs):
+        basis = self._denormalized_basis(n,field)
+        super(QuaternionHermitianEJA,self).__init__(field, basis, n, **kwargs)
+
 
-class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra):
+class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra, KnownRankEJA):
     """
     The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
     with the usual inner product and jordan product ``x*y =