X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_element.py;h=e1ea609efce334a9800ce162ba74d1342662da41;hb=31f14e8f8c51d34823ca26aaa9568f924304d08b;hp=b5f661bf037667d4d9210b6af65a3c1f95120732;hpb=958623afecbcfa2f0a7daa0bdaab458814e8e4a6;p=sage.d.git diff --git a/mjo/eja/eja_element.py b/mjo/eja/eja_element.py index b5f661b..e1ea609 100644 --- a/mjo/eja/eja_element.py +++ b/mjo/eja/eja_element.py @@ -4,7 +4,8 @@ from sage.modules.free_module import VectorSpace from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement from mjo.eja.eja_operator import FiniteDimensionalEJAOperator -from mjo.eja.eja_utils import _mat2vec +from mjo.eja.eja_utils import _scale + class FiniteDimensionalEJAElement(IndexedFreeModuleElement): """ @@ -42,14 +43,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The definition of `x^2` is the unambiguous `x*x`:: - sage: set_random_seed() sage: x = random_eja().random_element() sage: x*x == (x^2) True A few examples of power-associativity:: - sage: set_random_seed() sage: x = random_eja().random_element() sage: x*(x*x)*(x*x) == x^5 True @@ -59,7 +58,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): We also know that powers operator-commute (Koecher, Chapter III, Corollary 1):: - sage: set_random_seed() sage: x = random_eja().random_element() sage: m = ZZ.random_element(0,10) sage: n = ZZ.random_element(0,10) @@ -106,7 +104,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): We should always get back an element of the algebra:: - sage: set_random_seed() sage: p = PolynomialRing(AA, 't').random_element() sage: J = random_eja() sage: x = J.random_element() @@ -131,7 +128,8 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): SETUP:: - sage: from mjo.eja.eja_algebra import HadamardEJA + sage: from mjo.eja.eja_algebra import (random_eja, + ....: HadamardEJA) EXAMPLES: @@ -155,11 +153,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The characteristic polynomial of an element should evaluate to zero on that element:: - sage: set_random_seed() - sage: x = HadamardEJA(3).random_element() + sage: x = random_eja().random_element() sage: p = x.characteristic_polynomial() - sage: x.apply_univariate_polynomial(p) - 0 + sage: x.apply_univariate_polynomial(p).is_zero() + True The characteristic polynomials of the zero and unit elements should be what we think they are in a subalgebra, too:: @@ -167,8 +164,8 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: J = HadamardEJA(3) sage: p1 = J.one().characteristic_polynomial() sage: q1 = J.zero().characteristic_polynomial() - sage: e0,e1,e2 = J.gens() - sage: A = (e0 + 2*e1 + 3*e2).subalgebra_generated_by() # dim 3 + sage: b0,b1,b2 = J.gens() + sage: A = (b0 + 2*b1 + 3*b2).subalgebra_generated_by() # dim 3 sage: p2 = A.one().characteristic_polynomial() sage: q2 = A.zero().characteristic_polynomial() sage: p1 == p2 @@ -237,7 +234,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Ensure that we can always compute an inner product, and that it gives us back a real number:: - sage: set_random_seed() sage: J = random_eja() sage: x,y = J.random_elements(2) sage: x.inner_product(y) in RLF @@ -265,7 +261,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The definition of a Jordan algebra says that any element operator-commutes with its square:: - sage: set_random_seed() sage: x = random_eja().random_element() sage: x.operator_commutes_with(x^2) True @@ -274,7 +269,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Test Lemma 1 from Chapter III of Koecher:: - sage: set_random_seed() sage: u,v = random_eja().random_elements(2) sage: lhs = u.operator_commutes_with(u*v) sage: rhs = v.operator_commutes_with(u^2) @@ -284,7 +278,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Test the first polarization identity from my notes, Koecher Chapter III, or from Baes (2.3):: - sage: set_random_seed() sage: x,y = random_eja().random_elements(2) sage: Lx = x.operator() sage: Ly = y.operator() @@ -296,32 +289,32 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Test the second polarization identity from my notes or from Baes (2.4):: - sage: set_random_seed() - sage: x,y,z = random_eja().random_elements(3) - sage: Lx = x.operator() - sage: Ly = y.operator() - sage: Lz = z.operator() - sage: Lzy = (z*y).operator() - sage: Lxy = (x*y).operator() - sage: Lxz = (x*z).operator() - sage: bool(Lx*Lzy + Lz*Lxy + Ly*Lxz == Lzy*Lx + Lxy*Lz + Lxz*Ly) + sage: x,y,z = random_eja().random_elements(3) # long time + sage: Lx = x.operator() # long time + sage: Ly = y.operator() # long time + sage: Lz = z.operator() # long time + sage: Lzy = (z*y).operator() # long time + sage: Lxy = (x*y).operator() # long time + sage: Lxz = (x*z).operator() # long time + sage: lhs = Lx*Lzy + Lz*Lxy + Ly*Lxz # long time + sage: rhs = Lzy*Lx + Lxy*Lz + Lxz*Ly # long time + sage: bool(lhs == rhs) # long time True Test the third polarization identity from my notes or from Baes (2.5):: - sage: set_random_seed() - sage: u,y,z = random_eja().random_elements(3) - sage: Lu = u.operator() - sage: Ly = y.operator() - sage: Lz = z.operator() - sage: Lzy = (z*y).operator() - sage: Luy = (u*y).operator() - sage: Luz = (u*z).operator() - sage: Luyz = (u*(y*z)).operator() - sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz - sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly - sage: bool(lhs == rhs) + sage: u,y,z = random_eja().random_elements(3) # long time + sage: Lu = u.operator() # long time + sage: Ly = y.operator() # long time + sage: Lz = z.operator() # long time + sage: Lzy = (z*y).operator() # long time + sage: Luy = (u*y).operator() # long time + sage: Luz = (u*z).operator() # long time + sage: Luyz = (u*(y*z)).operator() # long time + sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz # long time + sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly # long time + sage: bool(lhs == rhs) # long time True """ @@ -348,7 +341,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): EXAMPLES:: sage: J = JordanSpinEJA(2) - sage: e0,e1 = J.gens() sage: x = sum( J.gens() ) sage: x.det() 0 @@ -356,7 +348,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): :: sage: J = JordanSpinEJA(3) - sage: e0,e1,e2 = J.gens() sage: x = sum( J.gens() ) sage: x.det() -1 @@ -377,7 +368,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): An element is invertible if and only if its determinant is non-zero:: - sage: set_random_seed() sage: x = random_eja().random_element() sage: x.is_invertible() == (x.det() != 0) True @@ -385,15 +375,14 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Ensure that the determinant is multiplicative on an associative subalgebra as in Faraut and Korányi's Proposition II.2.2:: - sage: set_random_seed() - sage: J = random_eja().random_element().subalgebra_generated_by() + sage: x0 = random_eja().random_element() + sage: J = x0.subalgebra_generated_by(orthonormalize=False) sage: x,y = J.random_elements(2) sage: (x*y).det() == x.det()*y.det() True - The determinant in matrix algebras is just the usual determinant:: + The determinant in real matrix algebras is the usual determinant:: - sage: set_random_seed() sage: X = matrix.random(QQ,3) sage: X = X + X.T sage: J1 = RealSymmetricEJA(3) @@ -406,21 +395,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: actual2 == expected True - :: - - sage: set_random_seed() - sage: J1 = ComplexHermitianEJA(2) - sage: J2 = ComplexHermitianEJA(2,field=QQ,orthonormalize=False) - sage: X = matrix.random(GaussianIntegers(), 2) - sage: X = X + X.H - sage: expected = AA(X.det()) - sage: actual1 = J1(J1.real_embed(X)).det() - sage: actual2 = J2(J2.real_embed(X)).det() - sage: expected == actual1 - True - sage: expected == actual2 - True - """ P = self.parent() r = P.rank() @@ -465,7 +439,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The inverse in the spin factor algebra is given in Alizadeh's Example 11.11:: - sage: set_random_seed() sage: J = JordanSpinEJA.random_instance() sage: x = J.random_element() sage: while not x.is_invertible(): @@ -490,14 +463,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The identity element is its own inverse:: - sage: set_random_seed() sage: J = random_eja() sage: J.one().inverse() == J.one() True If an element has an inverse, it acts like one:: - sage: set_random_seed() sage: J = random_eja() sage: x = J.random_element() sage: (not x.is_invertible()) or (x.inverse()*x == J.one()) @@ -505,7 +476,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The inverse of the inverse is what we started with:: - sage: set_random_seed() sage: J = random_eja() sage: x = J.random_element() sage: (not x.is_invertible()) or (x.inverse().inverse() == x) @@ -515,17 +485,17 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): of an element is the inverse of its left-multiplication operator applied to the algebra's identity, when that inverse exists:: - sage: set_random_seed() - sage: J = random_eja() - sage: x = J.random_element() - sage: (not x.operator().is_invertible()) or ( - ....: x.operator().inverse()(J.one()) == x.inverse() ) + sage: J = random_eja() # long time + sage: x = J.random_element() # long time + sage: (not x.operator().is_invertible()) or ( # long time + ....: x.operator().inverse()(J.one()) # long time + ....: == # long time + ....: x.inverse() ) # long time True Check that the fast (cached) and slow algorithms give the same answer:: - sage: set_random_seed() # long time sage: J = random_eja(field=QQ, orthonormalize=False) # long time sage: x = J.random_element() # long time sage: while not x.is_invertible(): # long time @@ -537,15 +507,18 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): True """ not_invertible_msg = "element is not invertible" - if self.parent()._charpoly_coefficients.is_in_cache(): + + algebra = self.parent() + if algebra._charpoly_coefficients.is_in_cache(): # We can invert using our charpoly if it will be fast to # compute. If the coefficients are cached, our rank had # better be too! if self.det().is_zero(): raise ZeroDivisionError(not_invertible_msg) - r = self.parent().rank() + r = algebra.rank() a = self.characteristic_polynomial().coefficients(sparse=False) - return (-1)**(r+1)*sum(a[i+1]*self**i for i in range(r))/self.det() + return (-1)**(r+1)*algebra.sum(a[i+1]*self**i + for i in range(r))/self.det() try: inv = (~self.quadratic_representation())(self) @@ -577,14 +550,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The identity element is always invertible:: - sage: set_random_seed() sage: J = random_eja() sage: J.one().is_invertible() True The zero element is never invertible in a non-trivial algebra:: - sage: set_random_seed() sage: J = random_eja() sage: (not J.is_trivial()) and J.zero().is_invertible() False @@ -592,7 +563,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Test that the fast (cached) and slow algorithms give the same answer:: - sage: set_random_seed() # long time sage: J = random_eja(field=QQ, orthonormalize=False) # long time sage: x = J.random_element() # long time sage: slow = x.is_invertible() # long time @@ -664,7 +634,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): element should always be in terms of minimal idempotents:: sage: J = JordanSpinEJA(4) - sage: x = sum( i*J.gens()[i] for i in range(len(J.gens())) ) + sage: x = sum( i*J.monomial(i) for i in range(len(J.gens())) ) sage: x.is_regular() True sage: [ c.is_primitive_idempotent() @@ -675,14 +645,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The identity element is minimal only in an EJA of rank one:: - sage: set_random_seed() sage: J = random_eja() sage: J.rank() == 1 or not J.one().is_primitive_idempotent() True A non-idempotent cannot be a minimal idempotent:: - sage: set_random_seed() sage: J = JordanSpinEJA(4) sage: x = J.random_element() sage: (not x.is_idempotent()) and x.is_primitive_idempotent() @@ -692,7 +660,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): idempotent if and only if it's idempotent with trace equal to unity:: - sage: set_random_seed() sage: J = JordanSpinEJA(4) sage: x = J.random_element() sage: expected = (x.is_idempotent() and x.trace() == 1) @@ -702,7 +669,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Primitive idempotents must be non-zero:: - sage: set_random_seed() sage: J = random_eja() sage: J.zero().is_idempotent() True @@ -759,14 +725,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The identity element is never nilpotent, except in a trivial EJA:: - sage: set_random_seed() sage: J = random_eja() sage: J.one().is_nilpotent() and not J.is_trivial() False The additive identity is always nilpotent:: - sage: set_random_seed() sage: random_eja().zero().is_nilpotent() True @@ -793,7 +757,9 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: J = JordanSpinEJA(5) sage: J.one().is_regular() False - sage: e0, e1, e2, e3, e4 = J.gens() # e0 is the identity + sage: b0, b1, b2, b3, b4 = J.gens() + sage: b0 == J.one() + True sage: for x in J.gens(): ....: (J.one() + x).is_regular() False @@ -807,7 +773,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The zero element should never be regular, unless the parent algebra has dimension less than or equal to one:: - sage: set_random_seed() sage: J = random_eja() sage: J.dimension() <= 1 or not J.zero().is_regular() True @@ -815,7 +780,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The unit element isn't regular unless the algebra happens to consist of only its scalar multiples:: - sage: set_random_seed() sage: J = random_eja() sage: J.dimension() <= 1 or not J.one().is_regular() True @@ -843,14 +807,13 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: J = JordanSpinEJA(4) sage: J.one().degree() 1 - sage: e0,e1,e2,e3 = J.gens() - sage: (e0 - e1).degree() + sage: b0,b1,b2,b3 = J.gens() + sage: (b0 - b1).degree() 2 In the spin factor algebra (of rank two), all elements that aren't multiples of the identity are regular:: - sage: set_random_seed() sage: J = JordanSpinEJA.random_instance() sage: n = J.dimension() sage: x = J.random_element() @@ -862,7 +825,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The zero and unit elements are both of degree one in nontrivial algebras:: - sage: set_random_seed() sage: J = random_eja() sage: d = J.zero().degree() sage: (J.is_trivial() and d == 0) or d == 1 @@ -873,7 +835,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Our implementation agrees with the definition:: - sage: set_random_seed() sage: x = random_eja().random_element() sage: x.degree() == x.minimal_polynomial().degree() True @@ -982,7 +943,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): always the same, except in trivial algebras where the minimal polynomial of the unit/zero element is ``1``:: - sage: set_random_seed() sage: J = random_eja() sage: mu = J.one().minimal_polynomial() sage: t = mu.parent().gen() @@ -996,7 +956,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The degree of an element is (by one definition) the degree of its minimal polynomial:: - sage: set_random_seed() sage: x = random_eja().random_element() sage: x.degree() == x.minimal_polynomial().degree() True @@ -1007,9 +966,8 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): identity. We require the dimension of the algebra to be at least two here so that said elements actually exist:: - sage: set_random_seed() - sage: n_max = max(2, JordanSpinEJA._max_random_instance_size()) - sage: n = ZZ.random_element(2, n_max) + sage: d_max = JordanSpinEJA._max_random_instance_dimension() + sage: n = ZZ.random_element(2, max(2,d_max)) sage: J = JordanSpinEJA(n) sage: y = J.random_element() sage: while y == y.coefficient(0)*J.one(): @@ -1024,18 +982,17 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The minimal polynomial should always kill its element:: - sage: set_random_seed() - sage: x = random_eja().random_element() - sage: p = x.minimal_polynomial() - sage: x.apply_univariate_polynomial(p) + sage: x = random_eja().random_element() # long time + sage: p = x.minimal_polynomial() # long time + sage: x.apply_univariate_polynomial(p) # long time 0 The minimal polynomial is invariant under a change of basis, and in particular, a re-scaling of the basis:: - sage: set_random_seed() - sage: n_max = RealSymmetricEJA._max_random_instance_size() - sage: n = ZZ.random_element(1, n_max) + sage: d_max = RealSymmetricEJA._max_random_instance_dimension() + sage: d = ZZ.random_element(1, d_max) + sage: n = RealSymmetricEJA._max_random_instance_size(d) sage: J1 = RealSymmetricEJA(n) sage: J2 = RealSymmetricEJA(n,orthonormalize=False) sage: X = random_matrix(AA,n) @@ -1047,19 +1004,30 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): """ if self.is_zero(): - # We would generate a zero-dimensional subalgebra - # where the minimal polynomial would be constant. - # That might be correct, but only if *this* algebra - # is trivial too. - if not self.parent().is_trivial(): - # Pretty sure we know what the minimal polynomial of - # the zero operator is going to be. This ensures - # consistency of e.g. the polynomial variable returned - # in the "normal" case without us having to think about it. - return self.operator().minimal_polynomial() - + # Pretty sure we know what the minimal polynomial of + # the zero operator is going to be. This ensures + # consistency of e.g. the polynomial variable returned + # in the "normal" case without us having to think about it. + return self.operator().minimal_polynomial() + + # If we don't orthonormalize the subalgebra's basis, then the + # first two monomials in the subalgebra will be self^0 and + # self^1... assuming that self^1 is not a scalar multiple of + # self^0 (the unit element). We special case these to avoid + # having to solve a system to coerce self into the subalgebra. A = self.subalgebra_generated_by(orthonormalize=False) - return A(self).operator().minimal_polynomial() + + if A.dimension() == 1: + # Does a solve to find the scalar multiple alpha such that + # alpha*unit = self. We have to do this because the basis + # for the subalgebra will be [ self^0 ], and not [ self^1 ]! + unit = self.parent().one() + alpha = self.to_vector() / unit.to_vector() + return (unit.operator()*alpha).minimal_polynomial() + else: + # If the dimension of the subalgebra is >= 2, then we just + # use the second basis element. + return A.monomial(1).operator().minimal_polynomial() @@ -1085,29 +1053,27 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: J = ComplexHermitianEJA(3) sage: J.one() - e0 + e3 + e8 + b0 + b3 + b8 sage: J.one().to_matrix() - [1 0 0 0 0 0] - [0 1 0 0 0 0] - [0 0 1 0 0 0] - [0 0 0 1 0 0] - [0 0 0 0 1 0] - [0 0 0 0 0 1] + +---+---+---+ + | 1 | 0 | 0 | + +---+---+---+ + | 0 | 1 | 0 | + +---+---+---+ + | 0 | 0 | 1 | + +---+---+---+ :: sage: J = QuaternionHermitianEJA(2) sage: J.one() - e0 + e5 + b0 + b5 sage: J.one().to_matrix() - [1 0 0 0 0 0 0 0] - [0 1 0 0 0 0 0 0] - [0 0 1 0 0 0 0 0] - [0 0 0 1 0 0 0 0] - [0 0 0 0 1 0 0 0] - [0 0 0 0 0 1 0 0] - [0 0 0 0 0 0 1 0] - [0 0 0 0 0 0 0 1] + +---+---+ + | 1 | 0 | + +---+---+ + | 0 | 1 | + +---+---+ This also works in Cartesian product algebras:: @@ -1125,19 +1091,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): B = self.parent().matrix_basis() W = self.parent().matrix_space() - if self.parent()._matrix_basis_is_cartesian: - # Aaaaand linear combinations don't work in Cartesian - # product spaces, even though they provide a method - # with that name. - pairs = zip(B, self.to_vector()) - return sum( ( W(tuple(alpha*b_i for b_i in b)) - for (b,alpha) in pairs ), - W.zero()) - else: - # This is just a manual "from_vector()", but of course - # matrix spaces aren't vector spaces in sage, so they - # don't have a from_vector() method. - return W.linear_combination( zip(B, self.to_vector()) ) + # This is just a manual "from_vector()", but of course + # matrix spaces aren't vector spaces in sage, so they + # don't have a from_vector() method. + return W.linear_combination( zip(B, self.to_vector()) ) @@ -1179,7 +1136,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): TESTS:: - sage: set_random_seed() sage: J = random_eja() sage: x,y = J.random_elements(2) sage: x.operator()(y) == x*y @@ -1208,7 +1164,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The explicit form in the spin factor algebra is given by Alizadeh's Example 11.12:: - sage: set_random_seed() sage: x = JordanSpinEJA.random_instance().random_element() sage: x_vec = x.to_vector() sage: Q = matrix.identity(x.base_ring(), 0) @@ -1228,7 +1183,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Test all of the properties from Theorem 11.2 in Alizadeh:: - sage: set_random_seed() sage: J = random_eja() sage: x,y = J.random_elements(2) sage: Lx = x.operator() @@ -1345,11 +1299,11 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: J = RealSymmetricEJA(3) sage: J.one() - e0 + e2 + e5 + b0 + b2 + b5 sage: J.one().spectral_decomposition() - [(1, e0 + e2 + e5)] + [(1, b0 + b2 + b5)] sage: J.zero().spectral_decomposition() - [(0, e0 + e2 + e5)] + [(0, b0 + b2 + b5)] TESTS:: @@ -1374,13 +1328,13 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The spectral decomposition should work in subalgebras, too:: sage: J = RealSymmetricEJA(4) - sage: (e0, e1, e2, e3, e4, e5, e6, e7, e8, e9) = J.gens() - sage: A = 2*e5 - 2*e8 + sage: (b0, b1, b2, b3, b4, b5, b6, b7, b8, b9) = J.gens() + sage: A = 2*b5 - 2*b8 sage: (lambda1, c1) = A.spectral_decomposition()[1] sage: (J0, J5, J1) = J.peirce_decomposition(c1) sage: (f0, f1, f2) = J1.gens() sage: f0.spectral_decomposition() - [(0, f2), (1, f0)] + [(0, 1.000000000000000?*c2), (1, 1.000000000000000?*c0)] """ A = self.subalgebra_generated_by(orthonormalize=True) @@ -1422,9 +1376,8 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): This subalgebra, being composed of only powers, is associative:: - sage: set_random_seed() sage: x0 = random_eja().random_element() - sage: A = x0.subalgebra_generated_by() + sage: A = x0.subalgebra_generated_by(orthonormalize=False) sage: x,y,z = A.random_elements(3) sage: (x*y)*z == x*(y*z) True @@ -1432,9 +1385,8 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): Squaring in the subalgebra should work the same as in the superalgebra:: - sage: set_random_seed() sage: x = random_eja().random_element() - sage: A = x.subalgebra_generated_by() + sage: A = x.subalgebra_generated_by(orthonormalize=False) sage: A(x^2) == A(x)*A(x) True @@ -1443,7 +1395,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): element... unless the original algebra was trivial, in which case the subalgebra is trivial too:: - sage: set_random_seed() sage: A = random_eja().zero().subalgebra_generated_by() sage: (A.is_trivial() and A.dimension() == 0) or A.dimension() == 1 True @@ -1474,8 +1425,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): where there are non-nilpotent elements, or that we get the dumb solution in the trivial algebra:: - sage: set_random_seed() - sage: J = random_eja() + sage: J = random_eja(field=QQ, orthonormalize=False) sage: x = J.random_element() sage: while x.is_nilpotent() and not J.is_trivial(): ....: x = J.random_element() @@ -1558,20 +1508,24 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The trace of an element is a real number:: - sage: set_random_seed() sage: J = random_eja() sage: J.random_element().trace() in RLF True The trace is linear:: - sage: set_random_seed() sage: J = random_eja() sage: x,y = J.random_elements(2) sage: alpha = J.base_ring().random_element() sage: (alpha*x + y).trace() == alpha*x.trace() + y.trace() True + The trace of a square is nonnegative:: + + sage: x = random_eja().random_element() + sage: (x*x).trace() >= 0 + True + """ P = self.parent() r = P.rank() @@ -1588,6 +1542,102 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): # we want the negative of THAT for the trace. return -p(*self.to_vector()) + def operator_inner_product(self, other): + r""" + Return the operator inner product of myself and ``other``. + + The "operator inner product," whose name is not standard, is + defined be the usual linear-algebraic trace of the + ``(x*y).operator()``. + + Proposition III.1.5 in Faraut and Korányi shows that on any + Euclidean Jordan algebra, this is another associative inner + product under which the cone of squares is symmetric. + + This *probably* works even if the basis hasn't been + orthonormalized because the eigenvalues of the corresponding + matrix don't change when the basis does (they're preserved by + any similarity transformation). + + SETUP:: + + sage: from mjo.eja.eja_algebra import (JordanSpinEJA, + ....: RealSymmetricEJA, + ....: ComplexHermitianEJA, + ....: random_eja) + + EXAMPLES: + + Proposition III.4.2 of Faraut and Korányi shows that on a + simple algebra of rank `r` and dimension `n`, this inner + product is `n/r` times the canonical + :meth:`trace_inner_product`:: + + sage: J = JordanSpinEJA(4, field=QQ) + sage: x,y = J.random_elements(2) + sage: n = J.dimension() + sage: r = J.rank() + sage: actual = x.operator_inner_product(y) + sage: expected = (n/r)*x.trace_inner_product(y) + sage: actual == expected + True + + :: + + sage: J = RealSymmetricEJA(3) + sage: x,y = J.random_elements(2) + sage: n = J.dimension() + sage: r = J.rank() + sage: actual = x.operator_inner_product(y) + sage: expected = (n/r)*x.trace_inner_product(y) + sage: actual == expected + True + + :: + + sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False) + sage: x,y = J.random_elements(2) + sage: n = J.dimension() + sage: r = J.rank() + sage: actual = x.operator_inner_product(y) + sage: expected = (n/r)*x.trace_inner_product(y) + sage: actual == expected + True + + TESTS: + + The operator inner product is commutative, bilinear, and + associative:: + + sage: J = random_eja() + sage: x,y,z = J.random_elements(3) + sage: # commutative + sage: x.operator_inner_product(y) == y.operator_inner_product(x) + True + sage: # bilinear + sage: a = J.base_ring().random_element() + sage: actual = (a*(x+z)).operator_inner_product(y) + sage: expected = ( a*x.operator_inner_product(y) + + ....: a*z.operator_inner_product(y) ) + sage: actual == expected + True + sage: actual = x.operator_inner_product(a*(y+z)) + sage: expected = ( a*x.operator_inner_product(y) + + ....: a*x.operator_inner_product(z) ) + sage: actual == expected + True + sage: # associative + sage: actual = (x*y).operator_inner_product(z) + sage: expected = y.operator_inner_product(x*z) + sage: actual == expected + True + + """ + if not other in self.parent(): + raise TypeError("'other' must live in the same algebra") + + return (self*other).operator().matrix().trace() + def trace_inner_product(self, other): """ @@ -1601,14 +1651,13 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): The trace inner product is commutative, bilinear, and associative:: - sage: set_random_seed() sage: J = random_eja() sage: x,y,z = J.random_elements(3) sage: # commutative sage: x.trace_inner_product(y) == y.trace_inner_product(x) True sage: # bilinear - sage: a = J.base_ring().random_element(); + sage: a = J.base_ring().random_element() sage: actual = (a*(x+z)).trace_inner_product(y) sage: expected = ( a*x.trace_inner_product(y) + ....: a*z.trace_inner_product(y) ) @@ -1655,3 +1704,29 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): """ return self.trace_inner_product(self).sqrt() + + +class CartesianProductEJAElement(FiniteDimensionalEJAElement): + def det(self): + r""" + Compute the determinant of this product-element using the + determianants of its factors. + + This result Follows from the spectral decomposition of (say) + the pair `(x,y)` in terms of the Jordan frame `\left\{ (c_1, + 0),(c_2, 0),...,(0,d_1),(0,d_2),... \right\}. + """ + from sage.misc.misc_c import prod + return prod( f.det() for f in self.cartesian_factors() ) + + def to_matrix(self): + # An override is necessary to call our custom _scale(). + B = self.parent().matrix_basis() + W = self.parent().matrix_space() + + # Aaaaand linear combinations don't work in Cartesian + # product spaces, even though they provide a method with + # that name. This is hidden behind an "if" because the + # _scale() function is slow. + pairs = zip(B, self.to_vector()) + return W.sum( _scale(b, alpha) for (b,alpha) in pairs )