X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_element.py;h=f4d5995ce9302d2fdc518dac2b0ca20f0fadf6d1;hb=HEAD;hp=e1ea609efce334a9800ce162ba74d1342662da41;hpb=31f14e8f8c51d34823ca26aaa9568f924304d08b;p=sage.d.git diff --git a/mjo/eja/eja_element.py b/mjo/eja/eja_element.py index e1ea609..f4d5995 100644 --- a/mjo/eja/eja_element.py +++ b/mjo/eja/eja_element.py @@ -3,11 +3,11 @@ from sage.misc.cachefunc import cached_method 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_operator import EJAOperator from mjo.eja.eja_utils import _scale -class FiniteDimensionalEJAElement(IndexedFreeModuleElement): +class EJAElement(IndexedFreeModuleElement): """ An element of a Euclidean Jordan algebra. """ @@ -332,7 +332,8 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): SETUP:: - sage: from mjo.eja.eja_algebra import (JordanSpinEJA, + sage: from mjo.eja.eja_algebra import (AlbertEJA, + ....: JordanSpinEJA, ....: TrivialEJA, ....: RealSymmetricEJA, ....: ComplexHermitianEJA, @@ -395,6 +396,22 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: actual2 == expected True + There's a formula for the determinant of the Albert algebra + (Yokota, Section 2.1):: + + sage: def albert_det(x): + ....: X = x.to_matrix() + ....: res = X[0,0]*X[1,1]*X[2,2] + ....: res += 2*(X[1,2]*X[2,0]*X[0,1]).real() + ....: res -= X[0,0]*X[1,2]*X[2,1] + ....: res -= X[1,1]*X[2,0]*X[0,2] + ....: res -= X[2,2]*X[0,1]*X[1,0] + ....: return res.leading_coefficient() + sage: J = AlbertEJA(field=QQ, orthonormalize=False) + sage: xs = J.random_elements(10) + sage: all( albert_det(x) == x.det() for x in xs ) + True + """ P = self.parent() r = P.rank() @@ -538,7 +555,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): If computing my determinant will be fast, we do so and compare with zero (Proposition II.2.4 in Faraut and - Koranyi). Otherwise, Proposition II.3.2 in Faraut and Koranyi + Korányi). Otherwise, Proposition II.3.2 in Faraut and Korányi reduces the problem to the invertibility of my quadratic representation. @@ -795,7 +812,23 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): ALGORITHM: - ......... + First we handle the special cases where the algebra is + trivial, this element is zero, or the dimension of the algebra + is one and this element is not zero. With those out of the + way, we may assume that ``self`` is nonzero, the algebra is + nontrivial, and that the dimension of the algebra is at least + two. + + Beginning with the algebra's unit element (power zero), we add + successive (basis representations of) powers of this element + to a matrix, row-reducing at each step. After row-reducing, we + check the rank of the matrix. If adding a row and row-reducing + does not increase the rank of the matrix at any point, the row + we've just added lives in the span of the previous ones; thus + the corresponding power of ``self`` lives in the span of its + lesser powers. When that happens, the degree of the minimal + polynomial is the rank of the matrix; if it never happens, the + degree must be the dimension of the entire space. SETUP:: @@ -838,7 +871,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: x = random_eja().random_element() sage: x.degree() == x.minimal_polynomial().degree() True - """ n = self.parent().dimension() @@ -1147,7 +1179,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): P = self.parent() left_mult_by_self = lambda y: self*y L = P.module_morphism(function=left_mult_by_self, codomain=P) - return FiniteDimensionalEJAOperator(P, P, L.matrix() ) + return EJAOperator(P, P, L.matrix() ) def quadratic_representation(self, other=None): @@ -1440,7 +1472,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): if self.is_nilpotent(): raise ValueError("this only works with non-nilpotent elements!") - J = self.subalgebra_generated_by() + # The subalgebra is transient (we return an element of the + # superalgebra, i.e. this algebra) so why bother + # orthonormalizing? + J = self.subalgebra_generated_by(orthonormalize=False) u = J(self) # The image of the matrix of left-u^m-multiplication @@ -1461,14 +1496,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): # subspace... or do we? Can't we just solve, knowing that # A(c) = u^(s+1) should have a solution in the big space, # too? - # - # Beware, solve_right() means that we're using COLUMN vectors. - # Our FiniteDimensionalAlgebraElement superclass uses rows. u_next = u**(s+1) A = u_next.operator().matrix() c = J.from_vector(A.solve_right(u_next.to_vector())) - # Now c is the idempotent we want, but it still lives in the subalgebra. + # Now c is the idempotent we want, but it still lives in + # the subalgebra. return c.superalgebra_element() @@ -1542,7 +1575,75 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): # we want the negative of THAT for the trace. return -p(*self.to_vector()) - def operator_inner_product(self, other): + + def trace_inner_product(self, other): + """ + Return the trace inner product of myself and ``other``. + + SETUP:: + + sage: from mjo.eja.eja_algebra import random_eja + + TESTS: + + The trace inner product is commutative, bilinear, and associative:: + + 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: actual = (a*(x+z)).trace_inner_product(y) + sage: expected = ( a*x.trace_inner_product(y) + + ....: a*z.trace_inner_product(y) ) + sage: actual == expected + True + sage: actual = x.trace_inner_product(a*(y+z)) + sage: expected = ( a*x.trace_inner_product(y) + + ....: a*x.trace_inner_product(z) ) + sage: actual == expected + True + sage: # associative + sage: (x*y).trace_inner_product(z) == y.trace_inner_product(x*z) + True + + """ + if not other in self.parent(): + raise TypeError("'other' must live in the same algebra") + + return (self*other).trace() + + + def trace_norm(self): + """ + The norm of this element with respect to :meth:`trace_inner_product`. + + SETUP:: + + sage: from mjo.eja.eja_algebra import (JordanSpinEJA, + ....: HadamardEJA) + + EXAMPLES:: + + sage: J = HadamardEJA(2) + sage: x = sum(J.gens()) + sage: x.trace_norm() + 1.414213562373095? + + :: + + sage: J = JordanSpinEJA(4) + sage: x = sum(J.gens()) + sage: x.trace_norm() + 2.828427124746190? + + """ + return self.trace_inner_product(self).sqrt() + + + def operator_trace_inner_product(self, other): r""" Return the operator inner product of myself and ``other``. @@ -1554,10 +1655,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): 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). + This 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:: @@ -1577,7 +1678,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: x,y = J.random_elements(2) sage: n = J.dimension() sage: r = J.rank() - sage: actual = x.operator_inner_product(y) + sage: actual = x.operator_trace_inner_product(y) sage: expected = (n/r)*x.trace_inner_product(y) sage: actual == expected True @@ -1588,7 +1689,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: x,y = J.random_elements(2) sage: n = J.dimension() sage: r = J.rank() - sage: actual = x.operator_inner_product(y) + sage: actual = x.operator_trace_inner_product(y) sage: expected = (n/r)*x.trace_inner_product(y) sage: actual == expected True @@ -1599,7 +1700,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): sage: x,y = J.random_elements(2) sage: n = J.dimension() sage: r = J.rank() - sage: actual = x.operator_inner_product(y) + sage: actual = x.operator_trace_inner_product(y) sage: expected = (n/r)*x.trace_inner_product(y) sage: actual == expected True @@ -1612,113 +1713,94 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): 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) + sage: actual = x.operator_trace_inner_product(y) + sage: expected = y.operator_trace_inner_product(x) + sage: actual == expected 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 = (a*(x+z)).operator_trace_inner_product(y) + sage: expected = ( a*x.operator_trace_inner_product(y) + + ....: a*z.operator_trace_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 = x.operator_trace_inner_product(a*(y+z)) + sage: expected = ( a*x.operator_trace_inner_product(y) + + ....: a*x.operator_trace_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 = (x*y).operator_trace_inner_product(z) + sage: expected = y.operator_trace_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() - + Despite the fact that the implementation uses a matrix representation, + the answer is independent of the basis used:: - def trace_inner_product(self, other): - """ - Return the trace inner product of myself and ``other``. - - SETUP:: - - sage: from mjo.eja.eja_algebra import random_eja - - TESTS: - - The trace inner product is commutative, bilinear, and associative:: - - 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: actual = (a*(x+z)).trace_inner_product(y) - sage: expected = ( a*x.trace_inner_product(y) + - ....: a*z.trace_inner_product(y) ) - sage: actual == expected - True - sage: actual = x.trace_inner_product(a*(y+z)) - sage: expected = ( a*x.trace_inner_product(y) + - ....: a*x.trace_inner_product(z) ) + sage: J = RealSymmetricEJA(3, field=QQ, orthonormalize=False) + sage: V = RealSymmetricEJA(3) + sage: x,y = J.random_elements(2) + sage: w = V(x.to_matrix()) + sage: z = V(y.to_matrix()) + sage: expected = x.operator_trace_inner_product(y) + sage: actual = w.operator_trace_inner_product(z) sage: actual == expected True - sage: # associative - sage: (x*y).trace_inner_product(z) == y.trace_inner_product(x*z) - True """ if not other in self.parent(): raise TypeError("'other' must live in the same algebra") - return (self*other).trace() + return (self*other).operator().matrix().trace() - def trace_norm(self): + def operator_trace_norm(self): """ - The norm of this element with respect to :meth:`trace_inner_product`. + The norm of this element with respect to + :meth:`operator_trace_inner_product`. SETUP:: sage: from mjo.eja.eja_algebra import (JordanSpinEJA, ....: HadamardEJA) - EXAMPLES:: + EXAMPLES: + + On a simple algebra, this will differ from :meth:`trace_norm` + by the scalar factor ``(n/r).sqrt()``, where `n` is the + dimension of the algebra and `r` its rank. This follows from + the corresponding result (Proposition III.4.2 of Faraut and + Korányi) for the trace inner product:: sage: J = HadamardEJA(2) sage: x = sum(J.gens()) - sage: x.trace_norm() + sage: x.operator_trace_norm() 1.414213562373095? :: sage: J = JordanSpinEJA(4) sage: x = sum(J.gens()) - sage: x.trace_norm() - 2.828427124746190? + sage: x.operator_trace_norm() + 4 """ - return self.trace_inner_product(self).sqrt() + return self.operator_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() ) +class CartesianProductParentEJAElement(EJAElement): + r""" + An intermediate class for elements that have a Cartesian + product as their parent algebra. + This is needed because the ``to_matrix`` method (which gives you a + representation from the superalgebra) needs to do special stuff + for Cartesian products. Specifically, an EJA subalgebra of a + Cartesian product EJA will not itself be a Cartesian product (it + has its own basis) -- but we want ``to_matrix()`` to be able to + give us a Cartesian product representation. + """ def to_matrix(self): # An override is necessary to call our custom _scale(). B = self.parent().matrix_basis() @@ -1726,7 +1808,20 @@ class CartesianProductEJAElement(FiniteDimensionalEJAElement): # 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 + # that name. This is hidden in a subclass because the # _scale() function is slow. pairs = zip(B, self.to_vector()) return W.sum( _scale(b, alpha) for (b,alpha) in pairs ) + +class CartesianProductEJAElement(CartesianProductParentEJAElement): + 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() )