X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_algebra.py;h=b155ed757e9e855ea9d18cdf5c0d1f1cb0878d4a;hb=e46ec3e9e875c496eeaad87585f0dd893c6ac71a;hp=f14218151d69db6fb90bdf2bec240a9506569f33;hpb=bbae1a93ab44c80896181ad8983393aa16a87868;p=sage.d.git diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index f142181..b155ed7 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -20,7 +20,9 @@ from itertools import repeat from sage.algebras.quatalg.quaternion_algebra import QuaternionAlgebra from sage.categories.magmatic_algebras import MagmaticAlgebras -from sage.combinat.free_module import CombinatorialFreeModule +from sage.categories.sets_cat import cartesian_product +from sage.combinat.free_module import (CombinatorialFreeModule, + CombinatorialFreeModule_CartesianProduct) from sage.matrix.constructor import matrix from sage.matrix.matrix_space import MatrixSpace from sage.misc.cachefunc import cached_method @@ -133,7 +135,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): deortho_vector_basis = tuple( V(b.list()) for b in basis ) from mjo.eja.eja_utils import gram_schmidt - basis = gram_schmidt(basis, inner_product) + basis = tuple(gram_schmidt(basis, inner_product)) # Save the (possibly orthonormalized) matrix basis for # later... @@ -158,7 +160,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # Now we actually compute the multiplication and inner-product # tables/matrices using the possibly-orthonormalized basis. - self._inner_product_matrix = matrix.zero(field, n) + self._inner_product_matrix = matrix.identity(field, n) self._multiplication_table = [ [0 for j in range(i+1)] for i in range(n) ] @@ -171,15 +173,20 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): q_i = basis[i] q_j = basis[j] - elt = jordan_product(q_i, q_j) - ip = inner_product(q_i, q_j) - # The jordan product returns a matrixy answer, so we # have to convert it to the algebra coordinates. + elt = jordan_product(q_i, q_j) elt = W.coordinate_vector(V(elt.list())) self._multiplication_table[i][j] = self.from_vector(elt) - self._inner_product_matrix[i,j] = ip - self._inner_product_matrix[j,i] = ip + + if not orthonormalize: + # If we're orthonormalizing the basis with respect + # to an inner-product, then the inner-product + # matrix with respect to the resulting basis is + # just going to be the identity. + ip = inner_product(q_i, q_j) + self._inner_product_matrix[i,j] = ip + self._inner_product_matrix[j,i] = ip self._inner_product_matrix._cache = {'hermitian': True} self._inner_product_matrix.set_immutable() @@ -315,7 +322,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): This method should of course always return ``True``, unless this algebra was constructed with ``check_axioms=False`` and - passed an invalid multiplication table. + passed an invalid Jordan or inner-product. """ # Used to check whether or not something is zero in an inexact @@ -682,7 +689,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): Why implement this for non-matrix algebras? Avoiding special cases for the :class:`BilinearFormEJA` pays with simplicity in its own right. But mainly, we would like to be able to assume - that elements of a :class:`DirectSumEJA` can be displayed + that elements of a :class:`CartesianProductEJA` can be displayed nicely, without having to have special classes for direct sums one of whose components was a matrix algebra. @@ -732,7 +739,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): if self.is_trivial(): return MatrixSpace(self.base_ring(), 0) else: - return self._matrix_basis[0].matrix_space() + return self.matrix_basis()[0].parent() @cached_method @@ -745,23 +752,57 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: from mjo.eja.eja_algebra import (HadamardEJA, ....: random_eja) - EXAMPLES:: + EXAMPLES: + + We can compute unit element in the Hadamard EJA:: sage: J = HadamardEJA(5) sage: J.one() e0 + e1 + e2 + e3 + e4 + 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 + sage: x = sum(J.gens()) + sage: A = x.subalgebra_generated_by(orthonormalize=False) + sage: A.one() + f0 + sage: A.one().superalgebra_element() + e0 + e1 + e2 + e3 + e4 + TESTS: - The identity element acts like the identity:: + The identity element acts like the identity, regardless of + whether or not we orthonormalize:: sage: set_random_seed() sage: J = random_eja() sage: x = J.random_element() sage: J.one()*x == x and x*J.one() == x True + sage: A = x.subalgebra_generated_by() + sage: y = A.random_element() + sage: A.one()*y == y and y*A.one() == y + True + + :: + + sage: set_random_seed() + sage: J = random_eja(field=QQ, orthonormalize=False) + sage: x = J.random_element() + sage: J.one()*x == x and x*J.one() == x + True + sage: A = x.subalgebra_generated_by(orthonormalize=False) + sage: y = A.random_element() + sage: A.one()*y == y and y*A.one() == y + True - The matrix of the unit element's operator is the identity:: + The matrix of the unit element's operator is the identity, + regardless of the base field and whether or not we + orthonormalize:: sage: set_random_seed() sage: J = random_eja() @@ -769,6 +810,27 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: expected = matrix.identity(J.base_ring(), J.dimension()) sage: actual == expected True + sage: x = J.random_element() + sage: A = x.subalgebra_generated_by() + sage: actual = A.one().operator().matrix() + sage: expected = matrix.identity(A.base_ring(), A.dimension()) + sage: actual == expected + True + + :: + + sage: set_random_seed() + sage: J = random_eja(field=QQ, orthonormalize=False) + sage: actual = J.one().operator().matrix() + sage: expected = matrix.identity(J.base_ring(), J.dimension()) + sage: actual == expected + True + sage: x = J.random_element() + sage: A = x.subalgebra_generated_by(orthonormalize=False) + sage: actual = A.one().operator().matrix() + sage: expected = matrix.identity(A.base_ring(), A.dimension()) + sage: actual == expected + True Ensure that the cached unit element (often precomputed by hand) agrees with the computed one:: @@ -780,6 +842,15 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: J.one() == cached True + :: + + sage: set_random_seed() + sage: J = random_eja(field=QQ, orthonormalize=False) + sage: cached = J.one() + sage: J.one.clear_cache() + sage: J.one() == cached + True + """ # We can brute-force compute the matrices of the operators # that correspond to the basis elements of this algebra. @@ -1032,6 +1103,21 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): r""" The `r` polynomial coefficients of the "characteristic polynomial of" function. + + SETUP:: + + sage: from mjo.eja.eja_algebra import random_eja + + TESTS: + + The theory shows that these are all homogeneous polynomials of + a known degree:: + + sage: set_random_seed() + sage: J = random_eja() + sage: all(p.is_homogeneous() for p in J._charpoly_coefficients()) + True + """ n = self.dimension() R = self.coordinate_polynomial_ring() @@ -1067,10 +1153,17 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): # The theory says that only the first "r" coefficients are # nonzero, and they actually live in the original polynomial - # ring and not the fraction field. We negate them because - # in the actual characteristic polynomial, they get moved - # to the other side where x^r lives. - return -A_rref.solve_right(E*b).change_ring(R)[:r] + # ring and not the fraction field. We negate them because in + # the actual characteristic polynomial, they get moved to the + # other side where x^r lives. We don't bother to trim A_rref + # down to a square matrix and solve the resulting system, + # because the upper-left r-by-r portion of A_rref is + # guaranteed to be the identity matrix, so e.g. + # + # A_rref.solve_right(Y) + # + # would just be returning Y. + return (-E*b)[:r].change_ring(R) @cached_method def rank(self): @@ -1131,7 +1224,7 @@ class FiniteDimensionalEJA(CombinatorialFreeModule): sage: set_random_seed() # long time sage: J = random_eja() # long time - sage: caches = J.rank() # long time + sage: cached = J.rank() # long time sage: J.rank.clear_cache() # long time sage: J.rank() == cached # long time True @@ -1187,9 +1280,7 @@ class RationalBasisEJA(FiniteDimensionalEJA): jordan_product, inner_product, field=AA, - orthonormalize=True, check_field=True, - check_axioms=True, **kwargs): if check_field: @@ -1198,6 +1289,7 @@ class RationalBasisEJA(FiniteDimensionalEJA): if not all( all(b_i in QQ for b_i in b.list()) for b in basis ): raise TypeError("basis not rational") + self._rational_algebra = None if field is not QQ: # There's no point in constructing the extra algebra if this # one is already rational. @@ -1212,15 +1304,13 @@ class RationalBasisEJA(FiniteDimensionalEJA): field=QQ, orthonormalize=False, check_field=False, - check_axioms=False, - **kwargs) + check_axioms=False) super().__init__(basis, jordan_product, inner_product, field=field, check_field=check_field, - check_axioms=check_axioms, **kwargs) @cached_method @@ -1260,7 +1350,14 @@ class RationalBasisEJA(FiniteDimensionalEJA): a = ( a_i.change_ring(self.base_ring()) for a_i in self._rational_algebra._charpoly_coefficients() ) - # Now convert the coordinate variables back to the + if self._deortho_matrix is None: + # This can happen if our base ring was, say, AA and we + # chose not to (or didn't need to) orthonormalize. It's + # still faster to do the computations over QQ even if + # the numbers in the boxes stay the same. + return tuple(a) + + # Otherwise, convert the coordinate variables back to the # deorthonormalized ones. R = self.coordinate_polynomial_ring() from sage.modules.free_module_element import vector @@ -1594,6 +1691,38 @@ class RealSymmetricEJA(ConcreteEJA, RealMatrixEJA): 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 @@ -1648,9 +1777,10 @@ class ComplexMatrixEJA(MatrixEJA): blocks = [] for z in M.list(): - a = z.list()[0] # real part, I guess - b = z.list()[1] # imag part, I guess - blocks.append(matrix(field, 2, [[a,b],[-b,a]])) + a = z.real() + b = z.imag() + blocks.append(matrix(field, 2, [ [ a, b], + [-b, a] ])) return matrix.block(field, n, blocks) @@ -1689,26 +1819,7 @@ class ComplexMatrixEJA(MatrixEJA): super(ComplexMatrixEJA,cls).real_unembed(M) n = ZZ(M.nrows()) d = cls.dimension_over_reals() - - # If "M" was normalized, its base ring might have roots - # adjoined and they can stick around after unembedding. - field = M.base_ring() - R = PolynomialRing(field, 'z') - z = R.gen() - - # 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. - F = field.extension(z**2 + 1, 'I', embedding=CLF(-1).sqrt()) + F = cls.complex_extension(M.base_ring()) i = F.gen() # Go top-left to bottom-right (reading order), converting every @@ -1804,7 +1915,6 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA): sage: set_random_seed() sage: n = ZZ.random_element(1,5) - sage: field = QuadraticField(2, 'sqrt2') sage: B = ComplexHermitianEJA._denormalized_basis(n) sage: all( M.is_symmetric() for M in B) True @@ -1822,18 +1932,27 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA): # * 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): - Eij = matrix(F, n, lambda k,l: k==i and l==j) + # "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. - Sij_real = cls.real_embed(Eij + Eij.transpose()) + Eij[j,i] = 1 # Eij = Eij + Eij.transpose() + Sij_real = cls.real_embed(Eij) S.append(Sij_real) - Sij_imag = cls.real_embed(I*Eij - I*Eij.transpose()) + # 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". @@ -1869,6 +1988,25 @@ class ComplexHermitianEJA(ConcreteEJA, ComplexMatrixEJA): 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 @@ -1973,8 +2111,7 @@ class QuaternionMatrixEJA(MatrixEJA): # Use the base ring of the matrix to ensure that its entries can be # multiplied by elements of the quaternion algebra. - field = M.base_ring() - Q = QuaternionAlgebra(field,-1,-1) + Q = cls.quaternion_extension(M.base_ring()) i,j,k = Q.gens() # Go top-left to bottom-right (reading order), converting every @@ -2089,23 +2226,39 @@ class QuaternionHermitianEJA(ConcreteEJA, QuaternionMatrixEJA): # * 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): - Eij = matrix(Q, n, lambda k,l: k==i and l==j) + # "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. - Sij_real = cls.real_embed(Eij + Eij.transpose()) + # Eij = Eij + Eij.transpose() + Eij[j,i] = 1 + Sij_real = cls.real_embed(Eij) S.append(Sij_real) - Sij_I = cls.real_embed(I*Eij - I*Eij.transpose()) + # Eij = I*(Eij - Eij.transpose()) + Eij[i,j] = I + Eij[j,i] = -I + Sij_I = cls.real_embed(Eij) S.append(Sij_I) - Sij_J = cls.real_embed(J*Eij - J*Eij.transpose()) + # Eij = J*(Eij - Eij.transpose()) + Eij[i,j] = J + Eij[j,i] = -J + Sij_J = cls.real_embed(Eij) S.append(Sij_J) - Sij_K = cls.real_embed(K*Eij - K*Eij.transpose()) + # 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". @@ -2311,31 +2464,38 @@ class BilinearFormEJA(ConcreteEJA): ....: for j in range(n-1) ] sage: actual == expected True + """ def __init__(self, B, **kwargs): - if not B.is_positive_definite(): - raise ValueError("bilinear form is not positive-definite") + # 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... + if ("check_axioms" not in kwargs) or kwargs["check_axioms"]: + if not B.is_positive_definite(): + raise ValueError("bilinear form is not positive-definite") + + # However, all of the other data for this EJA is computed + # by us in manner that guarantees the axioms are + # satisfied. So, again, unless we are specifically asked to + # verify things, we'll skip the rest of the checks. + if "check_axioms" not in kwargs: kwargs["check_axioms"] = False def inner_product(x,y): - return (B*x).inner_product(y) + return (y.T*B*x)[0,0] def jordan_product(x,y): P = x.parent() - x0 = x[0] - xbar = x[1:] - y0 = y[0] - ybar = y[1:] - z0 = inner_product(x,y) + x0 = x[0,0] + xbar = x[1:,0] + y0 = y[0,0] + ybar = y[1:,0] + z0 = inner_product(y,x) zbar = y0*xbar + x0*ybar - return P((z0,) + tuple(zbar)) - - # 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 + return P([z0] + zbar.list()) n = B.nrows() - standard_basis = FreeModule(ZZ, n).basis() - super(BilinearFormEJA, self).__init__(standard_basis, + column_basis = tuple( b.column() for b in FreeModule(ZZ, n).basis() ) + super(BilinearFormEJA, self).__init__(column_basis, jordan_product, inner_product, **kwargs) @@ -2516,100 +2676,214 @@ class TrivialEJA(ConcreteEJA): # inappropriate for us. return cls(**kwargs) -# class DirectSumEJA(ConcreteEJA): -# r""" -# The external (orthogonal) direct sum of two other Euclidean Jordan -# algebras. Essentially the Cartesian product of its two factors. -# Every Euclidean Jordan algebra decomposes into an orthogonal -# direct sum of simple Euclidean Jordan algebras, so no generality -# is lost by providing only this construction. - -# SETUP:: - -# sage: from mjo.eja.eja_algebra import (random_eja, -# ....: HadamardEJA, -# ....: RealSymmetricEJA, -# ....: DirectSumEJA) - -# EXAMPLES:: - -# sage: J1 = HadamardEJA(2) -# sage: J2 = RealSymmetricEJA(3) -# sage: J = DirectSumEJA(J1,J2) -# sage: J.dimension() -# 8 -# sage: J.rank() -# 5 - -# TESTS: - -# The external direct sum construction is only valid when the two factors -# have the same base ring; an error is raised otherwise:: - -# sage: set_random_seed() -# sage: J1 = random_eja(field=AA) -# sage: J2 = random_eja(field=QQ,orthonormalize=False) -# sage: J = DirectSumEJA(J1,J2) -# Traceback (most recent call last): -# ... -# ValueError: algebras must share the same base field - -# """ -# def __init__(self, J1, J2, **kwargs): -# if J1.base_ring() != J2.base_ring(): -# raise ValueError("algebras must share the same base field") -# field = J1.base_ring() - -# self._factors = (J1, J2) -# n1 = J1.dimension() -# n2 = J2.dimension() -# n = n1+n2 -# V = VectorSpace(field, n) -# mult_table = [ [ V.zero() for j in range(i+1) ] -# for i in range(n) ] -# for i in range(n1): -# for j in range(i+1): -# p = (J1.monomial(i)*J1.monomial(j)).to_vector() -# mult_table[i][j] = V(p.list() + [field.zero()]*n2) - -# for i in range(n2): -# for j in range(i+1): -# p = (J2.monomial(i)*J2.monomial(j)).to_vector() -# mult_table[n1+i][n1+j] = V([field.zero()]*n1 + p.list()) - -# # TODO: build the IP table here from the two constituent IP -# # matrices (it'll be block diagonal, I think). -# ip_table = [ [ field.zero() for j in range(i+1) ] -# for i in range(n) ] -# super(DirectSumEJA, self).__init__(field, -# mult_table, -# ip_table, -# check_axioms=False, -# **kwargs) -# self.rank.set_cache(J1.rank() + J2.rank()) - - -# def factors(self): -# r""" -# Return the pair of this algebra's factors. -# SETUP:: +class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct, + FiniteDimensionalEJA): + r""" + The external (orthogonal) direct sum of two or more Euclidean + Jordan algebras. Every Euclidean Jordan algebra decomposes into an + orthogonal direct sum of simple Euclidean Jordan algebras which is + then isometric to a Cartesian product, so no generality is lost by + providing only this construction. -# sage: from mjo.eja.eja_algebra import (HadamardEJA, -# ....: JordanSpinEJA, -# ....: DirectSumEJA) + SETUP:: -# EXAMPLES:: + sage: from mjo.eja.eja_algebra import (CartesianProductEJA, + ....: HadamardEJA, + ....: JordanSpinEJA, + ....: RealSymmetricEJA) -# sage: J1 = HadamardEJA(2, field=QQ) -# sage: J2 = JordanSpinEJA(3, field=QQ) -# sage: J = DirectSumEJA(J1,J2) -# sage: J.factors() -# (Euclidean Jordan algebra of dimension 2 over Rational Field, -# Euclidean Jordan algebra of dimension 3 over Rational Field) + EXAMPLES: + + The Jordan product is inherited from our factors and implemented by + our CombinatorialFreeModule Cartesian product superclass:: + + sage: set_random_seed() + sage: J1 = HadamardEJA(2) + sage: J2 = RealSymmetricEJA(2) + sage: J = cartesian_product([J1,J2]) + sage: x,y = J.random_elements(2) + sage: x*y in J + True + + The ability to retrieve the original factors is implemented by our + CombinatorialFreeModule Cartesian product superclass:: + + sage: J1 = HadamardEJA(2, field=QQ) + sage: J2 = JordanSpinEJA(3, field=QQ) + sage: J = cartesian_product([J1,J2]) + sage: J.cartesian_factors() + (Euclidean Jordan algebra of dimension 2 over Rational Field, + Euclidean Jordan algebra of dimension 3 over Rational Field) + + TESTS: + + All factors must share the same base field:: + + sage: J1 = HadamardEJA(2, field=QQ) + sage: J2 = RealSymmetricEJA(2) + sage: CartesianProductEJA((J1,J2)) + Traceback (most recent call last): + ... + ValueError: all factors must share the same base field + """ + def __init__(self, modules, **kwargs): + CombinatorialFreeModule_CartesianProduct.__init__(self, modules, **kwargs) + field = modules[0].base_ring() + if not all( J.base_ring() == field for J in modules ): + raise ValueError("all factors must share the same base field") + + M = cartesian_product( [J.matrix_space() for J in modules] ) + + m = len(modules) + W = VectorSpace(field,m) + self._matrix_basis = [] + for k in range(m): + for a in modules[k].matrix_basis(): + v = W.zero().list() + v[k] = a + self._matrix_basis.append(M(v)) + + self._matrix_basis = tuple(self._matrix_basis) + + n = len(self._matrix_basis) + # TODO: + # + # Initialize the FDEJA class, too. Does this override the + # initialization that we did for the + # CombinatorialFreeModule_CartesianProduct class? If not, we + # will probably have to duplicate some of the work (i.e. one + # of the constructors). Since the CartesianProduct one is + # smaller, that makes the most sense to copy/paste if it comes + # down to that. + # + + self.rank.set_cache(sum(J.rank() for J in modules)) + + @cached_method + def cartesian_projection(self, i): + r""" + SETUP:: + + sage: from mjo.eja.eja_algebra import (random_eja, + ....: HadamardEJA, + ....: RealSymmetricEJA) + + EXAMPLES:: + + sage: J1 = HadamardEJA(2) + sage: J2 = RealSymmetricEJA(2) + sage: J = cartesian_product([J1,J2]) + sage: J.cartesian_projection(0) + Linear operator between finite-dimensional Euclidean Jordan + algebras represented by the matrix: + [1 0 0 0 0] + [0 1 0 0 0] + Domain: Euclidean Jordan algebra of dimension 2 over Algebraic + Real Field (+) Euclidean Jordan algebra of dimension 3 over + Algebraic Real Field + Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic + Real Field + sage: J.cartesian_projection(1) + Linear operator between finite-dimensional Euclidean Jordan + algebras represented by the matrix: + [0 0 1 0 0] + [0 0 0 1 0] + [0 0 0 0 1] + Domain: Euclidean Jordan algebra of dimension 2 over Algebraic + Real Field (+) Euclidean Jordan algebra of dimension 3 over + Algebraic Real Field + Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic + Real Field + + TESTS: + + The answer never changes:: + + sage: set_random_seed() + sage: J1 = random_eja() + sage: J2 = random_eja() + sage: J = cartesian_product([J1,J2]) + sage: P0 = J.cartesian_projection(0) + sage: P1 = J.cartesian_projection(0) + sage: P0 == P1 + True + + """ + Ji = self.cartesian_factors()[i] + # We reimplement the CombinatorialFreeModule superclass method + # because if we don't, something gets messed up with the caching + # and the answer changes the second time you run it. See the TESTS. + Pi = self._module_morphism(lambda j_t: Ji.monomial(j_t[1]) + if i == j_t[0] else Ji.zero(), + codomain=Ji) + return FiniteDimensionalEJAOperator(self,Ji,Pi.matrix()) + + @cached_method + def cartesian_embedding(self, i): + r""" + SETUP:: + + sage: from mjo.eja.eja_algebra import (random_eja, + ....: HadamardEJA, + ....: RealSymmetricEJA) + + EXAMPLES:: + + sage: J1 = HadamardEJA(2) + sage: J2 = RealSymmetricEJA(2) + sage: J = cartesian_product([J1,J2]) + sage: J.cartesian_embedding(0) + Linear operator between finite-dimensional Euclidean Jordan + algebras represented by the matrix: + [1 0] + [0 1] + [0 0] + [0 0] + [0 0] + Domain: Euclidean Jordan algebra of dimension 2 over + Algebraic Real Field + Codomain: Euclidean Jordan algebra of dimension 2 over + Algebraic Real Field (+) Euclidean Jordan algebra of + dimension 3 over Algebraic Real Field + sage: J.cartesian_embedding(1) + Linear operator between finite-dimensional Euclidean Jordan + algebras represented by the matrix: + [0 0 0] + [0 0 0] + [1 0 0] + [0 1 0] + [0 0 1] + Domain: Euclidean Jordan algebra of dimension 3 over + Algebraic Real Field + Codomain: Euclidean Jordan algebra of dimension 2 over + Algebraic Real Field (+) Euclidean Jordan algebra of + dimension 3 over Algebraic Real Field + + TESTS: + + The answer never changes:: + + sage: set_random_seed() + sage: J1 = random_eja() + sage: J2 = random_eja() + sage: J = cartesian_product([J1,J2]) + sage: E0 = J.cartesian_embedding(0) + sage: E1 = J.cartesian_embedding(0) + sage: E0 == E1 + True + + """ + Ji = self.cartesian_factors()[i] + # We reimplement the CombinatorialFreeModule superclass method + # because if we don't, something gets messed up with the caching + # and the answer changes the second time you run it. See the TESTS. + Ei = Ji._module_morphism(lambda t: self.monomial((i, t)), codomain=self) + return FiniteDimensionalEJAOperator(Ji,self,Ei.matrix()) + + +FiniteDimensionalEJA.CartesianProduct = CartesianProductEJA -# """ -# return self._factors # def projections(self): # r"""