]> gitweb.michael.orlitzky.com - sage.d.git/commitdiff
eja: don't real-embed quaternion matrices.
authorMichael Orlitzky <michael@orlitzky.com>
Mon, 8 Mar 2021 21:00:31 +0000 (16:00 -0500)
committerMichael Orlitzky <michael@orlitzky.com>
Mon, 8 Mar 2021 21:00:31 +0000 (16:00 -0500)
mjo/all.py
mjo/eja/eja_algebra.py
mjo/eja/eja_utils.py

index f80276732c33e8cf3efd2384d3b950d52cfb758a..78550cccfb67854b8f1f7397b617ec9c5a8b5b80 100644 (file)
@@ -8,7 +8,7 @@ from mjo.cone.all import *
 from mjo.eja.all import *
 from mjo.interpolation import *
 from mjo.misc import *
-from mjo.octonions import Octonions
+from mjo.hurwitz import Octonions
 from mjo.orthogonal_polynomials import *
 from mjo.polynomial import *
 from mjo.symbol_sequence import *
index 2ecf8cd274bbee0be8c31e25c8a35a6e9561123a..827d40214cac92cb6e1ee46eb17f4285a44398f6 100644 (file)
@@ -1019,7 +1019,9 @@ class FiniteDimensionalEJA(CombinatorialFreeModule):
             Full MatrixSpace of 4 by 4 dense matrices over Rational Field
             sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
             sage: J.matrix_space()
-            Full MatrixSpace of 4 by 4 dense matrices over Rational Field
+            Module of 1 by 1 matrices with entries in Quaternion
+            Algebra (-1, -1) with base ring Rational Field over
+            the scalar ring Rational Field
 
         """
         if self.is_trivial():
@@ -1796,8 +1798,7 @@ class RealEmbeddedMatrixEJA(MatrixEJA):
 
         SETUP::
 
-            sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
-            ....:                                  QuaternionHermitianEJA)
+            sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
 
         EXAMPLES::
 
@@ -1813,20 +1814,6 @@ class RealEmbeddedMatrixEJA(MatrixEJA):
             sage: actual == expected
             True
 
-        ::
-
-            sage: set_random_seed()
-            sage: J = QuaternionHermitianEJA.random_instance()
-            sage: x,y = J.random_elements(2)
-            sage: Xe = x.to_matrix()
-            sage: Ye = y.to_matrix()
-            sage: X = J.real_unembed(Xe)
-            sage: Y = J.real_unembed(Ye)
-            sage: expected = (X*Y).trace().coefficient_tuple()[0]
-            sage: actual = J.trace_inner_product(Xe,Ye)
-            sage: actual == expected
-            True
-
         """
         # This does in fact compute the real part of the trace.
         # If we compute the trace of e.g. a complex matrix M,
@@ -2272,157 +2259,8 @@ class ComplexHermitianEJA(RationalBasisEJA, ConcreteEJA, ComplexMatrixEJA):
         n = ZZ.random_element(cls._max_random_instance_size() + 1)
         return cls(n, **kwargs)
 
-class QuaternionMatrixEJA(RealEmbeddedMatrixEJA):
-
-    # A manual dictionary-cache for the quaternion_extension() method,
-    # since apparently @classmethods can't also be @cached_methods.
-    _quaternion_extension = {}
-
-    @classmethod
-    def quaternion_extension(cls,field):
-        r"""
-        The quaternion field that we embed/unembed, as an extension
-        of the given ``field``.
-        """
-        if field in cls._quaternion_extension:
-            return cls._quaternion_extension[field]
-
-        Q = QuaternionAlgebra(field,-1,-1)
-
-        cls._quaternion_extension[field] = Q
-        return Q
-
-    @staticmethod
-    def dimension_over_reals():
-        return 4
-
-    @classmethod
-    def real_embed(cls,M):
-        """
-        Embed the n-by-n quaternion matrix ``M`` into the space of real
-        matrices of size 4n-by-4n by first sending each quaternion entry `z
-        = a + bi + cj + dk` to the block-complex matrix ``[[a + bi,
-        c+di],[-c + di, a-bi]]`, and then embedding those into a real
-        matrix.
-
-        SETUP::
-
-            sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
-
-        EXAMPLES::
-
-            sage: Q = QuaternionAlgebra(QQ,-1,-1)
-            sage: i,j,k = Q.gens()
-            sage: x = 1 + 2*i + 3*j + 4*k
-            sage: M = matrix(Q, 1, [[x]])
-            sage: QuaternionMatrixEJA.real_embed(M)
-            [ 1  2  3  4]
-            [-2  1 -4  3]
-            [-3  4  1 -2]
-            [-4 -3  2  1]
-
-        Embedding is a homomorphism (isomorphism, in fact)::
-
-            sage: set_random_seed()
-            sage: n = ZZ.random_element(2)
-            sage: Q = QuaternionAlgebra(QQ,-1,-1)
-            sage: X = random_matrix(Q, n)
-            sage: Y = random_matrix(Q, n)
-            sage: Xe = QuaternionMatrixEJA.real_embed(X)
-            sage: Ye = QuaternionMatrixEJA.real_embed(Y)
-            sage: XYe = QuaternionMatrixEJA.real_embed(X*Y)
-            sage: Xe*Ye == XYe
-            True
-
-        """
-        super().real_embed(M)
-        quaternions = M.base_ring()
-        n = M.nrows()
-
-        F = QuadraticField(-1, 'I')
-        i = F.gen()
-
-        blocks = []
-        for z in M.list():
-            t = z.coefficient_tuple()
-            a = t[0]
-            b = t[1]
-            c = t[2]
-            d = t[3]
-            cplxM = matrix(F, 2, [[ a + b*i, c + d*i],
-                                 [-c + d*i, a - b*i]])
-            realM = ComplexMatrixEJA.real_embed(cplxM)
-            blocks.append(realM)
-
-        # We should have real entries by now, so use the realest field
-        # we've got for the return value.
-        return matrix.block(quaternions.base_ring(), n, blocks)
-
-
-
-    @classmethod
-    def real_unembed(cls,M):
-        """
-        The inverse of _embed_quaternion_matrix().
-
-        SETUP::
-
-            sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
-
-        EXAMPLES::
-
-            sage: M = matrix(QQ, [[ 1,  2,  3,  4],
-            ....:                 [-2,  1, -4,  3],
-            ....:                 [-3,  4,  1, -2],
-            ....:                 [-4, -3,  2,  1]])
-            sage: QuaternionMatrixEJA.real_unembed(M)
-            [1 + 2*i + 3*j + 4*k]
-
-        TESTS:
-
-        Unembedding is the inverse of embedding::
-
-            sage: set_random_seed()
-            sage: Q = QuaternionAlgebra(QQ, -1, -1)
-            sage: M = random_matrix(Q, 3)
-            sage: Me = QuaternionMatrixEJA.real_embed(M)
-            sage: QuaternionMatrixEJA.real_unembed(Me) == M
-            True
-
-        """
-        super().real_unembed(M)
-        n = ZZ(M.nrows())
-        d = cls.dimension_over_reals()
-
-        # Use the base ring of the matrix to ensure that its entries can be
-        # multiplied by elements of the quaternion algebra.
-        Q = cls.quaternion_extension(M.base_ring())
-        i,j,k = Q.gens()
-
-        # Go top-left to bottom-right (reading order), converting every
-        # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
-        # quaternion block.
-        elements = []
-        for l in range(n/d):
-            for m in range(n/d):
-                submat = ComplexMatrixEJA.real_unembed(
-                    M[d*l:d*l+d,d*m:d*m+d] )
-                if submat[0,0] != submat[1,1].conjugate():
-                    raise ValueError('bad on-diagonal submatrix')
-                if submat[0,1] != -submat[1,0].conjugate():
-                    raise ValueError('bad off-diagonal submatrix')
-                z  = submat[0,0].real()
-                z += submat[0,0].imag()*i
-                z += submat[0,1].real()*j
-                z += submat[0,1].imag()*k
-                elements.append(z)
-
-        return matrix(Q, n/d, elements)
-
 
-class QuaternionHermitianEJA(RationalBasisEJA,
-                             ConcreteEJA,
-                             QuaternionMatrixEJA):
+class QuaternionHermitianEJA(RationalBasisEJA, ConcreteEJA, MatrixEJA):
     r"""
     The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
     matrices, the usual symmetric Jordan product, and the
@@ -2498,59 +2336,56 @@ class QuaternionHermitianEJA(RationalBasisEJA,
 
             sage: set_random_seed()
             sage: n = ZZ.random_element(1,5)
-            sage: B = QuaternionHermitianEJA._denormalized_basis(n,ZZ)
-            sage: all( M.is_symmetric() for M in B )
+            sage: B = QuaternionHermitianEJA._denormalized_basis(n,QQ)
+            sage: all( M.is_hermitian() for M in B )
             True
 
         """
-        Q = QuaternionAlgebra(QQ,-1,-1)
-        I,J,K = Q.gens()
+        from mjo.hurwitz import QuaternionMatrixAlgebra
+        A = QuaternionMatrixAlgebra(n, scalars=field)
+        es = A.entry_algebra_gens()
 
-        # This is like the symmetric case, but we need to be careful:
-        #
-        #   * We want conjugate-symmetry, not just symmetry.
-        #   * The diagonal will (as a result) be real.
-        #
-        S = []
-        Eij = matrix.zero(Q,n)
+        basis = []
         for i in range(n):
             for j in range(i+1):
-                # "build" E_ij
-                Eij[i,j] = 1
                 if i == j:
-                    Sij = cls.real_embed(Eij)
-                    S.append(Sij)
+                    E_ii = A.monomial( (i,j,es[0]) )
+                    basis.append(E_ii)
                 else:
-                    # The second, third, and fourth ones have a minus
-                    # because they're conjugated.
-                    # Eij = Eij + Eij.transpose()
-                    Eij[j,i] = 1
-                    Sij_real = cls.real_embed(Eij)
-                    S.append(Sij_real)
-                    # Eij = I*(Eij - Eij.transpose())
-                    Eij[i,j] = I
-                    Eij[j,i] = -I
-                    Sij_I = cls.real_embed(Eij)
-                    S.append(Sij_I)
-                    # Eij = J*(Eij - Eij.transpose())
-                    Eij[i,j] = J
-                    Eij[j,i] = -J
-                    Sij_J = cls.real_embed(Eij)
-                    S.append(Sij_J)
-                    # Eij = K*(Eij - Eij.transpose())
-                    Eij[i,j] = K
-                    Eij[j,i] = -K
-                    Sij_K = cls.real_embed(Eij)
-                    S.append(Sij_K)
-                    Eij[j,i] = 0
-                # "erase" E_ij
-                Eij[i,j] = 0
+                    for e in es:
+                        E_ij  = A.monomial( (i,j,e)             )
+                        ec = e.conjugate()
+                        # If the conjugate has a negative sign in front
+                        # of it, (j,i,ec) won't be a monomial!
+                        if (j,i,ec) in A.indices():
+                            E_ij += A.monomial( (j,i,ec) )
+                        else:
+                            E_ij -= A.monomial( (j,i,-ec) )
+                        basis.append(E_ij)
 
-        # Since we embedded the entries, we can drop back to the
-        # desired real "field" instead of the quaternion algebra "Q".
-        return tuple( s.change_ring(field) for s in S )
+        return tuple( basis )
 
 
+    @staticmethod
+    def trace_inner_product(X,Y):
+        r"""
+        Overload the superclass method because the quaternions are weird
+        and we need to use ``coefficient_tuple()`` to get the realpart.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
+
+        TESTS::
+
+            sage: J = QuaternionHermitianEJA(2,field=QQ,orthonormalize=False)
+            sage: I = J.one().to_matrix()
+            sage: J.trace_inner_product(I, -I)
+            -2
+
+        """
+        return (X*Y).trace().coefficient_tuple()[0]
+
     def __init__(self, n, field=AA, **kwargs):
         # We know this is a valid EJA, but will double-check
         # if the user passes check_axioms=True.
@@ -2571,7 +2406,7 @@ class QuaternionHermitianEJA(RationalBasisEJA,
         # because the MatrixEJA is not presently a subclass of the
         # FDEJA class that defines rank() and one().
         self.rank.set_cache(n)
-        idV = matrix.identity(ZZ, self.dimension_over_reals()*n)
+        idV = self.matrix_space().one()
         self.one.set_cache(self(idV))
 
 
@@ -2733,25 +2568,25 @@ class OctonionHermitianEJA(RationalBasisEJA, ConcreteEJA, MatrixEJA):
 
         """
         from mjo.hurwitz import OctonionMatrixAlgebra
-        MS = OctonionMatrixAlgebra(n, scalars=field)
-        es = MS.entry_algebra().gens()
+        A = OctonionMatrixAlgebra(n, scalars=field)
+        es = A.entry_algebra_gens()
 
         basis = []
         for i in range(n):
             for j in range(i+1):
                 if i == j:
-                    E_ii = MS.monomial( (i,j,es[0]) )
+                    E_ii = A.monomial( (i,j,es[0]) )
                     basis.append(E_ii)
                 else:
                     for e in es:
-                        E_ij  = MS.monomial( (i,j,e)             )
+                        E_ij  = A.monomial( (i,j,e)             )
                         ec = e.conjugate()
                         # If the conjugate has a negative sign in front
                         # of it, (j,i,ec) won't be a monomial!
-                        if (j,i,ec) in MS.indices():
-                            E_ij += MS.monomial( (j,i,ec) )
+                        if (j,i,ec) in A.indices():
+                            E_ij += A.monomial( (j,i,ec) )
                         else:
-                            E_ij -= MS.monomial( (j,i,-ec) )
+                            E_ij -= A.monomial( (j,i,-ec) )
                         basis.append(E_ij)
 
         return tuple( basis )
@@ -2774,7 +2609,7 @@ class OctonionHermitianEJA(RationalBasisEJA, ConcreteEJA, MatrixEJA):
             -2
 
         """
-        return (X*Y).trace().real().coefficient(0)
+        return (X*Y).trace().coefficient(0)
 
 
 class AlbertEJA(OctonionHermitianEJA):
index 0b2d2a315989949c2431641c8f82dea9b576f9b8..6f8cab6d8019dcbba0be1e81c3872a7ba738f807 100644 (file)
@@ -54,7 +54,9 @@ def _all2list(x):
     SETUP::
 
         sage: from mjo.eja.eja_utils import _all2list
-        sage: from mjo.octonions import Octonions, OctonionMatrixAlgebra
+        sage: from mjo.hurwitz import (QuaternionMatrixAlgebra,
+        ....:                          Octonions,
+        ....:                          OctonionMatrixAlgebra)
 
     EXAMPLES::
 
@@ -86,6 +88,13 @@ def _all2list(x):
         sage: _all2list(OctonionMatrixAlgebra(1).one())
         [1, 0, 0, 0, 0, 0, 0, 0]
 
+    ::
+
+        sage: _all2list(QuaternionAlgebra(QQ, -1, -1).one())
+        [1, 0, 0, 0]
+        sage: _all2list(QuaternionMatrixAlgebra(1).one())
+        [1, 0, 0, 0]
+
     ::
 
         sage: V1 = VectorSpace(QQ,2)
@@ -97,6 +106,13 @@ def _all2list(x):
         [3, 4, 1, 0, 0, 0, 0, 0, 0, 0]
 
     """
+    if hasattr(x, 'list') and hasattr(x, 'to_vector'):
+        # This avoids calling to_vector() on a matrix algebra with
+        # e.g. quaternions where the returned vector is of the wrong
+        # length (three instead of four) because the quaternions don't
+        # know how many generators they have.
+        return _all2list(x.list())
+
     if hasattr(x, 'to_vector'):
         # This works on matrices of e.g. octonions directly, without
         # first needing to convert them to a list of octonions and