2 Euclidean Jordan Algebras. These are formally-real Jordan Algebras;
3 specifically those where u^2 + v^2 = 0 implies that u = v = 0. They
4 are used in optimization, and have some additional nice methods beyond
5 what can be supported in a general Jordan Algebra.
8 from sage
.categories
.magmatic_algebras
import MagmaticAlgebras
9 from sage
.structure
.element
import is_Matrix
10 from sage
.structure
.category_object
import normalize_names
12 from sage
.algebras
.finite_dimensional_algebras
.finite_dimensional_algebra
import FiniteDimensionalAlgebra
13 from sage
.algebras
.finite_dimensional_algebras
.finite_dimensional_algebra_element
import FiniteDimensionalAlgebraElement
15 class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra
):
17 def __classcall_private__(cls
,
21 assume_associative
=False,
26 mult_table
= [b
.base_extend(field
) for b
in mult_table
]
29 if not (is_Matrix(b
) and b
.dimensions() == (n
, n
)):
30 raise ValueError("input is not a multiplication table")
31 mult_table
= tuple(mult_table
)
33 cat
= MagmaticAlgebras(field
).FiniteDimensional().WithBasis()
34 cat
.or_subcategory(category
)
35 if assume_associative
:
36 cat
= cat
.Associative()
38 names
= normalize_names(n
, names
)
40 fda
= super(FiniteDimensionalEuclideanJordanAlgebra
, cls
)
41 return fda
.__classcall
__(cls
,
44 assume_associative
=assume_associative
,
48 natural_basis
=natural_basis
)
55 assume_associative
=False,
62 By definition, Jordan multiplication commutes::
64 sage: set_random_seed()
65 sage: J = random_eja()
66 sage: x = J.random_element()
67 sage: y = J.random_element()
73 self
._natural
_basis
= natural_basis
74 self
._multiplication
_table
= mult_table
75 fda
= super(FiniteDimensionalEuclideanJordanAlgebra
, self
)
84 Return a string representation of ``self``.
86 fmt
= "Euclidean Jordan algebra of degree {} over {}"
87 return fmt
.format(self
.degree(), self
.base_ring())
92 def characteristic_polynomial(self
):
96 The characteristic polynomial in the spin algebra is given in
97 Alizadeh, Example 11.11::
99 sage: J = JordanSpinEJA(3)
100 sage: p = J.characteristic_polynomial(); p
101 X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
102 sage: xvec = J.one().vector()
110 # Now switch to the polynomial rings.
111 names
= ['X' + str(i
) for i
in range(1,n
+1)]
112 R
= PolynomialRing(self
.base_ring(), names
)
113 J
= FiniteDimensionalEuclideanJordanAlgebra(R
,
114 self
._multiplication
_table
,
118 # The coordinates of e_k with respect to the basis B.
119 # But, the e_k are elements of B...
120 return identity_matrix(J
.base_ring(), n
).column(k
-1).column()
122 # A matrix implementation 1
123 x
= J(vector(R
, R
.gens()))
124 l1
= [column_matrix((x
**k
).vector()) for k
in range(r
)]
125 l2
= [e(k
) for k
in range(r
+1, n
+1)]
126 A_of_x
= block_matrix(1, n
, (l1
+ l2
))
129 denominator
= A_of_x
.det() # This is constant
131 A_cols
= A_of_x
.columns()
133 numerator
= column_matrix(A_of_x
.base_ring(), A_cols
).det()
134 ai
= numerator
/denominator
137 # We go to a bit of trouble here to reorder the
138 # indeterminates, so that it's easier to evaluate the
139 # characteristic polynomial at x's coordinates and get back
140 # something in terms of t, which is what we want.
141 S
= PolynomialRing(self
.base_ring(),'t')
143 S
= PolynomialRing(S
, R
.variable_names())
146 # We're relying on the theory here to ensure that each entry
147 # a[i] is indeed back in R, and the added negative signs are
148 # to make the whole expression sum to zero.
149 a
= [R(-ai
) for ai
in a
] # corresponds to powerx x^0 through x^(r-1)
151 # Note: all entries past the rth should be zero. The
152 # coefficient of the highest power (x^r) is 1, but it doesn't
153 # appear in the solution vector which contains coefficients
154 # for the other powers (to make them sum to x^r).
156 a
[r
] = 1 # corresponds to x^r
158 # When the rank is equal to the dimension, trying to
159 # assign a[r] goes out-of-bounds.
160 a
.append(1) # corresponds to x^r
162 return sum( a
[k
]*(t
**k
) for k
in range(len(a
)) )
165 def inner_product(self
, x
, y
):
167 The inner product associated with this Euclidean Jordan algebra.
169 Defaults to the trace inner product, but can be overridden by
170 subclasses if they are sure that the necessary properties are
175 The inner product must satisfy its axiom for this algebra to truly
176 be a Euclidean Jordan Algebra::
178 sage: set_random_seed()
179 sage: J = random_eja()
180 sage: x = J.random_element()
181 sage: y = J.random_element()
182 sage: z = J.random_element()
183 sage: (x*y).inner_product(z) == y.inner_product(x*z)
187 if (not x
in self
) or (not y
in self
):
188 raise TypeError("arguments must live in this algebra")
189 return x
.trace_inner_product(y
)
192 def natural_basis(self
):
194 Return a more-natural representation of this algebra's basis.
196 Every finite-dimensional Euclidean Jordan Algebra is a direct
197 sum of five simple algebras, four of which comprise Hermitian
198 matrices. This method returns the original "natural" basis
199 for our underlying vector space. (Typically, the natural basis
200 is used to construct the multiplication table in the first place.)
202 Note that this will always return a matrix. The standard basis
203 in `R^n` will be returned as `n`-by-`1` column matrices.
207 sage: J = RealSymmetricEJA(2)
210 sage: J.natural_basis()
218 sage: J = JordanSpinEJA(2)
221 sage: J.natural_basis()
228 if self
._natural
_basis
is None:
229 return tuple( b
.vector().column() for b
in self
.basis() )
231 return self
._natural
_basis
236 Return the rank of this EJA.
238 if self
._rank
is None:
239 raise ValueError("no rank specified at genesis")
244 class Element(FiniteDimensionalAlgebraElement
):
246 An element of a Euclidean Jordan algebra.
249 def __init__(self
, A
, elt
=None):
253 The identity in `S^n` is converted to the identity in the EJA::
255 sage: J = RealSymmetricEJA(3)
256 sage: I = identity_matrix(QQ,3)
257 sage: J(I) == J.one()
260 This skew-symmetric matrix can't be represented in the EJA::
262 sage: J = RealSymmetricEJA(3)
263 sage: A = matrix(QQ,3, lambda i,j: i-j)
265 Traceback (most recent call last):
267 ArithmeticError: vector is not in free module
270 # Goal: if we're given a matrix, and if it lives in our
271 # parent algebra's "natural ambient space," convert it
272 # into an algebra element.
274 # The catch is, we make a recursive call after converting
275 # the given matrix into a vector that lives in the algebra.
276 # This we need to try the parent class initializer first,
277 # to avoid recursing forever if we're given something that
278 # already fits into the algebra, but also happens to live
279 # in the parent's "natural ambient space" (this happens with
282 FiniteDimensionalAlgebraElement
.__init
__(self
, A
, elt
)
284 natural_basis
= A
.natural_basis()
285 if elt
in natural_basis
[0].matrix_space():
286 # Thanks for nothing! Matrix spaces aren't vector
287 # spaces in Sage, so we have to figure out its
288 # natural-basis coordinates ourselves.
289 V
= VectorSpace(elt
.base_ring(), elt
.nrows()**2)
290 W
= V
.span( _mat2vec(s
) for s
in natural_basis
)
291 coords
= W
.coordinates(_mat2vec(elt
))
292 FiniteDimensionalAlgebraElement
.__init
__(self
, A
, coords
)
294 def __pow__(self
, n
):
296 Return ``self`` raised to the power ``n``.
298 Jordan algebras are always power-associative; see for
299 example Faraut and Koranyi, Proposition II.1.2 (ii).
303 We have to override this because our superclass uses row vectors
304 instead of column vectors! We, on the other hand, assume column
309 sage: set_random_seed()
310 sage: x = random_eja().random_element()
311 sage: x.operator_matrix()*x.vector() == (x^2).vector()
314 A few examples of power-associativity::
316 sage: set_random_seed()
317 sage: x = random_eja().random_element()
318 sage: x*(x*x)*(x*x) == x^5
320 sage: (x*x)*(x*x*x) == x^5
323 We also know that powers operator-commute (Koecher, Chapter
326 sage: set_random_seed()
327 sage: x = random_eja().random_element()
328 sage: m = ZZ.random_element(0,10)
329 sage: n = ZZ.random_element(0,10)
330 sage: Lxm = (x^m).operator_matrix()
331 sage: Lxn = (x^n).operator_matrix()
332 sage: Lxm*Lxn == Lxn*Lxm
342 return A( (self
.operator_matrix()**(n
-1))*self
.vector() )
345 def apply_univariate_polynomial(self
, p
):
347 Apply the univariate polynomial ``p`` to this element.
349 A priori, SageMath won't allow us to apply a univariate
350 polynomial to an element of an EJA, because we don't know
351 that EJAs are rings (they are usually not associative). Of
352 course, we know that EJAs are power-associative, so the
353 operation is ultimately kosher. This function sidesteps
354 the CAS to get the answer we want and expect.
358 sage: R = PolynomialRing(QQ, 't')
360 sage: p = t^4 - t^3 + 5*t - 2
361 sage: J = RealCartesianProductEJA(5)
362 sage: J.one().apply_univariate_polynomial(p) == 3*J.one()
367 We should always get back an element of the algebra::
369 sage: set_random_seed()
370 sage: p = PolynomialRing(QQ, 't').random_element()
371 sage: J = random_eja()
372 sage: x = J.random_element()
373 sage: x.apply_univariate_polynomial(p) in J
377 if len(p
.variables()) > 1:
378 raise ValueError("not a univariate polynomial")
381 # Convert the coeficcients to the parent's base ring,
382 # because a priori they might live in an (unnecessarily)
383 # larger ring for which P.sum() would fail below.
384 cs
= [ R(c
) for c
in p
.coefficients(sparse
=False) ]
385 return P
.sum( cs
[k
]*(self
**k
) for k
in range(len(cs
)) )
388 def characteristic_polynomial(self
):
390 Return the characteristic polynomial of this element.
394 The rank of `R^3` is three, and the minimal polynomial of
395 the identity element is `(t-1)` from which it follows that
396 the characteristic polynomial should be `(t-1)^3`::
398 sage: J = RealCartesianProductEJA(3)
399 sage: J.one().characteristic_polynomial()
400 t^3 - 3*t^2 + 3*t - 1
402 Likewise, the characteristic of the zero element in the
403 rank-three algebra `R^{n}` should be `t^{3}`::
405 sage: J = RealCartesianProductEJA(3)
406 sage: J.zero().characteristic_polynomial()
409 The characteristic polynomial of an element should evaluate
410 to zero on that element::
412 sage: set_random_seed()
413 sage: x = RealCartesianProductEJA(3).random_element()
414 sage: p = x.characteristic_polynomial()
415 sage: x.apply_univariate_polynomial(p)
419 p
= self
.parent().characteristic_polynomial()
420 return p(*self
.vector())
423 def inner_product(self
, other
):
425 Return the parent algebra's inner product of myself and ``other``.
429 The inner product in the Jordan spin algebra is the usual
430 inner product on `R^n` (this example only works because the
431 basis for the Jordan algebra is the standard basis in `R^n`)::
433 sage: J = JordanSpinEJA(3)
434 sage: x = vector(QQ,[1,2,3])
435 sage: y = vector(QQ,[4,5,6])
436 sage: x.inner_product(y)
438 sage: J(x).inner_product(J(y))
441 The inner product on `S^n` is `<X,Y> = trace(X*Y)`, where
442 multiplication is the usual matrix multiplication in `S^n`,
443 so the inner product of the identity matrix with itself
446 sage: J = RealSymmetricEJA(3)
447 sage: J.one().inner_product(J.one())
450 Likewise, the inner product on `C^n` is `<X,Y> =
451 Re(trace(X*Y))`, where we must necessarily take the real
452 part because the product of Hermitian matrices may not be
455 sage: J = ComplexHermitianEJA(3)
456 sage: J.one().inner_product(J.one())
459 Ditto for the quaternions::
461 sage: J = QuaternionHermitianEJA(3)
462 sage: J.one().inner_product(J.one())
467 Ensure that we can always compute an inner product, and that
468 it gives us back a real number::
470 sage: set_random_seed()
471 sage: J = random_eja()
472 sage: x = J.random_element()
473 sage: y = J.random_element()
474 sage: x.inner_product(y) in RR
480 raise TypeError("'other' must live in the same algebra")
482 return P
.inner_product(self
, other
)
485 def operator_commutes_with(self
, other
):
487 Return whether or not this element operator-commutes
492 The definition of a Jordan algebra says that any element
493 operator-commutes with its square::
495 sage: set_random_seed()
496 sage: x = random_eja().random_element()
497 sage: x.operator_commutes_with(x^2)
502 Test Lemma 1 from Chapter III of Koecher::
504 sage: set_random_seed()
505 sage: J = random_eja()
506 sage: u = J.random_element()
507 sage: v = J.random_element()
508 sage: lhs = u.operator_commutes_with(u*v)
509 sage: rhs = v.operator_commutes_with(u^2)
514 if not other
in self
.parent():
515 raise TypeError("'other' must live in the same algebra")
517 A
= self
.operator_matrix()
518 B
= other
.operator_matrix()
524 Return my determinant, the product of my eigenvalues.
528 sage: J = JordanSpinEJA(2)
529 sage: e0,e1 = J.gens()
533 sage: J = JordanSpinEJA(3)
534 sage: e0,e1,e2 = J.gens()
535 sage: x = e0 + e1 + e2
540 cs
= self
.characteristic_polynomial().coefficients(sparse
=False)
543 return cs
[0] * (-1)**r
545 raise ValueError('charpoly had no coefficients')
550 Return the Jordan-multiplicative inverse of this element.
552 We can't use the superclass method because it relies on the
553 algebra being associative.
557 The inverse in the spin factor algebra is given in Alizadeh's
560 sage: set_random_seed()
561 sage: n = ZZ.random_element(1,10)
562 sage: J = JordanSpinEJA(n)
563 sage: x = J.random_element()
564 sage: while x.is_zero():
565 ....: x = J.random_element()
566 sage: x_vec = x.vector()
568 sage: x_bar = x_vec[1:]
569 sage: coeff = 1/(x0^2 - x_bar.inner_product(x_bar))
570 sage: inv_vec = x_vec.parent()([x0] + (-x_bar).list())
571 sage: x_inverse = coeff*inv_vec
572 sage: x.inverse() == J(x_inverse)
577 The identity element is its own inverse::
579 sage: set_random_seed()
580 sage: J = random_eja()
581 sage: J.one().inverse() == J.one()
584 If an element has an inverse, it acts like one. TODO: this
585 can be a lot less ugly once ``is_invertible`` doesn't crash
586 on irregular elements::
588 sage: set_random_seed()
589 sage: J = random_eja()
590 sage: x = J.random_element()
592 ....: x.inverse()*x == J.one()
598 if self
.parent().is_associative():
599 elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
602 # TODO: we can do better once the call to is_invertible()
603 # doesn't crash on irregular elements.
604 #if not self.is_invertible():
605 # raise ValueError('element is not invertible')
607 # We do this a little different than the usual recursive
608 # call to a finite-dimensional algebra element, because we
609 # wind up with an inverse that lives in the subalgebra and
610 # we need information about the parent to convert it back.
611 V
= self
.span_of_powers()
612 assoc_subalg
= self
.subalgebra_generated_by()
613 # Mis-design warning: the basis used for span_of_powers()
614 # and subalgebra_generated_by() must be the same, and in
616 elt
= assoc_subalg(V
.coordinates(self
.vector()))
618 # This will be in the subalgebra's coordinates...
619 fda_elt
= FiniteDimensionalAlgebraElement(assoc_subalg
, elt
)
620 subalg_inverse
= fda_elt
.inverse()
622 # So we have to convert back...
623 basis
= [ self
.parent(v
) for v
in V
.basis() ]
624 pairs
= zip(subalg_inverse
.vector(), basis
)
625 return self
.parent().linear_combination(pairs
)
628 def is_invertible(self
):
630 Return whether or not this element is invertible.
632 We can't use the superclass method because it relies on
633 the algebra being associative.
637 The usual way to do this is to check if the determinant is
638 zero, but we need the characteristic polynomial for the
639 determinant. The minimal polynomial is a lot easier to get,
640 so we use Corollary 2 in Chapter V of Koecher to check
641 whether or not the paren't algebra's zero element is a root
642 of this element's minimal polynomial.
646 The identity element is always invertible::
648 sage: set_random_seed()
649 sage: J = random_eja()
650 sage: J.one().is_invertible()
653 The zero element is never invertible::
655 sage: set_random_seed()
656 sage: J = random_eja()
657 sage: J.zero().is_invertible()
661 zero
= self
.parent().zero()
662 p
= self
.minimal_polynomial()
663 return not (p(zero
) == zero
)
666 def is_nilpotent(self
):
668 Return whether or not some power of this element is zero.
670 The superclass method won't work unless we're in an
671 associative algebra, and we aren't. However, we generate
672 an assocoative subalgebra and we're nilpotent there if and
673 only if we're nilpotent here (probably).
677 The identity element is never nilpotent::
679 sage: set_random_seed()
680 sage: random_eja().one().is_nilpotent()
683 The additive identity is always nilpotent::
685 sage: set_random_seed()
686 sage: random_eja().zero().is_nilpotent()
690 # The element we're going to call "is_nilpotent()" on.
691 # Either myself, interpreted as an element of a finite-
692 # dimensional algebra, or an element of an associative
696 if self
.parent().is_associative():
697 elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
699 V
= self
.span_of_powers()
700 assoc_subalg
= self
.subalgebra_generated_by()
701 # Mis-design warning: the basis used for span_of_powers()
702 # and subalgebra_generated_by() must be the same, and in
704 elt
= assoc_subalg(V
.coordinates(self
.vector()))
706 # Recursive call, but should work since elt lives in an
707 # associative algebra.
708 return elt
.is_nilpotent()
711 def is_regular(self
):
713 Return whether or not this is a regular element.
717 The identity element always has degree one, but any element
718 linearly-independent from it is regular::
720 sage: J = JordanSpinEJA(5)
721 sage: J.one().is_regular()
723 sage: e0, e1, e2, e3, e4 = J.gens() # e0 is the identity
724 sage: for x in J.gens():
725 ....: (J.one() + x).is_regular()
733 return self
.degree() == self
.parent().rank()
738 Compute the degree of this element the straightforward way
739 according to the definition; by appending powers to a list
740 and figuring out its dimension (that is, whether or not
741 they're linearly dependent).
745 sage: J = JordanSpinEJA(4)
746 sage: J.one().degree()
748 sage: e0,e1,e2,e3 = J.gens()
749 sage: (e0 - e1).degree()
752 In the spin factor algebra (of rank two), all elements that
753 aren't multiples of the identity are regular::
755 sage: set_random_seed()
756 sage: n = ZZ.random_element(1,10)
757 sage: J = JordanSpinEJA(n)
758 sage: x = J.random_element()
759 sage: x == x.coefficient(0)*J.one() or x.degree() == 2
763 return self
.span_of_powers().dimension()
766 def minimal_polynomial(self
):
768 Return the minimal polynomial of this element,
769 as a function of the variable `t`.
773 We restrict ourselves to the associative subalgebra
774 generated by this element, and then return the minimal
775 polynomial of this element's operator matrix (in that
776 subalgebra). This works by Baes Proposition 2.3.16.
780 The minimal polynomial of the identity and zero elements are
783 sage: set_random_seed()
784 sage: J = random_eja()
785 sage: J.one().minimal_polynomial()
787 sage: J.zero().minimal_polynomial()
790 The degree of an element is (by one definition) the degree
791 of its minimal polynomial::
793 sage: set_random_seed()
794 sage: x = random_eja().random_element()
795 sage: x.degree() == x.minimal_polynomial().degree()
798 The minimal polynomial and the characteristic polynomial coincide
799 and are known (see Alizadeh, Example 11.11) for all elements of
800 the spin factor algebra that aren't scalar multiples of the
803 sage: set_random_seed()
804 sage: n = ZZ.random_element(2,10)
805 sage: J = JordanSpinEJA(n)
806 sage: y = J.random_element()
807 sage: while y == y.coefficient(0)*J.one():
808 ....: y = J.random_element()
809 sage: y0 = y.vector()[0]
810 sage: y_bar = y.vector()[1:]
811 sage: actual = y.minimal_polynomial()
812 sage: t = PolynomialRing(J.base_ring(),'t').gen(0)
813 sage: expected = t^2 - 2*y0*t + (y0^2 - norm(y_bar)^2)
814 sage: bool(actual == expected)
817 The minimal polynomial should always kill its element::
819 sage: set_random_seed()
820 sage: x = random_eja().random_element()
821 sage: p = x.minimal_polynomial()
822 sage: x.apply_univariate_polynomial(p)
826 V
= self
.span_of_powers()
827 assoc_subalg
= self
.subalgebra_generated_by()
828 # Mis-design warning: the basis used for span_of_powers()
829 # and subalgebra_generated_by() must be the same, and in
831 elt
= assoc_subalg(V
.coordinates(self
.vector()))
833 # We get back a symbolic polynomial in 'x' but want a real
835 p_of_x
= elt
.operator_matrix().minimal_polynomial()
836 return p_of_x
.change_variable_name('t')
839 def natural_representation(self
):
841 Return a more-natural representation of this element.
843 Every finite-dimensional Euclidean Jordan Algebra is a
844 direct sum of five simple algebras, four of which comprise
845 Hermitian matrices. This method returns the original
846 "natural" representation of this element as a Hermitian
847 matrix, if it has one. If not, you get the usual representation.
851 sage: J = ComplexHermitianEJA(3)
854 sage: J.one().natural_representation()
864 sage: J = QuaternionHermitianEJA(3)
867 sage: J.one().natural_representation()
868 [1 0 0 0 0 0 0 0 0 0 0 0]
869 [0 1 0 0 0 0 0 0 0 0 0 0]
870 [0 0 1 0 0 0 0 0 0 0 0 0]
871 [0 0 0 1 0 0 0 0 0 0 0 0]
872 [0 0 0 0 1 0 0 0 0 0 0 0]
873 [0 0 0 0 0 1 0 0 0 0 0 0]
874 [0 0 0 0 0 0 1 0 0 0 0 0]
875 [0 0 0 0 0 0 0 1 0 0 0 0]
876 [0 0 0 0 0 0 0 0 1 0 0 0]
877 [0 0 0 0 0 0 0 0 0 1 0 0]
878 [0 0 0 0 0 0 0 0 0 0 1 0]
879 [0 0 0 0 0 0 0 0 0 0 0 1]
882 B
= self
.parent().natural_basis()
883 W
= B
[0].matrix_space()
884 return W
.linear_combination(zip(self
.vector(), B
))
887 def operator_matrix(self
):
889 Return the matrix that represents left- (or right-)
890 multiplication by this element in the parent algebra.
892 We have to override this because the superclass method
893 returns a matrix that acts on row vectors (that is, on
898 Test the first polarization identity from my notes, Koecher Chapter
899 III, or from Baes (2.3)::
901 sage: set_random_seed()
902 sage: J = random_eja()
903 sage: x = J.random_element()
904 sage: y = J.random_element()
905 sage: Lx = x.operator_matrix()
906 sage: Ly = y.operator_matrix()
907 sage: Lxx = (x*x).operator_matrix()
908 sage: Lxy = (x*y).operator_matrix()
909 sage: bool(2*Lx*Lxy + Ly*Lxx == 2*Lxy*Lx + Lxx*Ly)
912 Test the second polarization identity from my notes or from
915 sage: set_random_seed()
916 sage: J = random_eja()
917 sage: x = J.random_element()
918 sage: y = J.random_element()
919 sage: z = J.random_element()
920 sage: Lx = x.operator_matrix()
921 sage: Ly = y.operator_matrix()
922 sage: Lz = z.operator_matrix()
923 sage: Lzy = (z*y).operator_matrix()
924 sage: Lxy = (x*y).operator_matrix()
925 sage: Lxz = (x*z).operator_matrix()
926 sage: bool(Lx*Lzy + Lz*Lxy + Ly*Lxz == Lzy*Lx + Lxy*Lz + Lxz*Ly)
929 Test the third polarization identity from my notes or from
932 sage: set_random_seed()
933 sage: J = random_eja()
934 sage: u = J.random_element()
935 sage: y = J.random_element()
936 sage: z = J.random_element()
937 sage: Lu = u.operator_matrix()
938 sage: Ly = y.operator_matrix()
939 sage: Lz = z.operator_matrix()
940 sage: Lzy = (z*y).operator_matrix()
941 sage: Luy = (u*y).operator_matrix()
942 sage: Luz = (u*z).operator_matrix()
943 sage: Luyz = (u*(y*z)).operator_matrix()
944 sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz
945 sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly
946 sage: bool(lhs == rhs)
950 fda_elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
951 return fda_elt
.matrix().transpose()
954 def quadratic_representation(self
, other
=None):
956 Return the quadratic representation of this element.
960 The explicit form in the spin factor algebra is given by
961 Alizadeh's Example 11.12::
963 sage: set_random_seed()
964 sage: n = ZZ.random_element(1,10)
965 sage: J = JordanSpinEJA(n)
966 sage: x = J.random_element()
967 sage: x_vec = x.vector()
969 sage: x_bar = x_vec[1:]
970 sage: A = matrix(QQ, 1, [x_vec.inner_product(x_vec)])
971 sage: B = 2*x0*x_bar.row()
972 sage: C = 2*x0*x_bar.column()
973 sage: D = identity_matrix(QQ, n-1)
974 sage: D = (x0^2 - x_bar.inner_product(x_bar))*D
975 sage: D = D + 2*x_bar.tensor_product(x_bar)
976 sage: Q = block_matrix(2,2,[A,B,C,D])
977 sage: Q == x.quadratic_representation()
980 Test all of the properties from Theorem 11.2 in Alizadeh::
982 sage: set_random_seed()
983 sage: J = random_eja()
984 sage: x = J.random_element()
985 sage: y = J.random_element()
989 sage: actual = x.quadratic_representation(y)
990 sage: expected = ( (x+y).quadratic_representation()
991 ....: -x.quadratic_representation()
992 ....: -y.quadratic_representation() ) / 2
993 sage: actual == expected
998 sage: alpha = QQ.random_element()
999 sage: actual = (alpha*x).quadratic_representation()
1000 sage: expected = (alpha^2)*x.quadratic_representation()
1001 sage: actual == expected
1006 sage: Qy = y.quadratic_representation()
1007 sage: actual = J(Qy*x.vector()).quadratic_representation()
1008 sage: expected = Qy*x.quadratic_representation()*Qy
1009 sage: actual == expected
1014 sage: k = ZZ.random_element(1,10)
1015 sage: actual = (x^k).quadratic_representation()
1016 sage: expected = (x.quadratic_representation())^k
1017 sage: actual == expected
1023 elif not other
in self
.parent():
1024 raise TypeError("'other' must live in the same algebra")
1026 L
= self
.operator_matrix()
1027 M
= other
.operator_matrix()
1028 return ( L
*M
+ M
*L
- (self
*other
).operator_matrix() )
1031 def span_of_powers(self
):
1033 Return the vector space spanned by successive powers of
1036 # The dimension of the subalgebra can't be greater than
1037 # the big algebra, so just put everything into a list
1038 # and let span() get rid of the excess.
1040 # We do the extra ambient_vector_space() in case we're messing
1041 # with polynomials and the direct parent is a module.
1042 V
= self
.vector().parent().ambient_vector_space()
1043 return V
.span( (self
**d
).vector() for d
in xrange(V
.dimension()) )
1046 def subalgebra_generated_by(self
):
1048 Return the associative subalgebra of the parent EJA generated
1053 sage: set_random_seed()
1054 sage: x = random_eja().random_element()
1055 sage: x.subalgebra_generated_by().is_associative()
1058 Squaring in the subalgebra should be the same thing as
1059 squaring in the superalgebra::
1061 sage: set_random_seed()
1062 sage: x = random_eja().random_element()
1063 sage: u = x.subalgebra_generated_by().random_element()
1064 sage: u.operator_matrix()*u.vector() == (u**2).vector()
1068 # First get the subspace spanned by the powers of myself...
1069 V
= self
.span_of_powers()
1070 F
= self
.base_ring()
1072 # Now figure out the entries of the right-multiplication
1073 # matrix for the successive basis elements b0, b1,... of
1076 for b_right
in V
.basis():
1077 eja_b_right
= self
.parent()(b_right
)
1079 # The first row of the right-multiplication matrix by
1080 # b1 is what we get if we apply that matrix to b1. The
1081 # second row of the right multiplication matrix by b1
1082 # is what we get when we apply that matrix to b2...
1084 # IMPORTANT: this assumes that all vectors are COLUMN
1085 # vectors, unlike our superclass (which uses row vectors).
1086 for b_left
in V
.basis():
1087 eja_b_left
= self
.parent()(b_left
)
1088 # Multiply in the original EJA, but then get the
1089 # coordinates from the subalgebra in terms of its
1091 this_row
= V
.coordinates((eja_b_left
*eja_b_right
).vector())
1092 b_right_rows
.append(this_row
)
1093 b_right_matrix
= matrix(F
, b_right_rows
)
1094 mats
.append(b_right_matrix
)
1096 # It's an algebra of polynomials in one element, and EJAs
1097 # are power-associative.
1099 # TODO: choose generator names intelligently.
1100 return FiniteDimensionalEuclideanJordanAlgebra(F
, mats
, assume_associative
=True, names
='f')
1103 def subalgebra_idempotent(self
):
1105 Find an idempotent in the associative subalgebra I generate
1106 using Proposition 2.3.5 in Baes.
1110 sage: set_random_seed()
1111 sage: J = RealCartesianProductEJA(5)
1112 sage: c = J.random_element().subalgebra_idempotent()
1115 sage: J = JordanSpinEJA(5)
1116 sage: c = J.random_element().subalgebra_idempotent()
1121 if self
.is_nilpotent():
1122 raise ValueError("this only works with non-nilpotent elements!")
1124 V
= self
.span_of_powers()
1125 J
= self
.subalgebra_generated_by()
1126 # Mis-design warning: the basis used for span_of_powers()
1127 # and subalgebra_generated_by() must be the same, and in
1129 u
= J(V
.coordinates(self
.vector()))
1131 # The image of the matrix of left-u^m-multiplication
1132 # will be minimal for some natural number s...
1134 minimal_dim
= V
.dimension()
1135 for i
in xrange(1, V
.dimension()):
1136 this_dim
= (u
**i
).operator_matrix().image().dimension()
1137 if this_dim
< minimal_dim
:
1138 minimal_dim
= this_dim
1141 # Now minimal_matrix should correspond to the smallest
1142 # non-zero subspace in Baes's (or really, Koecher's)
1145 # However, we need to restrict the matrix to work on the
1146 # subspace... or do we? Can't we just solve, knowing that
1147 # A(c) = u^(s+1) should have a solution in the big space,
1150 # Beware, solve_right() means that we're using COLUMN vectors.
1151 # Our FiniteDimensionalAlgebraElement superclass uses rows.
1153 A
= u_next
.operator_matrix()
1154 c_coordinates
= A
.solve_right(u_next
.vector())
1156 # Now c_coordinates is the idempotent we want, but it's in
1157 # the coordinate system of the subalgebra.
1159 # We need the basis for J, but as elements of the parent algebra.
1161 basis
= [self
.parent(v
) for v
in V
.basis()]
1162 return self
.parent().linear_combination(zip(c_coordinates
, basis
))
1167 Return my trace, the sum of my eigenvalues.
1171 sage: J = JordanSpinEJA(3)
1172 sage: e0,e1,e2 = J.gens()
1173 sage: x = e0 + e1 + e2
1178 cs
= self
.characteristic_polynomial().coefficients(sparse
=False)
1182 raise ValueError('charpoly had fewer than 2 coefficients')
1185 def trace_inner_product(self
, other
):
1187 Return the trace inner product of myself and ``other``.
1189 if not other
in self
.parent():
1190 raise TypeError("'other' must live in the same algebra")
1192 return (self
*other
).trace()
1195 class RealCartesianProductEJA(FiniteDimensionalEuclideanJordanAlgebra
):
1197 Return the Euclidean Jordan Algebra corresponding to the set
1198 `R^n` under the Hadamard product.
1200 Note: this is nothing more than the Cartesian product of ``n``
1201 copies of the spin algebra. Once Cartesian product algebras
1202 are implemented, this can go.
1206 This multiplication table can be verified by hand::
1208 sage: J = RealCartesianProductEJA(3)
1209 sage: e0,e1,e2 = J.gens()
1225 def __classcall_private__(cls
, n
, field
=QQ
):
1226 # The FiniteDimensionalAlgebra constructor takes a list of
1227 # matrices, the ith representing right multiplication by the ith
1228 # basis element in the vector space. So if e_1 = (1,0,0), then
1229 # right (Hadamard) multiplication of x by e_1 picks out the first
1230 # component of x; and likewise for the ith basis element e_i.
1231 Qs
= [ matrix(field
, n
, n
, lambda k
,j
: 1*(k
== j
== i
))
1232 for i
in xrange(n
) ]
1234 fdeja
= super(RealCartesianProductEJA
, cls
)
1235 return fdeja
.__classcall
_private
__(cls
, field
, Qs
, rank
=n
)
1237 def inner_product(self
, x
, y
):
1238 return _usual_ip(x
,y
)
1243 Return a "random" finite-dimensional Euclidean Jordan Algebra.
1247 For now, we choose a random natural number ``n`` (greater than zero)
1248 and then give you back one of the following:
1250 * The cartesian product of the rational numbers ``n`` times; this is
1251 ``QQ^n`` with the Hadamard product.
1253 * The Jordan spin algebra on ``QQ^n``.
1255 * The ``n``-by-``n`` rational symmetric matrices with the symmetric
1258 * The ``n``-by-``n`` complex-rational Hermitian matrices embedded
1259 in the space of ``2n``-by-``2n`` real symmetric matrices.
1261 * The ``n``-by-``n`` quaternion-rational Hermitian matrices embedded
1262 in the space of ``4n``-by-``4n`` real symmetric matrices.
1264 Later this might be extended to return Cartesian products of the
1270 Euclidean Jordan algebra of degree...
1274 # The max_n component lets us choose different upper bounds on the
1275 # value "n" that gets passed to the constructor. This is needed
1276 # because e.g. R^{10} is reasonable to test, while the Hermitian
1277 # 10-by-10 quaternion matrices are not.
1278 (constructor
, max_n
) = choice([(RealCartesianProductEJA
, 6),
1280 (RealSymmetricEJA
, 5),
1281 (ComplexHermitianEJA
, 4),
1282 (QuaternionHermitianEJA
, 3)])
1283 n
= ZZ
.random_element(1, max_n
)
1284 return constructor(n
, field
=QQ
)
1288 def _real_symmetric_basis(n
, field
=QQ
):
1290 Return a basis for the space of real symmetric n-by-n matrices.
1292 # The basis of symmetric matrices, as matrices, in their R^(n-by-n)
1296 for j
in xrange(i
+1):
1297 Eij
= matrix(field
, n
, lambda k
,l
: k
==i
and l
==j
)
1301 # Beware, orthogonal but not normalized!
1302 Sij
= Eij
+ Eij
.transpose()
1307 def _complex_hermitian_basis(n
, field
=QQ
):
1309 Returns a basis for the space of complex Hermitian n-by-n matrices.
1313 sage: set_random_seed()
1314 sage: n = ZZ.random_element(1,5)
1315 sage: all( M.is_symmetric() for M in _complex_hermitian_basis(n) )
1319 F
= QuadraticField(-1, 'I')
1322 # This is like the symmetric case, but we need to be careful:
1324 # * We want conjugate-symmetry, not just symmetry.
1325 # * The diagonal will (as a result) be real.
1329 for j
in xrange(i
+1):
1330 Eij
= matrix(field
, n
, lambda k
,l
: k
==i
and l
==j
)
1332 Sij
= _embed_complex_matrix(Eij
)
1335 # Beware, orthogonal but not normalized! The second one
1336 # has a minus because it's conjugated.
1337 Sij_real
= _embed_complex_matrix(Eij
+ Eij
.transpose())
1339 Sij_imag
= _embed_complex_matrix(I
*Eij
- I
*Eij
.transpose())
1344 def _quaternion_hermitian_basis(n
, field
=QQ
):
1346 Returns a basis for the space of quaternion Hermitian n-by-n matrices.
1350 sage: set_random_seed()
1351 sage: n = ZZ.random_element(1,5)
1352 sage: all( M.is_symmetric() for M in _quaternion_hermitian_basis(n) )
1356 Q
= QuaternionAlgebra(QQ
,-1,-1)
1359 # This is like the symmetric case, but we need to be careful:
1361 # * We want conjugate-symmetry, not just symmetry.
1362 # * The diagonal will (as a result) be real.
1366 for j
in xrange(i
+1):
1367 Eij
= matrix(Q
, n
, lambda k
,l
: k
==i
and l
==j
)
1369 Sij
= _embed_quaternion_matrix(Eij
)
1372 # Beware, orthogonal but not normalized! The second,
1373 # third, and fourth ones have a minus because they're
1375 Sij_real
= _embed_quaternion_matrix(Eij
+ Eij
.transpose())
1377 Sij_I
= _embed_quaternion_matrix(I
*Eij
- I
*Eij
.transpose())
1379 Sij_J
= _embed_quaternion_matrix(J
*Eij
- J
*Eij
.transpose())
1381 Sij_K
= _embed_quaternion_matrix(K
*Eij
- K
*Eij
.transpose())
1387 return vector(m
.base_ring(), m
.list())
1390 return matrix(v
.base_ring(), sqrt(v
.degree()), v
.list())
1392 def _multiplication_table_from_matrix_basis(basis
):
1394 At least three of the five simple Euclidean Jordan algebras have the
1395 symmetric multiplication (A,B) |-> (AB + BA)/2, where the
1396 multiplication on the right is matrix multiplication. Given a basis
1397 for the underlying matrix space, this function returns a
1398 multiplication table (obtained by looping through the basis
1399 elements) for an algebra of those matrices. A reordered copy
1400 of the basis is also returned to work around the fact that
1401 the ``span()`` in this function will change the order of the basis
1402 from what we think it is, to... something else.
1404 # In S^2, for example, we nominally have four coordinates even
1405 # though the space is of dimension three only. The vector space V
1406 # is supposed to hold the entire long vector, and the subspace W
1407 # of V will be spanned by the vectors that arise from symmetric
1408 # matrices. Thus for S^2, dim(V) == 4 and dim(W) == 3.
1409 field
= basis
[0].base_ring()
1410 dimension
= basis
[0].nrows()
1412 V
= VectorSpace(field
, dimension
**2)
1413 W
= V
.span( _mat2vec(s
) for s
in basis
)
1415 # Taking the span above reorders our basis (thanks, jerk!) so we
1416 # need to put our "matrix basis" in the same order as the
1417 # (reordered) vector basis.
1418 S
= tuple( _vec2mat(b
) for b
in W
.basis() )
1422 # Brute force the multiplication-by-s matrix by looping
1423 # through all elements of the basis and doing the computation
1424 # to find out what the corresponding row should be. BEWARE:
1425 # these multiplication tables won't be symmetric! It therefore
1426 # becomes REALLY IMPORTANT that the underlying algebra
1427 # constructor uses ROW vectors and not COLUMN vectors. That's
1428 # why we're computing rows here and not columns.
1431 this_row
= _mat2vec((s
*t
+ t
*s
)/2)
1432 Q_rows
.append(W
.coordinates(this_row
))
1433 Q
= matrix(field
, W
.dimension(), Q_rows
)
1439 def _embed_complex_matrix(M
):
1441 Embed the n-by-n complex matrix ``M`` into the space of real
1442 matrices of size 2n-by-2n via the map the sends each entry `z = a +
1443 bi` to the block matrix ``[[a,b],[-b,a]]``.
1447 sage: F = QuadraticField(-1,'i')
1448 sage: x1 = F(4 - 2*i)
1449 sage: x2 = F(1 + 2*i)
1452 sage: M = matrix(F,2,[[x1,x2],[x3,x4]])
1453 sage: _embed_complex_matrix(M)
1462 Embedding is a homomorphism (isomorphism, in fact)::
1464 sage: set_random_seed()
1465 sage: n = ZZ.random_element(5)
1466 sage: F = QuadraticField(-1, 'i')
1467 sage: X = random_matrix(F, n)
1468 sage: Y = random_matrix(F, n)
1469 sage: actual = _embed_complex_matrix(X) * _embed_complex_matrix(Y)
1470 sage: expected = _embed_complex_matrix(X*Y)
1471 sage: actual == expected
1477 raise ValueError("the matrix 'M' must be square")
1478 field
= M
.base_ring()
1483 blocks
.append(matrix(field
, 2, [[a
,b
],[-b
,a
]]))
1485 # We can drop the imaginaries here.
1486 return block_matrix(field
.base_ring(), n
, blocks
)
1489 def _unembed_complex_matrix(M
):
1491 The inverse of _embed_complex_matrix().
1495 sage: A = matrix(QQ,[ [ 1, 2, 3, 4],
1496 ....: [-2, 1, -4, 3],
1497 ....: [ 9, 10, 11, 12],
1498 ....: [-10, 9, -12, 11] ])
1499 sage: _unembed_complex_matrix(A)
1501 [ 10*i + 9 12*i + 11]
1505 Unembedding is the inverse of embedding::
1507 sage: set_random_seed()
1508 sage: F = QuadraticField(-1, 'i')
1509 sage: M = random_matrix(F, 3)
1510 sage: _unembed_complex_matrix(_embed_complex_matrix(M)) == M
1516 raise ValueError("the matrix 'M' must be square")
1517 if not n
.mod(2).is_zero():
1518 raise ValueError("the matrix 'M' must be a complex embedding")
1520 F
= QuadraticField(-1, 'i')
1523 # Go top-left to bottom-right (reading order), converting every
1524 # 2-by-2 block we see to a single complex element.
1526 for k
in xrange(n
/2):
1527 for j
in xrange(n
/2):
1528 submat
= M
[2*k
:2*k
+2,2*j
:2*j
+2]
1529 if submat
[0,0] != submat
[1,1]:
1530 raise ValueError('bad on-diagonal submatrix')
1531 if submat
[0,1] != -submat
[1,0]:
1532 raise ValueError('bad off-diagonal submatrix')
1533 z
= submat
[0,0] + submat
[0,1]*i
1536 return matrix(F
, n
/2, elements
)
1539 def _embed_quaternion_matrix(M
):
1541 Embed the n-by-n quaternion matrix ``M`` into the space of real
1542 matrices of size 4n-by-4n by first sending each quaternion entry
1543 `z = a + bi + cj + dk` to the block-complex matrix
1544 ``[[a + bi, c+di],[-c + di, a-bi]]`, and then embedding those into
1549 sage: Q = QuaternionAlgebra(QQ,-1,-1)
1550 sage: i,j,k = Q.gens()
1551 sage: x = 1 + 2*i + 3*j + 4*k
1552 sage: M = matrix(Q, 1, [[x]])
1553 sage: _embed_quaternion_matrix(M)
1559 Embedding is a homomorphism (isomorphism, in fact)::
1561 sage: set_random_seed()
1562 sage: n = ZZ.random_element(5)
1563 sage: Q = QuaternionAlgebra(QQ,-1,-1)
1564 sage: X = random_matrix(Q, n)
1565 sage: Y = random_matrix(Q, n)
1566 sage: actual = _embed_quaternion_matrix(X)*_embed_quaternion_matrix(Y)
1567 sage: expected = _embed_quaternion_matrix(X*Y)
1568 sage: actual == expected
1572 quaternions
= M
.base_ring()
1575 raise ValueError("the matrix 'M' must be square")
1577 F
= QuadraticField(-1, 'i')
1582 t
= z
.coefficient_tuple()
1587 cplx_matrix
= matrix(F
, 2, [[ a
+ b
*i
, c
+ d
*i
],
1588 [-c
+ d
*i
, a
- b
*i
]])
1589 blocks
.append(_embed_complex_matrix(cplx_matrix
))
1591 # We should have real entries by now, so use the realest field
1592 # we've got for the return value.
1593 return block_matrix(quaternions
.base_ring(), n
, blocks
)
1596 def _unembed_quaternion_matrix(M
):
1598 The inverse of _embed_quaternion_matrix().
1602 sage: M = matrix(QQ, [[ 1, 2, 3, 4],
1603 ....: [-2, 1, -4, 3],
1604 ....: [-3, 4, 1, -2],
1605 ....: [-4, -3, 2, 1]])
1606 sage: _unembed_quaternion_matrix(M)
1607 [1 + 2*i + 3*j + 4*k]
1611 Unembedding is the inverse of embedding::
1613 sage: set_random_seed()
1614 sage: Q = QuaternionAlgebra(QQ, -1, -1)
1615 sage: M = random_matrix(Q, 3)
1616 sage: _unembed_quaternion_matrix(_embed_quaternion_matrix(M)) == M
1622 raise ValueError("the matrix 'M' must be square")
1623 if not n
.mod(4).is_zero():
1624 raise ValueError("the matrix 'M' must be a complex embedding")
1626 Q
= QuaternionAlgebra(QQ
,-1,-1)
1629 # Go top-left to bottom-right (reading order), converting every
1630 # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
1633 for l
in xrange(n
/4):
1634 for m
in xrange(n
/4):
1635 submat
= _unembed_complex_matrix(M
[4*l
:4*l
+4,4*m
:4*m
+4])
1636 if submat
[0,0] != submat
[1,1].conjugate():
1637 raise ValueError('bad on-diagonal submatrix')
1638 if submat
[0,1] != -submat
[1,0].conjugate():
1639 raise ValueError('bad off-diagonal submatrix')
1640 z
= submat
[0,0].real() + submat
[0,0].imag()*i
1641 z
+= submat
[0,1].real()*j
+ submat
[0,1].imag()*k
1644 return matrix(Q
, n
/4, elements
)
1647 # The usual inner product on R^n.
1649 return x
.vector().inner_product(y
.vector())
1651 # The inner product used for the real symmetric simple EJA.
1652 # We keep it as a separate function because e.g. the complex
1653 # algebra uses the same inner product, except divided by 2.
1654 def _matrix_ip(X
,Y
):
1655 X_mat
= X
.natural_representation()
1656 Y_mat
= Y
.natural_representation()
1657 return (X_mat
*Y_mat
).trace()
1660 class RealSymmetricEJA(FiniteDimensionalEuclideanJordanAlgebra
):
1662 The rank-n simple EJA consisting of real symmetric n-by-n
1663 matrices, the usual symmetric Jordan product, and the trace inner
1664 product. It has dimension `(n^2 + n)/2` over the reals.
1668 sage: J = RealSymmetricEJA(2)
1669 sage: e0, e1, e2 = J.gens()
1679 The degree of this algebra is `(n^2 + n) / 2`::
1681 sage: set_random_seed()
1682 sage: n = ZZ.random_element(1,5)
1683 sage: J = RealSymmetricEJA(n)
1684 sage: J.degree() == (n^2 + n)/2
1687 The Jordan multiplication is what we think it is::
1689 sage: set_random_seed()
1690 sage: n = ZZ.random_element(1,5)
1691 sage: J = RealSymmetricEJA(n)
1692 sage: x = J.random_element()
1693 sage: y = J.random_element()
1694 sage: actual = (x*y).natural_representation()
1695 sage: X = x.natural_representation()
1696 sage: Y = y.natural_representation()
1697 sage: expected = (X*Y + Y*X)/2
1698 sage: actual == expected
1700 sage: J(expected) == x*y
1705 def __classcall_private__(cls
, n
, field
=QQ
):
1706 S
= _real_symmetric_basis(n
, field
=field
)
1707 (Qs
, T
) = _multiplication_table_from_matrix_basis(S
)
1709 fdeja
= super(RealSymmetricEJA
, cls
)
1710 return fdeja
.__classcall
_private
__(cls
,
1716 def inner_product(self
, x
, y
):
1717 return _matrix_ip(x
,y
)
1720 class ComplexHermitianEJA(FiniteDimensionalEuclideanJordanAlgebra
):
1722 The rank-n simple EJA consisting of complex Hermitian n-by-n
1723 matrices over the real numbers, the usual symmetric Jordan product,
1724 and the real-part-of-trace inner product. It has dimension `n^2` over
1729 The degree of this algebra is `n^2`::
1731 sage: set_random_seed()
1732 sage: n = ZZ.random_element(1,5)
1733 sage: J = ComplexHermitianEJA(n)
1734 sage: J.degree() == n^2
1737 The Jordan multiplication is what we think it is::
1739 sage: set_random_seed()
1740 sage: n = ZZ.random_element(1,5)
1741 sage: J = ComplexHermitianEJA(n)
1742 sage: x = J.random_element()
1743 sage: y = J.random_element()
1744 sage: actual = (x*y).natural_representation()
1745 sage: X = x.natural_representation()
1746 sage: Y = y.natural_representation()
1747 sage: expected = (X*Y + Y*X)/2
1748 sage: actual == expected
1750 sage: J(expected) == x*y
1755 def __classcall_private__(cls
, n
, field
=QQ
):
1756 S
= _complex_hermitian_basis(n
)
1757 (Qs
, T
) = _multiplication_table_from_matrix_basis(S
)
1759 fdeja
= super(ComplexHermitianEJA
, cls
)
1760 return fdeja
.__classcall
_private
__(cls
,
1766 def inner_product(self
, x
, y
):
1767 # Since a+bi on the diagonal is represented as
1772 # we'll double-count the "a" entries if we take the trace of
1774 return _matrix_ip(x
,y
)/2
1777 class QuaternionHermitianEJA(FiniteDimensionalEuclideanJordanAlgebra
):
1779 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
1780 matrices, the usual symmetric Jordan product, and the
1781 real-part-of-trace inner product. It has dimension `2n^2 - n` over
1786 The degree of this algebra is `n^2`::
1788 sage: set_random_seed()
1789 sage: n = ZZ.random_element(1,5)
1790 sage: J = QuaternionHermitianEJA(n)
1791 sage: J.degree() == 2*(n^2) - n
1794 The Jordan multiplication is what we think it is::
1796 sage: set_random_seed()
1797 sage: n = ZZ.random_element(1,5)
1798 sage: J = QuaternionHermitianEJA(n)
1799 sage: x = J.random_element()
1800 sage: y = J.random_element()
1801 sage: actual = (x*y).natural_representation()
1802 sage: X = x.natural_representation()
1803 sage: Y = y.natural_representation()
1804 sage: expected = (X*Y + Y*X)/2
1805 sage: actual == expected
1807 sage: J(expected) == x*y
1812 def __classcall_private__(cls
, n
, field
=QQ
):
1813 S
= _quaternion_hermitian_basis(n
)
1814 (Qs
, T
) = _multiplication_table_from_matrix_basis(S
)
1816 fdeja
= super(QuaternionHermitianEJA
, cls
)
1817 return fdeja
.__classcall
_private
__(cls
,
1823 def inner_product(self
, x
, y
):
1824 # Since a+bi+cj+dk on the diagonal is represented as
1826 # a + bi +cj + dk = [ a b c d]
1831 # we'll quadruple-count the "a" entries if we take the trace of
1833 return _matrix_ip(x
,y
)/4
1836 class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra
):
1838 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
1839 with the usual inner product and jordan product ``x*y =
1840 (<x_bar,y_bar>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
1845 This multiplication table can be verified by hand::
1847 sage: J = JordanSpinEJA(4)
1848 sage: e0,e1,e2,e3 = J.gens()
1866 def __classcall_private__(cls
, n
, field
=QQ
):
1868 id_matrix
= identity_matrix(field
, n
)
1870 ei
= id_matrix
.column(i
)
1871 Qi
= zero_matrix(field
, n
)
1873 Qi
.set_column(0, ei
)
1874 Qi
+= diagonal_matrix(n
, [ei
[0]]*n
)
1875 # The addition of the diagonal matrix adds an extra ei[0] in the
1876 # upper-left corner of the matrix.
1877 Qi
[0,0] = Qi
[0,0] * ~
field(2)
1880 # The rank of the spin algebra is two, unless we're in a
1881 # one-dimensional ambient space (because the rank is bounded by
1882 # the ambient dimension).
1883 fdeja
= super(JordanSpinEJA
, cls
)
1884 return fdeja
.__classcall
_private
__(cls
, field
, Qs
, rank
=min(n
,2))
1886 def inner_product(self
, x
, y
):
1887 return _usual_ip(x
,y
)