X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_algebra.py;h=106a0cddec06355a952e71d75699677b25dd7da9;hb=fc29add6cf1d9ff4e8a240b0f8f2ca6672d4ea57;hp=799490044dbb4d0834577f435afd2dbc89429baf;hpb=78bee5c30c1cd2828d9834fc7d652db21331d4fe;p=sage.d.git diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index 7994900..106a0cd 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -32,22 +32,25 @@ for these simple algebras: * :class:`RealSymmetricEJA` * :class:`ComplexHermitianEJA` * :class:`QuaternionHermitianEJA` + * :class:`OctonionHermitianEJA` -Missing from this list is the algebra of three-by-three octononion -Hermitian matrices, as there is (as of yet) no implementation of the -octonions in SageMath. In addition to these, we provide two other -example constructions, +In addition to these, we provide two other example constructions, + * :class:`JordanSpinEJA` * :class:`HadamardEJA` + * :class:`AlbertEJA` * :class:`TrivialEJA` The Jordan spin algebra is a bilinear form algebra where the bilinear form is the identity. The Hadamard EJA is simply a Cartesian product -of one-dimensional spin algebras. And last but not least, the trivial -EJA is exactly what you think. Cartesian products of these are also -supported using the usual ``cartesian_product()`` function; as a -result, we support (up to isomorphism) all Euclidean Jordan algebras -that don't involve octonions. +of one-dimensional spin algebras. The Albert EJA is simply a special +case of the :class:`OctonionHermitianEJA` where the matrices are +three-by-three and the resulting space has dimension 27. And +last/least, the trivial EJA is exactly what you think it is; it could +also be obtained by constructing a dimension-zero instance of any of +the other algebras. Cartesian products of these are also supported +using the usual ``cartesian_product()`` function; as a result, we +support (up to isomorphism) all Euclidean Jordan algebras. SETUP:: @@ -59,13 +62,10 @@ EXAMPLES:: Euclidean Jordan algebra of dimension... """ -from itertools import repeat - from sage.algebras.quatalg.quaternion_algebra import QuaternionAlgebra from sage.categories.magmatic_algebras import MagmaticAlgebras from sage.categories.sets_cat import cartesian_product -from sage.combinat.free_module import (CombinatorialFreeModule, - CombinatorialFreeModule_CartesianProduct) +from sage.combinat.free_module import CombinatorialFreeModule from sage.matrix.constructor import matrix from sage.matrix.matrix_space import MatrixSpace from sage.misc.cachefunc import cached_method @@ -142,7 +142,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): cartesian_product=False, check_field=True, check_axioms=True, - prefix='e'): + prefix="b"): n = len(basis) @@ -153,11 +153,6 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # we've specified a real embedding. raise ValueError("scalar field is not real") - from mjo.eja.eja_utils import _change_ring - # If the basis given to us wasn't over the field that it's - # supposed to be over, fix that. Or, you know, crash. - basis = tuple( _change_ring(b, field) for b in basis ) - if check_axioms: # Check commutativity of the Jordan and inner-products. # This has to be done before we build the multiplication @@ -177,6 +172,11 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): category = MagmaticAlgebras(field).FiniteDimensional() category = category.WithBasis().Unital().Commutative() + if n <= 1: + # All zero- and one-dimensional algebras are just the real + # numbers with (some positive multiples of) the usual + # multiplication as its Jordan and inner-product. + associative = True 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() @@ -211,7 +211,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # ambient vector space V that our (vectorized) basis lives in, # as well as a subspace W of V spanned by those (vectorized) # basis elements. The W-coordinates are the coefficients that - # we see in things like x = 1*e1 + 2*e2. + # we see in things like x = 1*b1 + 2*b2. vector_basis = basis degree = 0 @@ -338,16 +338,16 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: set_random_seed() sage: J = random_eja() sage: n = J.dimension() - sage: ei = J.zero() - sage: ej = J.zero() - sage: ei_ej = J.zero()*J.zero() + sage: bi = J.zero() + sage: bj = J.zero() + sage: bi_bj = J.zero()*J.zero() sage: if n > 0: ....: i = ZZ.random_element(n) ....: j = ZZ.random_element(n) - ....: ei = J.monomial(i) - ....: ej = J.monomial(j) - ....: ei_ej = J.product_on_basis(i,j) - sage: ei*ej == ei_ej + ....: bi = J.monomial(i) + ....: bj = J.monomial(j) + ....: bi_bj = J.product_on_basis(i,j) + sage: bi*bj == bi_bj True """ @@ -619,7 +619,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: J2 = RealSymmetricEJA(2) sage: J = cartesian_product([J1,J2]) sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) ) - e1 + e5 + b1 + b5 TESTS: @@ -894,15 +894,15 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: J = JordanSpinEJA(4) sage: J.multiplication_table() +----++----+----+----+----+ - | * || e0 | e1 | e2 | e3 | + | * || b0 | b1 | b2 | b3 | +====++====+====+====+====+ - | e0 || e0 | e1 | e2 | e3 | + | b0 || b0 | b1 | b2 | b3 | +----++----+----+----+----+ - | e1 || e1 | e0 | 0 | 0 | + | b1 || b1 | b0 | 0 | 0 | +----++----+----+----+----+ - | e2 || e2 | 0 | e0 | 0 | + | b2 || b2 | 0 | b0 | 0 | +----++----+----+----+----+ - | e3 || e3 | 0 | 0 | e0 | + | b3 || b3 | 0 | 0 | b0 | +----++----+----+----+----+ """ @@ -956,7 +956,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: J = RealSymmetricEJA(2) sage: J.basis() - Finite family {0: e0, 1: e1, 2: e2} + Finite family {0: b0, 1: b1, 2: b2} sage: J.matrix_basis() ( [1 0] [ 0 0.7071067811865475?] [0 0] @@ -967,7 +967,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: J = JordanSpinEJA(2) sage: J.basis() - Finite family {0: e0, 1: e1} + Finite family {0: b0, 1: b1} sage: J.matrix_basis() ( [1] [0] @@ -1024,7 +1024,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(): @@ -1049,20 +1051,20 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: J = HadamardEJA(5) sage: J.one() - e0 + e1 + e2 + e3 + e4 + b0 + b1 + b2 + b3 + b4 The unit element in the Hadamard EJA is inherited in the subalgebras generated by its elements:: sage: J = HadamardEJA(5) sage: J.one() - e0 + e1 + e2 + e3 + e4 + b0 + b1 + b2 + b3 + b4 sage: x = sum(J.gens()) sage: A = x.subalgebra_generated_by(orthonormalize=False) sage: A.one() - f0 + c0 sage: A.one().superalgebra_element() - e0 + e1 + e2 + e3 + e4 + b0 + b1 + b2 + b3 + b4 TESTS: @@ -1549,7 +1551,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): class RationalBasisEJA(FiniteDimensionalEJA): r""" - New class for algebras whose supplied basis elements have all rational entries. + Algebras whose supplied basis elements have all rational entries. SETUP:: @@ -1580,7 +1582,11 @@ class RationalBasisEJA(FiniteDimensionalEJA): if check_field: # Abuse the check_field parameter to check that the entries of # out basis (in ambient coordinates) are in the field QQ. - if not all( all(b_i in QQ for b_i in b.list()) for b in basis ): + # Use _all2list to get the vector coordinates of octonion + # entries and not the octonions themselves (which are not + # rational). + if not all( all(b_i in QQ for b_i in _all2list(b)) + for b in basis ): raise TypeError("basis not rational") super().__init__(basis, @@ -1662,7 +1668,7 @@ class RationalBasisEJA(FiniteDimensionalEJA): subs_dict = { X[i]: BX[i] for i in range(len(X)) } return tuple( a_i.subs(subs_dict) for a_i in a ) -class ConcreteEJA(RationalBasisEJA): +class ConcreteEJA(FiniteDimensionalEJA): r""" A class for the Euclidean Jordan algebras that we know by name. @@ -1732,132 +1738,147 @@ class ConcreteEJA(RationalBasisEJA): class MatrixEJA: @staticmethod - def dimension_over_reals(): - r""" - The dimension of this matrix's base ring over the reals. + def _denormalized_basis(A): + """ + Returns a basis for the space of complex Hermitian n-by-n matrices. - The reals are dimension one over themselves, obviously; that's - just `\mathbb{R}^{1}`. Likewise, the complex numbers `a + bi` - have dimension two. Finally, the quaternions have dimension - four over the reals. + Why do we embed these? Basically, because all of numerical linear + algebra assumes that you're working with vectors consisting of `n` + entries from a field and scalars from the same field. There's no way + to tell SageMath that (for example) the vectors contain complex + numbers, while the scalar field is real. - This is used to determine the size of the matrix returned from - :meth:`real_embed`, among other things. - """ - raise NotImplementedError + SETUP:: - @classmethod - def real_embed(cls,M): - """ - Embed the matrix ``M`` into a space of real matrices. + sage: from mjo.hurwitz import (ComplexMatrixAlgebra, + ....: QuaternionMatrixAlgebra, + ....: OctonionMatrixAlgebra) + sage: from mjo.eja.eja_algebra import MatrixEJA + + TESTS:: + + sage: set_random_seed() + sage: n = ZZ.random_element(1,5) + sage: A = MatrixSpace(QQ, n) + sage: B = MatrixEJA._denormalized_basis(A) + sage: all( M.is_hermitian() for M in B) + True - 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) + sage: set_random_seed() + sage: n = ZZ.random_element(1,5) + sage: A = ComplexMatrixAlgebra(n, scalars=QQ) + sage: B = MatrixEJA._denormalized_basis(A) + sage: all( M.is_hermitian() for M in B) + True - """ - if M.ncols() != M.nrows(): - raise ValueError("the matrix 'M' must be square") - return M + :: + sage: set_random_seed() + sage: n = ZZ.random_element(1,5) + sage: A = QuaternionMatrixAlgebra(n, scalars=QQ) + sage: B = MatrixEJA._denormalized_basis(A) + sage: all( M.is_hermitian() for M in B ) + True + + :: + + sage: set_random_seed() + sage: n = ZZ.random_element(1,5) + sage: A = OctonionMatrixAlgebra(n, scalars=QQ) + sage: B = MatrixEJA._denormalized_basis(A) + sage: all( M.is_hermitian() for M in B ) + True - @classmethod - def real_unembed(cls,M): - """ - The inverse of :meth:`real_embed`. """ - if M.ncols() != M.nrows(): - raise ValueError("the matrix 'M' must be square") - if not ZZ(M.nrows()).mod(cls.dimension_over_reals()).is_zero(): - raise ValueError("the matrix 'M' must be a real embedding") - return M + # These work for real MatrixSpace, whose monomials only have + # two coordinates (because the last one would always be "1"). + es = A.base_ring().gens() + gen = lambda A,m: A.monomial(m[:2]) + + if hasattr(A, 'entry_algebra_gens'): + # We've got a MatrixAlgebra, and its monomials will have + # three coordinates. + es = A.entry_algebra_gens() + gen = lambda A,m: A.monomial(m) + + basis = [] + for i in range(A.nrows()): + for j in range(i+1): + if i == j: + E_ii = gen(A, (i,j,es[0])) + basis.append(E_ii) + else: + for e in es: + E_ij = gen(A, (i,j,e)) + E_ij += E_ij.conjugate_transpose() + basis.append(E_ij) + + return tuple( basis ) @staticmethod def jordan_product(X,Y): return (X*Y + Y*X)/2 - @classmethod - def trace_inner_product(cls,X,Y): + @staticmethod + def trace_inner_product(X,Y): r""" - Compute the trace inner-product of two real-embeddings. + A trace inner-product for matrices that aren't embedded in the + reals. It takes MATRICES as arguments, not EJA elements. SETUP:: sage: from mjo.eja.eja_algebra import (RealSymmetricEJA, ....: ComplexHermitianEJA, - ....: QuaternionHermitianEJA) + ....: QuaternionHermitianEJA, + ....: OctonionHermitianEJA) EXAMPLES:: - This gives the same answer as it would if we computed the trace - from the unembedded (original) matrices:: + sage: J = RealSymmetricEJA(2,field=QQ,orthonormalize=False) + sage: I = J.one().to_matrix() + sage: J.trace_inner_product(I, -I) + -2 - sage: set_random_seed() - sage: J = RealSymmetricEJA.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() - sage: actual = J.trace_inner_product(Xe,Ye) - sage: actual == expected - True + :: + + sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False) + sage: I = J.one().to_matrix() + sage: J.trace_inner_product(I, -I) + -2 :: - sage: set_random_seed() - sage: J = ComplexHermitianEJA.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().real() - sage: actual = J.trace_inner_product(Xe,Ye) - sage: actual == expected - True + sage: J = QuaternionHermitianEJA(2,field=QQ,orthonormalize=False) + sage: I = J.one().to_matrix() + sage: J.trace_inner_product(I, -I) + -2 :: - 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 + sage: J = OctonionHermitianEJA(2,field=QQ,orthonormalize=False) + sage: I = J.one().to_matrix() + sage: J.trace_inner_product(I, -I) + -2 """ - Xu = cls.real_unembed(X) - Yu = cls.real_unembed(Y) - tr = (Xu*Yu).trace() - - try: - # Works in QQ, AA, RDF, et cetera. - return tr.real() - except AttributeError: - # A quaternion doesn't have a real() method, but does - # have coefficient_tuple() method that returns the - # coefficients of 1, i, j, and k -- in that order. + tr = (X*Y).trace() + if hasattr(tr, 'coefficient'): + # Works for octonions, and has to come first because they + # also have a "real()" method that doesn't return an + # element of the scalar ring. + return tr.coefficient(0) + elif hasattr(tr, 'coefficient_tuple'): + # Works for quaternions. return tr.coefficient_tuple()[0] + # Works for real and complex numbers. + return tr.real() -class RealMatrixEJA(MatrixEJA): - @staticmethod - def dimension_over_reals(): - return 1 -class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA): +class RealSymmetricEJA(RationalBasisEJA, ConcreteEJA, MatrixEJA): """ The rank-n simple EJA consisting of real symmetric n-by-n matrices, the usual symmetric Jordan product, and the trace inner @@ -1870,13 +1891,13 @@ class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA): EXAMPLES:: sage: J = RealSymmetricEJA(2) - sage: e0, e1, e2 = J.gens() - sage: e0*e0 - e0 - sage: e1*e1 - 1/2*e0 + 1/2*e2 - sage: e2*e2 - e2 + sage: b0, b1, b2 = J.gens() + sage: b0*b0 + b0 + sage: b1*b1 + 1/2*b0 + 1/2*b2 + sage: b2*b2 + b2 In theory, our "field" can be any subfield of the reals:: @@ -1922,38 +1943,6 @@ class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA): Euclidean Jordan algebra of dimension 0 over Algebraic Real Field """ - @classmethod - def _denormalized_basis(cls, n): - """ - Return a basis for the space of real symmetric n-by-n matrices. - - SETUP:: - - sage: from mjo.eja.eja_algebra import RealSymmetricEJA - - TESTS:: - - sage: set_random_seed() - sage: n = ZZ.random_element(1,5) - sage: B = RealSymmetricEJA._denormalized_basis(n) - sage: all( M.is_symmetric() for M in B) - True - - """ - # The basis of symmetric matrices, as matrices, in their R^(n-by-n) - # coordinates. - S = [] - for i in range(n): - for j in range(i+1): - Eij = matrix(ZZ, n, lambda k,l: k==i and l==j) - if i == j: - Sij = Eij - else: - Sij = Eij + Eij.transpose() - S.append(Sij) - return tuple(S) - - @staticmethod def _max_random_instance_size(): return 4 # Dimension 10 @@ -1966,179 +1955,27 @@ class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA): n = ZZ.random_element(cls._max_random_instance_size() + 1) return cls(n, **kwargs) - def __init__(self, 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 - associative = False - if n <= 1: - associative = True - - super().__init__(self._denormalized_basis(n), + A = MatrixSpace(field, n) + super().__init__(self._denormalized_basis(A), self.jordan_product, self.trace_inner_product, - associative=associative, + field=field, **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(). self.rank.set_cache(n) - idV = matrix.identity(ZZ, self.dimension_over_reals()*n) - self.one.set_cache(self(idV)) - - - -class ComplexMatrixEJA(MatrixEJA): - # A manual dictionary-cache for the complex_extension() method, - # since apparently @classmethods can't also be @cached_methods. - _complex_extension = {} - - @classmethod - def complex_extension(cls,field): - r""" - The complex field that we embed/unembed, as an extension - of the given ``field``. - """ - if field in cls._complex_extension: - return cls._complex_extension[field] - - # Sage doesn't know how to adjoin the complex "i" (the root of - # x^2 + 1) to a field in a general way. Here, we just enumerate - # all of the cases that I have cared to support so far. - if field is AA: - # Sage doesn't know how to embed AA into QQbar, i.e. how - # to adjoin sqrt(-1) to AA. - F = QQbar - elif not field.is_exact(): - # RDF or RR - F = field.complex_field() - else: - # Works for QQ and... maybe some other fields. - R = PolynomialRing(field, 'z') - z = R.gen() - F = field.extension(z**2 + 1, 'I', embedding=CLF(-1).sqrt()) - - cls._complex_extension[field] = F - return F - - @staticmethod - def dimension_over_reals(): - return 2 - - @classmethod - def real_embed(cls,M): - """ - Embed the n-by-n complex matrix ``M`` into the space of real - matrices of size 2n-by-2n via the map the sends each entry `z = a + - bi` to the block matrix ``[[a,b],[-b,a]]``. - - SETUP:: - - sage: from mjo.eja.eja_algebra import ComplexMatrixEJA - - EXAMPLES:: - - sage: F = QuadraticField(-1, 'I') - sage: x1 = F(4 - 2*i) - sage: x2 = F(1 + 2*i) - sage: x3 = F(-i) - sage: x4 = F(6) - sage: M = matrix(F,2,[[x1,x2],[x3,x4]]) - sage: ComplexMatrixEJA.real_embed(M) - [ 4 -2| 1 2] - [ 2 4|-2 1] - [-----+-----] - [ 0 -1| 6 0] - [ 1 0| 0 6] - - TESTS: - - Embedding is a homomorphism (isomorphism, in fact):: - - sage: set_random_seed() - sage: n = ZZ.random_element(3) - sage: F = QuadraticField(-1, 'I') - sage: X = random_matrix(F, n) - sage: Y = random_matrix(F, n) - sage: Xe = ComplexMatrixEJA.real_embed(X) - sage: Ye = ComplexMatrixEJA.real_embed(Y) - sage: XYe = ComplexMatrixEJA.real_embed(X*Y) - sage: Xe*Ye == XYe - True - - """ - super().real_embed(M) - n = M.nrows() - - # We don't need any adjoined elements... - field = M.base_ring().base_ring() - - blocks = [] - for z in M.list(): - a = z.real() - b = z.imag() - blocks.append(matrix(field, 2, [ [ a, b], - [-b, a] ])) - - return matrix.block(field, n, blocks) + self.one.set_cache(self(A.one())) - @classmethod - def real_unembed(cls,M): - """ - The inverse of _embed_complex_matrix(). - - SETUP:: - - sage: from mjo.eja.eja_algebra import ComplexMatrixEJA - - EXAMPLES:: - - sage: A = matrix(QQ,[ [ 1, 2, 3, 4], - ....: [-2, 1, -4, 3], - ....: [ 9, 10, 11, 12], - ....: [-10, 9, -12, 11] ]) - sage: ComplexMatrixEJA.real_unembed(A) - [ 2*I + 1 4*I + 3] - [ 10*I + 9 12*I + 11] - - TESTS: - - Unembedding is the inverse of embedding:: - sage: set_random_seed() - sage: F = QuadraticField(-1, 'I') - sage: M = random_matrix(F, 3) - sage: Me = ComplexMatrixEJA.real_embed(M) - sage: ComplexMatrixEJA.real_unembed(Me) == M - True - - """ - super().real_unembed(M) - n = ZZ(M.nrows()) - d = cls.dimension_over_reals() - F = cls.complex_extension(M.base_ring()) - i = F.gen() - - # Go top-left to bottom-right (reading order), converting every - # 2-by-2 block we see to a single complex element. - elements = [] - for k in range(n/d): - for j in range(n/d): - submat = M[d*k:d*k+d,d*j:d*j+d] - if submat[0,0] != submat[1,1]: - raise ValueError('bad on-diagonal submatrix') - if submat[0,1] != -submat[1,0]: - raise ValueError('bad off-diagonal submatrix') - z = submat[0,0] + submat[0,1]*i - elements.append(z) - - return matrix(F, n/d, elements) - - -class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA): +class ComplexHermitianEJA(RationalBasisEJA, ConcreteEJA, MatrixEJA): """ The rank-n simple EJA consisting of complex Hermitian n-by-n matrices over the real numbers, the usual symmetric Jordan product, @@ -2195,90 +2032,24 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA): Euclidean Jordan algebra of dimension 0 over Algebraic Real Field """ - - @classmethod - def _denormalized_basis(cls, n): - """ - Returns a basis for the space of complex Hermitian n-by-n matrices. - - Why do we embed these? Basically, because all of numerical linear - algebra assumes that you're working with vectors consisting of `n` - entries from a field and scalars from the same field. There's no way - to tell SageMath that (for example) the vectors contain complex - numbers, while the scalar field is real. - - SETUP:: - - sage: from mjo.eja.eja_algebra import ComplexHermitianEJA - - TESTS:: - - sage: set_random_seed() - sage: n = ZZ.random_element(1,5) - sage: B = ComplexHermitianEJA._denormalized_basis(n) - sage: all( M.is_symmetric() for M in B) - True - - """ - field = ZZ - R = PolynomialRing(field, 'z') - z = R.gen() - F = field.extension(z**2 + 1, 'I') - I = F.gen(1) - - # 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(F,n) - 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) - else: - # The second one has a minus because it's conjugated. - Eij[j,i] = 1 # Eij = Eij + Eij.transpose() - Sij_real = cls.real_embed(Eij) - S.append(Sij_real) - # Eij = I*Eij - I*Eij.transpose() - Eij[i,j] = I - Eij[j,i] = -I - Sij_imag = cls.real_embed(Eij) - S.append(Sij_imag) - Eij[j,i] = 0 - # "erase" E_ij - Eij[i,j] = 0 - - # 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 ) - - - def __init__(self, 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 - associative = False - if n <= 1: - associative = True - - super().__init__(self._denormalized_basis(n), + from mjo.hurwitz import ComplexMatrixAlgebra + A = ComplexMatrixAlgebra(n, scalars=field) + super().__init__(self._denormalized_basis(A), self.jordan_product, self.trace_inner_product, - associative=associative, + field=field, **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(). self.rank.set_cache(n) - idV = matrix.identity(ZZ, self.dimension_over_reals()*n) - self.one.set_cache(self(idV)) + self.one.set_cache(self(A.one())) @staticmethod def _max_random_instance_size(): @@ -2292,155 +2063,8 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA): n = ZZ.random_element(cls._max_random_instance_size() + 1) return cls(n, **kwargs) -class QuaternionMatrixEJA(MatrixEJA): - - # 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(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 @@ -2497,100 +2121,24 @@ class QuaternionHermitianEJA(ConcreteEJA, QuaternionMatrixEJA): Euclidean Jordan algebra of dimension 0 over Algebraic Real Field """ - @classmethod - def _denormalized_basis(cls, n): - """ - Returns a basis for the space of quaternion Hermitian n-by-n matrices. - - Why do we embed these? Basically, because all of numerical - linear algebra assumes that you're working with vectors consisting - of `n` entries from a field and scalars from the same field. There's - no way to tell SageMath that (for example) the vectors contain - complex numbers, while the scalar field is real. - - SETUP:: - - sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA - - TESTS:: - - sage: set_random_seed() - sage: n = ZZ.random_element(1,5) - sage: B = QuaternionHermitianEJA._denormalized_basis(n) - sage: all( M.is_symmetric() for M in B ) - True - - """ - field = ZZ - Q = QuaternionAlgebra(QQ,-1,-1) - I,J,K = Q.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) - 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) - 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 - - # Since we embedded these, we can drop back to the "field" that we - # started with instead of the quaternion algebra "Q". - return tuple( s.change_ring(field) for s in S ) - - - def __init__(self, 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 - associative = False - if n <= 1: - associative = True - - super().__init__(self._denormalized_basis(n), + from mjo.hurwitz import QuaternionMatrixAlgebra + A = QuaternionMatrixAlgebra(n, scalars=field) + super().__init__(self._denormalized_basis(A), self.jordan_product, self.trace_inner_product, - associative=associative, + field=field, **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(). self.rank.set_cache(n) - idV = matrix.identity(ZZ, self.dimension_over_reals()*n) - self.one.set_cache(self(idV)) + self.one.set_cache(self(A.one())) @staticmethod @@ -2608,15 +2156,160 @@ class QuaternionHermitianEJA(ConcreteEJA, QuaternionMatrixEJA): n = ZZ.random_element(cls._max_random_instance_size() + 1) return cls(n, **kwargs) +class OctonionHermitianEJA(RationalBasisEJA, ConcreteEJA, MatrixEJA): + r""" + SETUP:: + + sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA, + ....: OctonionHermitianEJA) + + EXAMPLES: + + The 3-by-3 algebra satisfies the axioms of an EJA:: + + sage: OctonionHermitianEJA(3, # long time + ....: field=QQ, # long time + ....: orthonormalize=False, # long time + ....: check_axioms=True) # long time + Euclidean Jordan algebra of dimension 27 over Rational Field + + After a change-of-basis, the 2-by-2 algebra has the same + multiplication table as the ten-dimensional Jordan spin algebra:: + + sage: b = OctonionHermitianEJA._denormalized_basis(2,QQ) + sage: basis = (b[0] + b[9],) + b[1:9] + (b[0] - b[9],) + sage: jp = OctonionHermitianEJA.jordan_product + sage: ip = OctonionHermitianEJA.trace_inner_product + sage: J = FiniteDimensionalEJA(basis, + ....: jp, + ....: ip, + ....: field=QQ, + ....: orthonormalize=False) + sage: J.multiplication_table() + +----++----+----+----+----+----+----+----+----+----+----+ + | * || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 | + +====++====+====+====+====+====+====+====+====+====+====+ + | b0 || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b1 || b1 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b2 || b2 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b3 || b3 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b4 || b4 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b5 || b5 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b6 || b6 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b7 || b7 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b8 || b8 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | + +----++----+----+----+----+----+----+----+----+----+----+ + | b9 || b9 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | + +----++----+----+----+----+----+----+----+----+----+----+ + + TESTS: + + We can actually construct the 27-dimensional Albert algebra, + and we get the right unit element if we recompute it:: + + sage: J = OctonionHermitianEJA(3, # long time + ....: field=QQ, # long time + ....: orthonormalize=False) # long time + sage: J.one.clear_cache() # long time + sage: J.one() # long time + b0 + b9 + b26 + sage: J.one().to_matrix() # long time + +----+----+----+ + | e0 | 0 | 0 | + +----+----+----+ + | 0 | e0 | 0 | + +----+----+----+ + | 0 | 0 | e0 | + +----+----+----+ + + The 2-by-2 algebra is isomorphic to the ten-dimensional Jordan + spin algebra, but just to be sure, we recompute its rank:: + + sage: J = OctonionHermitianEJA(2, # long time + ....: field=QQ, # long time + ....: orthonormalize=False) # long time + sage: J.rank.clear_cache() # long time + sage: J.rank() # long time + 2 -class HadamardEJA(ConcreteEJA): """ - Return the Euclidean Jordan Algebra corresponding to the set - `R^n` under the Hadamard product. + @staticmethod + def _max_random_instance_size(): + r""" + The maximum rank of a random QuaternionHermitianEJA. + """ + return 1 # Dimension 1 + + @classmethod + def random_instance(cls, **kwargs): + """ + Return a random instance of this type of algebra. + """ + n = ZZ.random_element(cls._max_random_instance_size() + 1) + return cls(n, **kwargs) - Note: this is nothing more than the Cartesian product of ``n`` - copies of the spin algebra. Once Cartesian product algebras - are implemented, this can go. + def __init__(self, n, field=AA, **kwargs): + if n > 3: + # Otherwise we don't get an EJA. + raise ValueError("n cannot exceed 3") + + # 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 OctonionMatrixAlgebra + A = OctonionMatrixAlgebra(n, scalars=field) + super().__init__(self._denormalized_basis(A), + self.jordan_product, + self.trace_inner_product, + field=field, + **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(). + self.rank.set_cache(n) + self.one.set_cache(self(A.one())) + + +class AlbertEJA(OctonionHermitianEJA): + r""" + The Albert algebra is the algebra of three-by-three Hermitian + matrices whose entries are octonions. + + SETUP:: + + sage: from mjo.eja.eja_algebra import AlbertEJA + + EXAMPLES:: + + sage: AlbertEJA(field=QQ, orthonormalize=False) + Euclidean Jordan algebra of dimension 27 over Rational Field + sage: AlbertEJA() # long time + Euclidean Jordan algebra of dimension 27 over Algebraic Real Field + + """ + def __init__(self, *args, **kwargs): + super().__init__(3, *args, **kwargs) + + +class HadamardEJA(RationalBasisEJA, ConcreteEJA): + """ + Return the Euclidean Jordan algebra on `R^n` with the Hadamard + (pointwise real-number multiplication) Jordan product and the + usual inner-product. + + This is nothing more than the Cartesian product of ``n`` copies of + the one-dimensional Jordan spin algebra, and is the most common + example of a non-simple Euclidean Jordan algebra. SETUP:: @@ -2627,19 +2320,19 @@ class HadamardEJA(ConcreteEJA): This multiplication table can be verified by hand:: sage: J = HadamardEJA(3) - sage: e0,e1,e2 = J.gens() - sage: e0*e0 - e0 - sage: e0*e1 + sage: b0,b1,b2 = J.gens() + sage: b0*b0 + b0 + sage: b0*b1 0 - sage: e0*e2 + sage: b0*b2 0 - sage: e1*e1 - e1 - sage: e1*e2 + sage: b1*b1 + b1 + sage: b1*b2 0 - sage: e2*e2 - e2 + sage: b2*b2 + b2 TESTS: @@ -2647,9 +2340,8 @@ class HadamardEJA(ConcreteEJA): sage: HadamardEJA(3, prefix='r').gens() (r0, r1, r2) - """ - def __init__(self, n, **kwargs): + def __init__(self, n, field=AA, **kwargs): if n == 0: jordan_product = lambda x,y: x inner_product = lambda x,y: x @@ -2670,10 +2362,12 @@ class HadamardEJA(ConcreteEJA): if "orthonormalize" not in kwargs: kwargs["orthonormalize"] = False if "check_axioms" not in kwargs: kwargs["check_axioms"] = False - column_basis = tuple( b.column() for b in FreeModule(ZZ, n).basis() ) + column_basis = tuple( b.column() + for b in FreeModule(field, n).basis() ) super().__init__(column_basis, jordan_product, inner_product, + field=field, associative=True, **kwargs) self.rank.set_cache(n) @@ -2699,7 +2393,7 @@ class HadamardEJA(ConcreteEJA): return cls(n, **kwargs) -class BilinearFormEJA(ConcreteEJA): +class BilinearFormEJA(RationalBasisEJA, ConcreteEJA): r""" The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)`` with the half-trace inner product and jordan product ``x*y = @@ -2781,7 +2475,7 @@ class BilinearFormEJA(ConcreteEJA): True """ - def __init__(self, B, **kwargs): + def __init__(self, B, field=AA, **kwargs): # The matrix "B" is supplied by the user in most cases, # so it makes sense to check whether or not its positive- # definite unless we are specifically asked not to... @@ -2809,7 +2503,8 @@ class BilinearFormEJA(ConcreteEJA): return P([z0] + zbar.list()) n = B.nrows() - column_basis = tuple( b.column() for b in FreeModule(ZZ, n).basis() ) + column_basis = tuple( b.column() + for b in FreeModule(field, n).basis() ) # TODO: I haven't actually checked this, but it seems legit. associative = False @@ -2819,6 +2514,7 @@ class BilinearFormEJA(ConcreteEJA): super().__init__(column_basis, jordan_product, inner_product, + field=field, associative=associative, **kwargs) @@ -2880,20 +2576,20 @@ class JordanSpinEJA(BilinearFormEJA): This multiplication table can be verified by hand:: sage: J = JordanSpinEJA(4) - sage: e0,e1,e2,e3 = J.gens() - sage: e0*e0 - e0 - sage: e0*e1 - e1 - sage: e0*e2 - e2 - sage: e0*e3 - e3 - sage: e1*e2 + sage: b0,b1,b2,b3 = J.gens() + sage: b0*b0 + b0 + sage: b0*b1 + b1 + sage: b0*b2 + b2 + sage: b0*b3 + b3 + sage: b1*b2 0 - sage: e1*e3 + sage: b1*b3 0 - sage: e2*e3 + sage: b2*b3 0 We can change the generator prefix:: @@ -2914,7 +2610,7 @@ class JordanSpinEJA(BilinearFormEJA): True """ - def __init__(self, n, **kwargs): + def __init__(self, n, *args, **kwargs): # This is a special case of the BilinearFormEJA with the # identity matrix as its bilinear form. B = matrix.identity(ZZ, n) @@ -2925,7 +2621,7 @@ class JordanSpinEJA(BilinearFormEJA): # But also don't pass check_field=False here, because the user # can pass in a field! - super().__init__(B, **kwargs) + super().__init__(B, *args, **kwargs) @staticmethod def _max_random_instance_size(): @@ -2945,7 +2641,7 @@ class JordanSpinEJA(BilinearFormEJA): return cls(n, **kwargs) -class TrivialEJA(ConcreteEJA): +class TrivialEJA(RationalBasisEJA, ConcreteEJA): """ The trivial Euclidean Jordan algebra consisting of only a zero element. @@ -3105,23 +2801,23 @@ class CartesianProductEJA(FiniteDimensionalEJA): sage: J = cartesian_product([J1,cartesian_product([J2,J3])]) sage: J.multiplication_table() +----++----+----+----+ - | * || e0 | e1 | e2 | + | * || b0 | b1 | b2 | +====++====+====+====+ - | e0 || e0 | 0 | 0 | + | b0 || b0 | 0 | 0 | +----++----+----+----+ - | e1 || 0 | e1 | 0 | + | b1 || 0 | b1 | 0 | +----++----+----+----+ - | e2 || 0 | 0 | e2 | + | b2 || 0 | 0 | b2 | +----++----+----+----+ sage: HadamardEJA(3).multiplication_table() +----++----+----+----+ - | * || e0 | e1 | e2 | + | * || b0 | b1 | b2 | +====++====+====+====+ - | e0 || e0 | 0 | 0 | + | b0 || b0 | 0 | 0 | +----++----+----+----+ - | e1 || 0 | e1 | 0 | + | b1 || 0 | b1 | 0 | +----++----+----+----+ - | e2 || 0 | 0 | e2 | + | b2 || 0 | 0 | b2 | +----++----+----+----+ TESTS: @@ -3230,9 +2926,17 @@ class CartesianProductEJA(FiniteDimensionalEJA): 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 (HadamardEJA, + sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA, + ....: HadamardEJA, + ....: OctonionHermitianEJA, ....: RealSymmetricEJA) EXAMPLES:: @@ -3245,10 +2949,44 @@ class CartesianProductEJA(FiniteDimensionalEJA): 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 0] + [0 1] + sage: J.one().to_matrix()[1] + [1 0] + [0 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 | + +----+ + """ - from sage.categories.cartesian_product import cartesian_product - return cartesian_product( [J.matrix_space() - for J in self.cartesian_factors()] ) + scalars = self.cartesian_factor(0).base_ring() + + # This category isn't perfect, but is good enough for what we + # need to do. + cat = MagmaticAlgebras(scalars).FiniteDimensional().WithBasis() + cat = cat.Unital().CartesianProducts() + factors = tuple( J.matrix_space() for J in self.cartesian_factors() ) + + from sage.sets.cartesian_product import CartesianProduct + return CartesianProduct(factors, cat) + @cached_method def cartesian_projection(self, i): @@ -3450,7 +3188,9 @@ class RationalBasisCartesianProductEJA(CartesianProductEJA, SETUP:: - sage: from mjo.eja.eja_algebra import (JordanSpinEJA, + sage: from mjo.eja.eja_algebra import (HadamardEJA, + ....: JordanSpinEJA, + ....: OctonionHermitianEJA, ....: RealSymmetricEJA) EXAMPLES: @@ -3467,30 +3207,45 @@ class RationalBasisCartesianProductEJA(CartesianProductEJA, sage: J.rank() 5 + TESTS: + + 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 + """ def __init__(self, algebras, **kwargs): CartesianProductEJA.__init__(self, algebras, **kwargs) self._rational_algebra = None if self.vector_space().base_field() is not QQ: - self._rational_algebra = cartesian_product([ - r._rational_algebra for r in algebras - ]) + if all( hasattr(r, "_rational_algebra") for r in algebras ): + self._rational_algebra = cartesian_product([ + r._rational_algebra for r in algebras + ]) RationalBasisEJA.CartesianProduct = RationalBasisCartesianProductEJA -random_eja = ConcreteEJA.random_instance - -# def random_eja(*args, **kwargs): -# J1 = ConcreteEJA.random_instance(*args, **kwargs) - -# # This might make Cartesian products appear roughly as often as -# # any other ConcreteEJA. -# if ZZ.random_element(len(ConcreteEJA.__subclasses__()) + 1) == 0: -# # Use random_eja() again so we can get more than two factors. -# J2 = random_eja(*args, **kwargs) -# J = cartesian_product([J1,J2]) -# return J -# else: -# return J1 +def random_eja(*args, **kwargs): + J1 = ConcreteEJA.random_instance(*args, **kwargs) + + # This might make Cartesian products appear roughly as often as + # any other ConcreteEJA. + if ZZ.random_element(len(ConcreteEJA.__subclasses__()) + 1) == 0: + # Use random_eja() again so we can get more than two factors. + J2 = random_eja(*args, **kwargs) + J = cartesian_product([J1,J2]) + return J + else: + return J1