From fa3651ef28df68f4c37b16339074fbf86b8feff4 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Thu, 27 Jan 2022 21:01:22 -0500 Subject: [PATCH] eja: use the correct definition for the new algebra. --- mjo/eja/all.py | 2 +- mjo/eja/eja_algebra.py | 181 +++++++++++++++++++++++++++++++---------- 2 files changed, 141 insertions(+), 42 deletions(-) diff --git a/mjo/eja/all.py b/mjo/eja/all.py index 9fe9cb7..5161195 100644 --- a/mjo/eja/all.py +++ b/mjo/eja/all.py @@ -5,7 +5,7 @@ All imports from mjo.eja modules. from mjo.eja.eja_algebra import (AlbertEJA, BilinearFormEJA, ComplexHermitianEJA, - ComplexSkewHermitianEJA, + ComplexSkewSymmetricEJA, HadamardEJA, JordanSpinEJA, OctonionHermitianEJA, diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index fb01646..3671e36 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -3544,24 +3544,124 @@ def random_eja(max_dimension=None, *args, **kwargs): return cartesian_product([J1,J2]) -class ComplexSkewHermitianEJA(RationalBasisEJA): +class ComplexSkewSymmetricEJA(RationalBasisEJA, ConcreteEJA): r""" - The EJA described in Faraut and Koranyi's Exercise III.1.b. + The skew-symmetric EJA of order `2n` described in Faraut and + Koranyi's Exercise III.1.b. It has dimension `2n^2 - n`. + + It is (not obviously) isomorphic to the QuaternionHermitianEJA of + order `n`, as can be inferred by comparing rank/dimension or + explicitly from their "characteristic polynomial of" functions, + which just so happen to align nicely. + + SETUP:: + + sage: from mjo.eja.eja_algebra import (ComplexSkewSymmetricEJA, + ....: QuaternionHermitianEJA) + sage: from mjo.eja.eja_operator import FiniteDimensionalEJAOperator + + EXAMPLES: + + This EJA is isomorphic to the quaternions:: + + sage: J = ComplexSkewSymmetricEJA(2, field=QQ, orthonormalize=False) + sage: K = QuaternionHermitianEJA(2, field=QQ, orthonormalize=False) + sage: jordan_isom_matrix = matrix.diagonal(QQ,[-1,1,1,1,1,-1]) + sage: phi = FiniteDimensionalEJAOperator(J,K,jordan_isom_matrix) + sage: all( phi(x*y) == phi(x)*phi(y) + ....: for x in J.gens() + ....: for y in J.gens() ) + True + sage: x,y = J.random_elements(2) + sage: phi(x*y) == phi(x)*phi(y) + True + + TESTS: + + Random elements should satisfy the same conditions that the basis + elements do:: + + sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ, + ....: orthonormalize=False) + sage: x,y = K.random_elements(2) + sage: z = x*y + sage: x = x.to_matrix() + sage: y = y.to_matrix() + sage: z = z.to_matrix() + sage: all( e.is_skew_symmetric() for e in (x,y,z) ) + True + sage: J = -K.one().to_matrix() + sage: all( e*J == J*e.conjugate() for e in (x,y,z) ) + True + + The power law in Faraut & Koranyi's II.7.a is satisfied. + We're in a subalgebra of theirs, but powers are still + defined the same:: + + sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ, + ....: orthonormalize=False) + sage: x = K.random_element() + sage: k = ZZ.random_element(5) + sage: actual = x^k + sage: J = -K.one().to_matrix() + sage: expected = K(-J*(J*x.to_matrix())^k) + sage: actual == expected + True + """ + @staticmethod + def _max_random_instance_size(max_dimension): + # Obtained by solving d = 2n^2 - n, which comes from noticing + # that, in 2x2 block form, any element of this algebra has a + # free skew-symmetric top-left block, a Hermitian top-right + # block, and two bottom blocks that are determined by the top. + # The ZZ-int-ZZ thing is just "floor." + return ZZ(int(ZZ(8*max_dimension + 1).sqrt()/4 + 1/4)) + + @classmethod + def random_instance(cls, max_dimension=None, *args, **kwargs): + """ + Return a random instance of this type of algebra. + """ + class_max_d = cls._max_random_instance_dimension() + if (max_dimension is None or max_dimension > class_max_d): + max_dimension = class_max_d + max_size = cls._max_random_instance_size(max_dimension) + n = ZZ.random_element(max_size + 1) + return cls(n, **kwargs) + @staticmethod def _denormalized_basis(A): """ SETUP:: sage: from mjo.hurwitz import ComplexMatrixAlgebra - sage: from mjo.eja.eja_algebra import ComplexSkewHermitianEJA + sage: from mjo.eja.eja_algebra import ComplexSkewSymmetricEJA - TESTS:: + TESTS: + + The basis elements are all skew-Hermitian:: + + sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension() + sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max) + sage: n = ZZ.random_element(n_max + 1) + sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ) + sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A) + sage: all( M.is_skew_symmetric() for M in B) + True - sage: n = ZZ.random_element(1,2) + The basis elements ``b`` all satisfy ``b*J == J*b.conjugate()``, + as in the definition of the algebra:: + + sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension() + sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max) + sage: n = ZZ.random_element(n_max + 1) sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ) - sage: B = ComplexSkewHermitianEJA._denormalized_basis(A) - sage: all( M.is_skew_hermitian() for M in B) + sage: I_n = matrix.identity(ZZ, n) + sage: J = matrix.block(ZZ, 2, 2, (0, I_n, -I_n, 0), subdivide=False) + sage: J = A.from_list(J.rows()) + sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A) + sage: all( b*J == J*b.conjugate() for b in B ) True """ @@ -3573,53 +3673,44 @@ class ComplexSkewHermitianEJA(RationalBasisEJA): # The size of the blocks. We're going to treat these thing as # 2x2 block matrices, # - # [ x1 x2 ] - # [ -x2^* x1-conj ] + # [ x1 x2 ] + # [ -x2-conj x1-conj ] # - # where x1 is skew-Hermitian and x2 is symmetric. + # where x1 is skew-symmetric and x2 is Hermitian. # m = A.nrows()/2 # We only loop through the top half of the matrix, because the # bottom can be constructed from the top. for i in range(m): - - # First do the top-left block, which is skew-Hermitian. + # First do the top-left block, which is skew-symmetric. # We can compute the bottom-right block in the process. for j in range(i+1): - if i == j: - # Top-left block's entry. - E_ii = gen(A, (i,j,es[1])) - - # Bottom-right block's entry. - E_ii += gen(A, (i+m,j+m,es[1])).conjugate() - basis.append(E_ii) - else: + if i != j: + # Skew-symmetry implies zeros for (i == j). for e in es: # Top-left block's entry. E_ij = gen(A, (i,j,e)) - E_ij -= E_ij.conjugate_transpose() + E_ij -= gen(A, (j,i,e)) # Bottom-right block's entry. F_ij = gen(A, (i+m,j+m,e)).conjugate() - F_ij -= F_ij.conjugate_transpose() + F_ij -= gen(A, (j+m,i+m,e)).conjugate() basis.append(E_ij + F_ij) - # Now do the top-right block, which is symmetric, and compute + # Now do the top-right block, which is Hermitian, and compute # the bottom-left block along the way. for j in range(m,i+m+1): if (i+m) == j: - # A symmetric (not Hermitian!) complex matrix can - # have both real and complex entries on its - # diagonal. - for e in es: - # Top-right block's entry. - E_ii = gen(A, (i,j,e)) + # Hermitian matrices have real diagonal entries. + # Top-right block's entry. + E_ii = gen(A, (i,j,es[0])) - # Bottom-left block's entry. - E_ii -= gen(A, (i-m,j-m,e)).conjugate() - basis.append(E_ii) + # Bottom-left block's entry. Don't conjugate + # 'cause it's real. + E_ii -= gen(A, (i+m,j-m,es[0])) + basis.append(E_ii) else: for e in es: # Top-right block's entry. BEWARE! We're not @@ -3629,38 +3720,46 @@ class ComplexSkewHermitianEJA(RationalBasisEJA): E_ij = gen(A, (i,j,e)) # Shift it back to non-offset coords, transpose, - # and put it back: + # conjugate, and put it back: # # (i,j) -> (i,j-m) -> (j-m, i) -> (j-m, i+m) - E_ij += gen(A, (j-m,i+m,e)) + E_ij += gen(A, (j-m,i+m,e)).conjugate() # Bottom-left's block's below-diagonal entry. # Just shift the top-right coords down m and # left m. F_ij = -gen(A, (i+m,j-m,e)).conjugate() - F_ij += -gen(A, (j,i,e)).conjugate() + F_ij += -gen(A, (j,i,e)) # double-conjugate cancels basis.append(E_ij + F_ij) return tuple( basis ) + @staticmethod + @cached_method + def _J_matrix(matrix_space): + n = matrix_space.nrows() // 2 + F = matrix_space.base_ring() + I_n = matrix.identity(F, n) + J = matrix.block(F, 2, 2, (0, I_n, -I_n, 0), subdivide=False) + return matrix_space.from_list(J.rows()) + + def J_matrix(self): + return ComplexSkewSymmetricEJA._J_matrix(self.matrix_space()) def __init__(self, n, field=AA, **kwargs): # New code; always check the axioms. - if "check_axioms" not in kwargs: kwargs["check_axioms"] = False + #if "check_axioms" not in kwargs: kwargs["check_axioms"] = False from mjo.hurwitz import ComplexMatrixAlgebra A = ComplexMatrixAlgebra(2*n, scalars=field) - - I_n = matrix.identity(ZZ, n) - J = matrix.block(ZZ, 2, 2, (0, I_n, -I_n, 0), subdivide=False) - J = A.from_list(J.rows()) + J = ComplexSkewSymmetricEJA._J_matrix(A) def jordan_product(X,Y): return (X*J*Y + Y*J*X)/2 def inner_product(X,Y): - return (X*Y.conjugate_transpose()).trace().real() + return (X.conjugate_transpose()*Y).trace().real() super().__init__(self._denormalized_basis(A), jordan_product, -- 2.43.2