X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feuclidean_jordan_algebra.py;h=b25839815289ba6f60630a14b790dbfe29f8969b;hb=5978697af0b924d2d2b6174975bf06f1830746f3;hp=3124132b61e5f196268cb027ce1bd26028aaa942;hpb=4d6b173a3ef90c7243e15febf60b0dc45282fd14;p=sage.d.git diff --git a/mjo/eja/euclidean_jordan_algebra.py b/mjo/eja/euclidean_jordan_algebra.py index 3124132..b258398 100644 --- a/mjo/eja/euclidean_jordan_algebra.py +++ b/mjo/eja/euclidean_jordan_algebra.py @@ -21,7 +21,8 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): assume_associative=False, category=None, rank=None, - natural_basis=None): + natural_basis=None, + inner_product=None): n = len(mult_table) mult_table = [b.base_extend(field) for b in mult_table] for b in mult_table: @@ -45,7 +46,8 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): names=names, category=cat, rank=rank, - natural_basis=natural_basis) + natural_basis=natural_basis, + inner_product=inner_product) def __init__(self, field, @@ -54,7 +56,8 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): assume_associative=False, category=None, rank=None, - natural_basis=None): + natural_basis=None, + inner_product=None): """ EXAMPLES: @@ -70,6 +73,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): """ self._rank = rank self._natural_basis = natural_basis + self._inner_product = inner_product fda = super(FiniteDimensionalEuclideanJordanAlgebra, self) fda.__init__(field, mult_table, @@ -85,6 +89,20 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): return fmt.format(self.degree(), self.base_ring()) + def inner_product(self, x, y): + """ + The inner product associated with this Euclidean Jordan algebra. + + Will default to the trace inner product if nothing else. + """ + if (not x in self) or (not y in self): + raise TypeError("arguments must live in this algebra") + if self._inner_product is None: + return x.trace_inner_product(y) + else: + return self._inner_product(x,y) + + def natural_basis(self): """ Return a more-natural representation of this algebra's basis. @@ -142,6 +160,51 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): An element of a Euclidean Jordan algebra. """ + def __init__(self, A, elt=None): + """ + EXAMPLES: + + The identity in `S^n` is converted to the identity in the EJA:: + + sage: J = RealSymmetricSimpleEJA(3) + sage: I = identity_matrix(QQ,3) + sage: J(I) == J.one() + True + + This skew-symmetric matrix can't be represented in the EJA:: + + sage: J = RealSymmetricSimpleEJA(3) + sage: A = matrix(QQ,3, lambda i,j: i-j) + sage: J(A) + Traceback (most recent call last): + ... + ArithmeticError: vector is not in free module + + """ + # Goal: if we're given a matrix, and if it lives in our + # parent algebra's "natural ambient space," convert it + # into an algebra element. + # + # The catch is, we make a recursive call after converting + # the given matrix into a vector that lives in the algebra. + # This we need to try the parent class initializer first, + # to avoid recursing forever if we're given something that + # already fits into the algebra, but also happens to live + # in the parent's "natural ambient space" (this happens with + # vectors in R^n). + try: + FiniteDimensionalAlgebraElement.__init__(self, A, elt) + except ValueError: + natural_basis = A.natural_basis() + if elt in natural_basis[0].matrix_space(): + # Thanks for nothing! Matrix spaces aren't vector + # spaces in Sage, so we have to figure out its + # natural-basis coordinates ourselves. + V = VectorSpace(elt.base_ring(), elt.nrows()**2) + W = V.span( _mat2vec(s) for s in natural_basis ) + coords = W.coordinates(_mat2vec(elt)) + FiniteDimensionalAlgebraElement.__init__(self, A, coords) + def __pow__(self, n): """ Return ``self`` raised to the power ``n``. @@ -208,6 +271,62 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): raise NotImplementedError('irregular element') + def inner_product(self, other): + """ + Return the parent algebra's inner product of myself and ``other``. + + EXAMPLES: + + The inner product in the Jordan spin algebra is the usual + inner product on `R^n` (this example only works because the + basis for the Jordan algebra is the standard basis in `R^n`):: + + sage: J = JordanSpinSimpleEJA(3) + sage: x = vector(QQ,[1,2,3]) + sage: y = vector(QQ,[4,5,6]) + sage: x.inner_product(y) + 32 + sage: J(x).inner_product(J(y)) + 32 + + The inner product on `S^n` is ` = trace(X*Y)`, where + multiplication is the usual matrix multiplication in `S^n`, + so the inner product of the identity matrix with itself + should be the `n`:: + + sage: J = RealSymmetricSimpleEJA(3) + sage: J.one().inner_product(J.one()) + 3 + + Likewise, the inner product on `C^n` is ` = + Re(trace(X*Y))`, where we must necessarily take the real + part because the product of Hermitian matrices may not be + Hermitian:: + + sage: J = ComplexHermitianSimpleEJA(3) + sage: J.one().inner_product(J.one()) + 3 + + TESTS: + + 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 = J.random_element() + sage: y = J.random_element() + sage: x.inner_product(y) in RR + True + + """ + P = self.parent() + if not other in P: + raise TypeError("'other' must live in the same algebra") + + return P.inner_product(self, other) + + def operator_commutes_with(self, other): """ Return whether or not this element operator-commutes @@ -238,7 +357,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): """ if not other in self.parent(): - raise ArgumentError("'other' must live in the same algebra") + raise TypeError("'other' must live in the same algebra") A = self.operator_matrix() B = other.operator_matrix() @@ -328,7 +447,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): # TODO: we can do better once the call to is_invertible() # doesn't crash on irregular elements. #if not self.is_invertible(): - # raise ArgumentError('element is not invertible') + # raise ValueError('element is not invertible') # We do this a little different than the usual recursive # call to a finite-dimensional algebra element, because we @@ -461,6 +580,92 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): return self.span_of_powers().dimension() + def minimal_polynomial(self): + """ + EXAMPLES:: + + sage: set_random_seed() + sage: x = random_eja().random_element() + sage: x.degree() == x.minimal_polynomial().degree() + True + + :: + + sage: set_random_seed() + sage: x = random_eja().random_element() + sage: x.degree() == x.minimal_polynomial().degree() + True + + The minimal polynomial and the characteristic polynomial coincide + and are known (see Alizadeh, Example 11.11) for all elements of + the spin factor algebra that aren't scalar multiples of the + identity:: + + sage: set_random_seed() + sage: n = ZZ.random_element(2,10) + sage: J = JordanSpinSimpleEJA(n) + sage: y = J.random_element() + sage: while y == y.coefficient(0)*J.one(): + ....: y = J.random_element() + sage: y0 = y.vector()[0] + sage: y_bar = y.vector()[1:] + sage: actual = y.minimal_polynomial() + sage: x = SR.symbol('x', domain='real') + sage: expected = x^2 - 2*y0*x + (y0^2 - norm(y_bar)^2) + sage: bool(actual == expected) + True + + """ + # The element we're going to call "minimal_polynomial()" on. + # Either myself, interpreted as an element of a finite- + # dimensional algebra, or an element of an associative + # subalgebra. + elt = None + + if self.parent().is_associative(): + elt = FiniteDimensionalAlgebraElement(self.parent(), self) + else: + V = self.span_of_powers() + assoc_subalg = self.subalgebra_generated_by() + # Mis-design warning: the basis used for span_of_powers() + # and subalgebra_generated_by() must be the same, and in + # the same order! + elt = assoc_subalg(V.coordinates(self.vector())) + + # Recursive call, but should work since elt lives in an + # associative algebra. + return elt.minimal_polynomial() + + + def natural_representation(self): + """ + Return a more-natural representation of this element. + + Every finite-dimensional Euclidean Jordan Algebra is a + direct sum of five simple algebras, four of which comprise + Hermitian matrices. This method returns the original + "natural" representation of this element as a Hermitian + matrix, if it has one. If not, you get the usual representation. + + EXAMPLES:: + + sage: J = ComplexHermitianSimpleEJA(3) + sage: J.one() + e0 + e5 + e8 + sage: J.one().natural_representation() + [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] + + """ + B = self.parent().natural_basis() + W = B[0].matrix_space() + return W.linear_combination(zip(self.vector(), B)) + + def operator_matrix(self): """ Return the matrix that represents left- (or right-) @@ -528,64 +733,6 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): return fda_elt.matrix().transpose() - - def minimal_polynomial(self): - """ - EXAMPLES:: - - sage: set_random_seed() - sage: x = random_eja().random_element() - sage: x.degree() == x.minimal_polynomial().degree() - True - - :: - - sage: set_random_seed() - sage: x = random_eja().random_element() - sage: x.degree() == x.minimal_polynomial().degree() - True - - The minimal polynomial and the characteristic polynomial coincide - and are known (see Alizadeh, Example 11.11) for all elements of - the spin factor algebra that aren't scalar multiples of the - identity:: - - sage: set_random_seed() - sage: n = ZZ.random_element(2,10) - sage: J = JordanSpinSimpleEJA(n) - sage: y = J.random_element() - sage: while y == y.coefficient(0)*J.one(): - ....: y = J.random_element() - sage: y0 = y.vector()[0] - sage: y_bar = y.vector()[1:] - sage: actual = y.minimal_polynomial() - sage: x = SR.symbol('x', domain='real') - sage: expected = x^2 - 2*y0*x + (y0^2 - norm(y_bar)^2) - sage: bool(actual == expected) - True - - """ - # The element we're going to call "minimal_polynomial()" on. - # Either myself, interpreted as an element of a finite- - # dimensional algebra, or an element of an associative - # subalgebra. - elt = None - - if self.parent().is_associative(): - elt = FiniteDimensionalAlgebraElement(self.parent(), self) - else: - V = self.span_of_powers() - assoc_subalg = self.subalgebra_generated_by() - # Mis-design warning: the basis used for span_of_powers() - # and subalgebra_generated_by() must be the same, and in - # the same order! - elt = assoc_subalg(V.coordinates(self.vector())) - - # Recursive call, but should work since elt lives in an - # associative algebra. - return elt.minimal_polynomial() - - def quadratic_representation(self, other=None): """ Return the quadratic representation of this element. @@ -656,7 +803,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): if other is None: other=self elif not other in self.parent(): - raise ArgumentError("'other' must live in the same algebra") + raise TypeError("'other' must live in the same algebra") L = self.operator_matrix() M = other.operator_matrix() @@ -819,7 +966,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): Return the trace inner product of myself and ``other``. """ if not other in self.parent(): - raise ArgumentError("'other' must live in the same algebra") + raise TypeError("'other' must live in the same algebra") return (self*other).trace() @@ -857,7 +1004,10 @@ def eja_rn(dimension, field=QQ): Qs = [ matrix(field, dimension, dimension, lambda k,j: 1*(k == j == i)) for i in xrange(dimension) ] - return FiniteDimensionalEuclideanJordanAlgebra(field,Qs,rank=dimension) + return FiniteDimensionalEuclideanJordanAlgebra(field, + Qs, + rank=dimension, + inner_product=_usual_ip) @@ -952,6 +1102,12 @@ def _complex_hermitian_basis(n, field=QQ): return tuple(S) +def _mat2vec(m): + return vector(m.base_ring(), m.list()) + +def _vec2mat(v): + return matrix(v.base_ring(), sqrt(v.degree()), v.list()) + def _multiplication_table_from_matrix_basis(basis): """ At least three of the five simple Euclidean Jordan algebras have the @@ -972,19 +1128,13 @@ def _multiplication_table_from_matrix_basis(basis): field = basis[0].base_ring() dimension = basis[0].nrows() - def mat2vec(m): - return vector(field, m.list()) - - def vec2mat(v): - return matrix(field, dimension, v.list()) - V = VectorSpace(field, dimension**2) - W = V.span( mat2vec(s) for s in basis ) + W = V.span( _mat2vec(s) for s in basis ) # Taking the span above reorders our basis (thanks, jerk!) so we # need to put our "matrix basis" in the same order as the # (reordered) vector basis. - S = tuple( vec2mat(b) for b in W.basis() ) + S = tuple( _vec2mat(b) for b in W.basis() ) Qs = [] for s in S: @@ -997,7 +1147,7 @@ def _multiplication_table_from_matrix_basis(basis): # why we're computing rows here and not columns. Q_rows = [] for t in S: - this_row = mat2vec((s*t + t*s)/2) + this_row = _mat2vec((s*t + t*s)/2) Q_rows.append(W.coordinates(this_row)) Q = matrix(field, W.dimension(), Q_rows) Qs.append(Q) @@ -1029,7 +1179,7 @@ def _embed_complex_matrix(M): """ n = M.nrows() if M.ncols() != n: - raise ArgumentError("the matrix 'M' must be square") + raise ValueError("the matrix 'M' must be square") field = M.base_ring() blocks = [] for z in M.list(): @@ -1057,9 +1207,9 @@ def _unembed_complex_matrix(M): """ n = ZZ(M.nrows()) if M.ncols() != n: - raise ArgumentError("the matrix 'M' must be square") + raise ValueError("the matrix 'M' must be square") if not n.mod(2).is_zero(): - raise ArgumentError("the matrix 'M' must be a complex embedding") + raise ValueError("the matrix 'M' must be a complex embedding") F = QuadraticField(-1, 'i') i = F.gen() @@ -1071,14 +1221,26 @@ def _unembed_complex_matrix(M): for j in xrange(n/2): submat = M[2*k:2*k+2,2*j:2*j+2] if submat[0,0] != submat[1,1]: - raise ArgumentError('bad real submatrix') + raise ValueError('bad real submatrix') if submat[0,1] != -submat[1,0]: - raise ArgumentError('bad imag submatrix') + raise ValueError('bad imag submatrix') z = submat[0,0] + submat[1,0]*i elements.append(z) return matrix(F, n/2, elements) +# The usual inner product on R^n. +def _usual_ip(x,y): + return x.vector().inner_product(y.vector()) + +# The inner product used for the real symmetric simple EJA. +# We keep it as a separate function because e.g. the complex +# algebra uses the same inner product, except divided by 2. +def _matrix_ip(X,Y): + X_mat = X.natural_representation() + Y_mat = Y.natural_representation() + return (X_mat*Y_mat).trace() + def RealSymmetricSimpleEJA(n, field=QQ): """ @@ -1114,7 +1276,8 @@ def RealSymmetricSimpleEJA(n, field=QQ): return FiniteDimensionalEuclideanJordanAlgebra(field, Qs, rank=n, - natural_basis=T) + natural_basis=T, + inner_product=_matrix_ip) def ComplexHermitianSimpleEJA(n, field=QQ): @@ -1137,10 +1300,21 @@ def ComplexHermitianSimpleEJA(n, field=QQ): """ S = _complex_hermitian_basis(n) (Qs, T) = _multiplication_table_from_matrix_basis(S) + + # Since a+bi on the diagonal is represented as + # + # a + bi = [ a b ] + # [ -b a ], + # + # we'll double-count the "a" entries if we take the trace of + # the embedding. + ip = lambda X,Y: _matrix_ip(X,Y)/2 + return FiniteDimensionalEuclideanJordanAlgebra(field, Qs, rank=n, - natural_basis=T) + natural_basis=T, + inner_product=ip) def QuaternionHermitianSimpleEJA(n): @@ -1189,10 +1363,10 @@ def JordanSpinSimpleEJA(n, field=QQ): In one dimension, this is the reals under multiplication:: - sage: J1 = JordanSpinSimpleEJA(1) - sage: J2 = eja_rn(1) - sage: J1 == J2 - True + sage: J1 = JordanSpinSimpleEJA(1) + sage: J2 = eja_rn(1) + sage: J1 == J2 + True """ Qs = [] @@ -1211,4 +1385,7 @@ def JordanSpinSimpleEJA(n, field=QQ): # The rank of the spin factor algebra is two, UNLESS we're in a # one-dimensional ambient space (the rank is bounded by the # ambient dimension). - return FiniteDimensionalEuclideanJordanAlgebra(field, Qs, rank=min(n,2)) + return FiniteDimensionalEuclideanJordanAlgebra(field, + Qs, + rank=min(n,2), + inner_product=_usual_ip)