X-Git-Url: http://gitweb.michael.orlitzky.com/?p=sage.d.git;a=blobdiff_plain;f=mjo%2Feja%2Feja_algebra.py;h=2bad32c2f500193e4126b7c5e209c0acb3116ede;hp=8e522527413b89c0e704687583dc9c4636303170;hb=f0cabe7c6e37781e4f92c9ba0e0c7413a5f6b939;hpb=da6160ad4c8156d51290fd056aa7d5b3d609c8e9 diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index 8e52252..2bad32c 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -166,9 +166,21 @@ from sage.modules.free_module import FreeModule, VectorSpace from sage.rings.all import (ZZ, QQ, AA, QQbar, RR, RLF, CLF, PolynomialRing, QuadraticField) -from mjo.eja.eja_element import FiniteDimensionalEJAElement +from mjo.eja.eja_element import (CartesianProductEJAElement, + FiniteDimensionalEJAElement) from mjo.eja.eja_operator import FiniteDimensionalEJAOperator -from mjo.eja.eja_utils import _all2list, _mat2vec +from mjo.eja.eja_utils import _all2list + +def EuclideanJordanAlgebras(field): + r""" + The category of Euclidean Jordan algebras over ``field``, which + must be a subfield of the real numbers. For now this is just a + convenient wrapper around all of the other category axioms that + apply to all EJAs. + """ + category = MagmaticAlgebras(field).FiniteDimensional() + category = category.WithBasis().Unital().Commutative() + return category class FiniteDimensionalEJA(CombinatorialFreeModule): r""" @@ -228,6 +240,26 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): """ Element = FiniteDimensionalEJAElement + @staticmethod + def _check_input_field(field): + if not field.is_subring(RR): + # Note: this does return true for the real algebraic + # field, the rationals, and any quadratic field where + # we've specified a real embedding. + raise ValueError("scalar field is not real") + + @staticmethod + def _check_input_axioms(basis, jordan_product, inner_product): + if not all( jordan_product(bi,bj) == jordan_product(bj,bi) + for bi in basis + for bj in basis ): + raise ValueError("Jordan product is not commutative") + + if not all( inner_product(bi,bj) == inner_product(bj,bi) + for bi in basis + for bj in basis ): + raise ValueError("inner-product is not commutative") + def __init__(self, basis, jordan_product, @@ -236,7 +268,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): matrix_space=None, orthonormalize=True, associative=None, - cartesian_product=False, check_field=True, check_axioms=True, prefix="b"): @@ -244,30 +275,14 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): n = len(basis) if check_field: - if not field.is_subring(RR): - # Note: this does return true for the real algebraic - # field, the rationals, and any quadratic field where - # we've specified a real embedding. - raise ValueError("scalar field is not real") + self._check_input_field(field) if check_axioms: # Check commutativity of the Jordan and inner-products. # This has to be done before we build the multiplication # and inner-product tables/matrices, because we take # advantage of symmetry in the process. - if not all( jordan_product(bi,bj) == jordan_product(bj,bi) - for bi in basis - for bj in basis ): - raise ValueError("Jordan product is not commutative") - - if not all( inner_product(bi,bj) == inner_product(bj,bi) - for bi in basis - for bj in basis ): - raise ValueError("inner-product is not commutative") - - - category = MagmaticAlgebras(field).FiniteDimensional() - category = category.WithBasis().Unital().Commutative() + self._check_input_axioms(basis, jordan_product, inner_product) if n <= 1: # All zero- and one-dimensional algebras are just the real @@ -286,14 +301,11 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): for bj in basis for bk in basis) + category = EuclideanJordanAlgebras(field) + if associative: # Element subalgebras can take advantage of this. category = category.Associative() - if cartesian_product: - # 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. @@ -368,7 +380,8 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # Now we actually compute the multiplication and inner-product # tables/matrices using the possibly-orthonormalized basis. self._inner_product_matrix = matrix.identity(field, n) - self._multiplication_table = [ [0 for j in range(i+1)] + zed = self.zero() + self._multiplication_table = [ [zed for j in range(i+1)] for i in range(n) ] # Note: the Jordan and inner-products are defined in terms @@ -779,13 +792,11 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): if elt.parent().superalgebra() == self: return elt.superalgebra_element() - try: - # Try to convert a vector into a column-matrix... + if hasattr(elt, 'sparse_vector'): + # Convert a vector into a column-matrix. We check for + # "sparse_vector" and not "column" because matrices also + # have a "column" method. elt = elt.column() - except (AttributeError, TypeError): - # and ignore failure, because we weren't really expecting - # a vector as an argument anyway. - pass if elt not in self.matrix_space(): raise ValueError(msg) @@ -1271,7 +1282,9 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # # Of course, matrices aren't vectors in sage, so we have to # appeal to the "long vectors" isometry. - oper_vecs = [ _mat2vec(g.operator().matrix()) for g in self.gens() ] + + V = VectorSpace(self.base_ring(), self.dimension()**2) + oper_vecs = [ V(g.operator().matrix().list()) for g in self.gens() ] # Now we use basic linear algebra to find the coefficients, # of the matrices-as-vectors-linear-combination, which should @@ -1281,7 +1294,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # We used the isometry on the left-hand side already, but we # still need to do it for the right-hand side. Recall that we # wanted something that summed to the identity matrix. - b = _mat2vec( matrix.identity(self.base_ring(), self.dimension()) ) + b = V( matrix.identity(self.base_ring(), self.dimension()).list() ) # Now if there's an identity element in the algebra, this # should work. We solve on the left to avoid having to @@ -1411,7 +1424,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # corresponding to trivial spaces (e.g. it returns only the # eigenspace corresponding to lambda=1 if you take the # decomposition relative to the identity element). - trivial = self.subalgebra(()) + trivial = self.subalgebra((), check_axioms=False) J0 = trivial # eigenvalue zero J5 = VectorSpace(self.base_ring(), 0) # eigenvalue one-half J1 = trivial # eigenvalue one @@ -1733,6 +1746,15 @@ class RationalBasisEJA(FiniteDimensionalEJA): check_field=False, check_axioms=False) + def rational_algebra(self): + # Using None as a flag here (rather than just assigning "self" + # to self._rational_algebra by default) feels a little bit + # more sane to me in a garbage-collected environment. + if self._rational_algebra is None: + return self + else: + return self._rational_algebra + @cached_method def _charpoly_coefficients(self): r""" @@ -1757,18 +1779,15 @@ class RationalBasisEJA(FiniteDimensionalEJA): Algebraic Real Field """ - if self._rational_algebra is None: - # There's no need to construct *another* algebra over the - # rationals if this one is already over the - # rationals. Likewise, if we never orthonormalized our - # basis, we might as well just use the given one. + if self.rational_algebra() is self: + # Bypass the hijinks if they won't benefit us. return super()._charpoly_coefficients() # Do the computation over the rationals. The answer will be # the same, because all we've done is a change of basis. # Then, change back from QQ to our real base ring a = ( a_i.change_ring(self.base_ring()) - for a_i in self._rational_algebra._charpoly_coefficients() ) + for a_i in self.rational_algebra()._charpoly_coefficients() ) # Otherwise, convert the coordinate variables back to the # deorthonormalized ones. @@ -2029,7 +2048,6 @@ class MatrixEJA(FiniteDimensionalEJA): # if the user passes check_axioms=True. if "check_axioms" not in kwargs: kwargs["check_axioms"] = False - super().__init__(self._denormalized_basis(matrix_space), self.jordan_product, self.trace_inner_product, @@ -2124,20 +2142,13 @@ class RealSymmetricEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA): return cls(n, **kwargs) 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. - if "check_axioms" not in kwargs: kwargs["check_axioms"] = False - A = MatrixSpace(field, n) super().__init__(A, **kwargs) from mjo.eja.eja_cache import real_symmetric_eja_coeffs a = real_symmetric_eja_coeffs(self) if a is not None: - if self._rational_algebra is None: - self._charpoly_coefficients.set_cache(a) - else: - self._rational_algebra._charpoly_coefficients.set_cache(a) + self.rational_algebra()._charpoly_coefficients.set_cache(a) @@ -2214,10 +2225,6 @@ class ComplexHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA): """ 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. - if "check_axioms" not in kwargs: kwargs["check_axioms"] = False - from mjo.hurwitz import ComplexMatrixAlgebra A = ComplexMatrixAlgebra(n, scalars=field) super().__init__(A, **kwargs) @@ -2225,10 +2232,7 @@ class ComplexHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA): from mjo.eja.eja_cache import complex_hermitian_eja_coeffs a = complex_hermitian_eja_coeffs(self) if a is not None: - if self._rational_algebra is None: - self._charpoly_coefficients.set_cache(a) - else: - self._rational_algebra._charpoly_coefficients.set_cache(a) + self.rational_algebra()._charpoly_coefficients.set_cache(a) @staticmethod def _max_random_instance_size(max_dimension): @@ -2307,10 +2311,6 @@ class QuaternionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA): """ 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. - if "check_axioms" not in kwargs: kwargs["check_axioms"] = False - from mjo.hurwitz import QuaternionMatrixAlgebra A = QuaternionMatrixAlgebra(n, scalars=field) super().__init__(A, **kwargs) @@ -2318,10 +2318,7 @@ class QuaternionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA): from mjo.eja.eja_cache import quaternion_hermitian_eja_coeffs a = quaternion_hermitian_eja_coeffs(self) if a is not None: - if self._rational_algebra is None: - self._charpoly_coefficients.set_cache(a) - else: - self._rational_algebra._charpoly_coefficients.set_cache(a) + self.rational_algebra()._charpoly_coefficients.set_cache(a) @@ -2477,10 +2474,7 @@ class OctonionHermitianEJA(MatrixEJA, RationalBasisEJA, ConcreteEJA): from mjo.eja.eja_cache import octonion_hermitian_eja_coeffs a = octonion_hermitian_eja_coeffs(self) if a is not None: - if self._rational_algebra is None: - self._charpoly_coefficients.set_cache(a) - else: - self._rational_algebra._charpoly_coefficients.set_cache(a) + self.rational_algebra()._charpoly_coefficients.set_cache(a) class AlbertEJA(OctonionHermitianEJA): @@ -2939,6 +2933,7 @@ class CartesianProductEJA(FiniteDimensionalEJA): sage: from mjo.eja.eja_algebra import (random_eja, ....: CartesianProductEJA, + ....: ComplexHermitianEJA, ....: HadamardEJA, ....: JordanSpinEJA, ....: RealSymmetricEJA) @@ -3050,6 +3045,28 @@ class CartesianProductEJA(FiniteDimensionalEJA): | b2 || 0 | 0 | b2 | +----++----+----+----+ + The "matrix space" of a Cartesian product always consists of + ordered pairs (or triples, or...) whose components are the + matrix spaces of its factors:: + + sage: J1 = HadamardEJA(2) + sage: J2 = ComplexHermitianEJA(2) + sage: J = cartesian_product([J1,J2]) + sage: J.matrix_space() + The Cartesian product of (Full MatrixSpace of 2 by 1 dense + matrices over Algebraic Real Field, Module of 2 by 2 matrices + with entries in Algebraic Field over the scalar ring Algebraic + Real Field) + sage: J.one().to_matrix()[0] + [1] + [1] + sage: J.one().to_matrix()[1] + +---+---+ + | 1 | 0 | + +---+---+ + | 0 | 1 | + +---+---+ + TESTS: All factors must share the same base field:: @@ -3072,11 +3089,8 @@ class CartesianProductEJA(FiniteDimensionalEJA): sage: expected = J.one() # long time sage: actual == expected # long time True - """ - Element = FiniteDimensionalEJAElement - - + Element = CartesianProductEJAElement def __init__(self, factors, **kwargs): m = len(factors) if m == 0: @@ -3088,68 +3102,126 @@ class CartesianProductEJA(FiniteDimensionalEJA): if not all( J.base_ring() == field for J in factors ): raise ValueError("all factors must share the same base field") + # Figure out the category to use. associative = all( f.is_associative() for f in factors ) - - # Compute my matrix space. This category isn't perfect, but - # is good enough for what we need to do. + category = EuclideanJordanAlgebras(field) + if associative: category = category.Associative() + category = category.join([category, category.CartesianProducts()]) + + # Compute my matrix space. We don't simply use the + # ``cartesian_product()`` functor here because it acts + # differently on SageMath MatrixSpaces and our custom + # MatrixAlgebras, which are CombinatorialFreeModules. We + # always want the result to be represented (and indexed) as an + # ordered tuple. This category isn't perfect, but is good + # enough for what we need to do. MS_cat = MagmaticAlgebras(field).FiniteDimensional().WithBasis() MS_cat = MS_cat.Unital().CartesianProducts() MS_factors = tuple( J.matrix_space() for J in factors ) from sage.sets.cartesian_product import CartesianProduct - MS = CartesianProduct(MS_factors, MS_cat) + self._matrix_space = CartesianProduct(MS_factors, MS_cat) - basis = [] - zero = MS.zero() + self._matrix_basis = [] + zero = self._matrix_space.zero() for i in range(m): for b in factors[i].matrix_basis(): z = list(zero) z[i] = b - basis.append(z) + self._matrix_basis.append(z) - basis = tuple( MS(b) for b in basis ) + self._matrix_basis = tuple( self._matrix_space(b) + for b in self._matrix_basis ) + n = len(self._matrix_basis) - # Define jordan/inner products that operate on that matrix_basis. - def jordan_product(x,y): - return MS(tuple( - (factors[i](x[i])*factors[i](y[i])).to_matrix() - for i in range(m) - )) - - def inner_product(x, y): - return sum( - factors[i](x[i]).inner_product(factors[i](y[i])) - for i in range(m) - ) + # We already have what we need for the super-superclass constructor. + CombinatorialFreeModule.__init__(self, + field, + range(n), + prefix="b", + category=category, + bracket=False) - # There's no need to check the field since it already came - # from an EJA. Likewise the axioms are guaranteed to be - # satisfied, unless the guy writing this class sucks. - # - # If you want the basis to be orthonormalized, orthonormalize - # the factors. - FiniteDimensionalEJA.__init__(self, - basis, - jordan_product, - inner_product, - field=field, - matrix_space=MS, - orthonormalize=False, - associative=associative, - cartesian_product=True, - check_field=False, - check_axioms=False) + # Now create the vector space for the algebra, which will have + # its own set of non-ambient coordinates (in terms of the + # supplied basis). + degree = sum( f._matrix_span.ambient_vector_space().degree() + for f in factors ) + V = VectorSpace(field, degree) + vector_basis = tuple( V(_all2list(b)) for b in self._matrix_basis ) + + # Save the span of our matrix basis (when written out as long + # vectors) because otherwise we'll have to reconstruct it + # every time we want to coerce a matrix into the algebra. + self._matrix_span = V.span_of_basis( vector_basis, check=False) # Since we don't (re)orthonormalize the basis, the FDEJA # constructor is going to set self._deortho_matrix to the # identity matrix. Here we set it to the correct value using # the deortho matrices from our factors. - self._deortho_matrix = matrix.block_diagonal( [J._deortho_matrix - for J in factors] ) + self._deortho_matrix = matrix.block_diagonal( + [J._deortho_matrix for J in factors] + ) + + self._inner_product_matrix = matrix.block_diagonal( + [J._inner_product_matrix for J in factors] + ) + self._inner_product_matrix._cache = {'hermitian': True} + self._inner_product_matrix.set_immutable() + + # Building the multiplication table is a bit more tricky + # because we have to embed the entries of the factors' + # multiplication tables into the product EJA. + zed = self.zero() + self._multiplication_table = [ [zed for j in range(i+1)] + for i in range(n) ] + + # Keep track of an offset that tallies the dimensions of all + # previous factors. If the second factor is dim=2 and if the + # first one is dim=3, then we want to skip the first 3x3 block + # when copying the multiplication table for the second factor. + offset = 0 + for f in range(m): + phi_f = self.cartesian_embedding(f) + factor_dim = factors[f].dimension() + for i in range(factor_dim): + for j in range(i+1): + f_ij = factors[f]._multiplication_table[i][j] + e = phi_f(f_ij) + self._multiplication_table[offset+i][offset+j] = e + offset += factor_dim self.rank.set_cache(sum(J.rank() for J in factors)) ones = tuple(J.one().to_matrix() for J in factors) self.one.set_cache(self(ones)) + def _sets_keys(self): + r""" + + SETUP:: + + sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA, + ....: RealSymmetricEJA) + + TESTS: + + The superclass uses ``_sets_keys()`` to implement its + ``cartesian_factors()`` method:: + + sage: J1 = RealSymmetricEJA(2, + ....: field=QQ, + ....: orthonormalize=False, + ....: prefix="a") + sage: J2 = ComplexHermitianEJA(2,field=QQ,orthonormalize=False) + sage: J = cartesian_product([J1,J2]) + sage: x = sum(i*J.gens()[i] for i in range(len(J.gens()))) + sage: x.cartesian_factors() + (a1 + 2*a2, 3*b0 + 4*b1 + 5*b2 + 6*b3) + + """ + # Copy/pasted from CombinatorialFreeModule_CartesianProduct, + # but returning a tuple instead of a list. + return tuple(range(len(self.cartesian_factors()))) + def cartesian_factors(self): # Copy/pasted from CombinatorialFreeModule_CartesianProduct. return self._sets @@ -3166,65 +3238,6 @@ class CartesianProductEJA(FiniteDimensionalEJA): return cartesian_product.symbol.join("%s" % factor for factor in self._sets) - def matrix_space(self): - r""" - Return the space that our matrix basis lives in as a Cartesian - product. - - We don't simply use the ``cartesian_product()`` functor here - because it acts differently on SageMath MatrixSpaces and our - custom MatrixAlgebras, which are CombinatorialFreeModules. We - always want the result to be represented (and indexed) as - an ordered tuple. - - SETUP:: - - sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA, - ....: HadamardEJA, - ....: OctonionHermitianEJA, - ....: RealSymmetricEJA) - - EXAMPLES:: - - sage: J1 = HadamardEJA(1) - sage: J2 = RealSymmetricEJA(2) - sage: J = cartesian_product([J1,J2]) - sage: J.matrix_space() - The Cartesian product of (Full MatrixSpace of 1 by 1 dense - matrices over Algebraic Real Field, Full MatrixSpace of 2 - by 2 dense matrices over Algebraic Real Field) - - :: - - sage: J1 = ComplexHermitianEJA(1) - sage: J2 = ComplexHermitianEJA(1) - sage: J = cartesian_product([J1,J2]) - sage: J.one().to_matrix()[0] - +---+ - | 1 | - +---+ - sage: J.one().to_matrix()[1] - +---+ - | 1 | - +---+ - - :: - - sage: J1 = OctonionHermitianEJA(1) - sage: J2 = OctonionHermitianEJA(1) - sage: J = cartesian_product([J1,J2]) - sage: J.one().to_matrix()[0] - +----+ - | e0 | - +----+ - sage: J.one().to_matrix()[1] - +----+ - | e0 | - +----+ - - """ - return super().matrix_space() - @cached_method def cartesian_projection(self, i): @@ -3426,9 +3439,9 @@ class RationalBasisCartesianProductEJA(CartesianProductEJA, SETUP:: - sage: from mjo.eja.eja_algebra import (HadamardEJA, + sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA, + ....: HadamardEJA, ....: JordanSpinEJA, - ....: OctonionHermitianEJA, ....: RealSymmetricEJA) EXAMPLES: @@ -3449,28 +3462,38 @@ class RationalBasisCartesianProductEJA(CartesianProductEJA, The ``cartesian_product()`` function only uses the first factor to decide where the result will live; thus we have to be careful to - check that all factors do indeed have a `_rational_algebra` member - before we try to access it:: - - sage: J1 = OctonionHermitianEJA(1) # no rational basis - sage: J2 = HadamardEJA(2) - sage: cartesian_product([J1,J2]) - Euclidean Jordan algebra of dimension 1 over Algebraic Real Field - (+) Euclidean Jordan algebra of dimension 2 over Algebraic Real Field - sage: cartesian_product([J2,J1]) - Euclidean Jordan algebra of dimension 2 over Algebraic Real Field - (+) Euclidean Jordan algebra of dimension 1 over Algebraic Real Field + check that all factors do indeed have a ``rational_algebra()`` method + before we construct an algebra that claims to have a rational basis:: + + sage: J1 = HadamardEJA(2) + sage: jp = lambda X,Y: X*Y + sage: ip = lambda X,Y: X[0,0]*Y[0,0] + sage: b1 = matrix(QQ, [[1]]) + sage: J2 = FiniteDimensionalEJA((b1,), jp, ip) + sage: cartesian_product([J2,J1]) # factor one not RationalBasisEJA + Euclidean Jordan algebra of dimension 1 over Algebraic Real + Field (+) Euclidean Jordan algebra of dimension 2 over Algebraic + Real Field + sage: cartesian_product([J1,J2]) # factor one is RationalBasisEJA + Traceback (most recent call last): + ... + ValueError: factor not a RationalBasisEJA """ def __init__(self, algebras, **kwargs): + if not all( hasattr(r, "rational_algebra") for r in algebras ): + raise ValueError("factor not a RationalBasisEJA") + CartesianProductEJA.__init__(self, algebras, **kwargs) - self._rational_algebra = None - if self.vector_space().base_field() is not QQ: - if all( hasattr(r, "_rational_algebra") for r in algebras ): - self._rational_algebra = cartesian_product([ - r._rational_algebra for r in algebras - ]) + @cached_method + def rational_algebra(self): + if self.base_ring() is QQ: + return self + + return cartesian_product([ + r.rational_algebra() for r in self.cartesian_factors() + ]) RationalBasisEJA.CartesianProduct = RationalBasisCartesianProductEJA