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,
27 mult_table
= [b
.base_extend(field
) for b
in mult_table
]
30 if not (is_Matrix(b
) and b
.dimensions() == (n
, n
)):
31 raise ValueError("input is not a multiplication table")
32 mult_table
= tuple(mult_table
)
34 cat
= MagmaticAlgebras(field
).FiniteDimensional().WithBasis()
35 cat
.or_subcategory(category
)
36 if assume_associative
:
37 cat
= cat
.Associative()
39 names
= normalize_names(n
, names
)
41 fda
= super(FiniteDimensionalEuclideanJordanAlgebra
, cls
)
42 return fda
.__classcall
__(cls
,
45 assume_associative
=assume_associative
,
49 natural_basis
=natural_basis
,
50 inner_product
=inner_product
)
53 def __init__(self
, field
,
56 assume_associative
=False,
64 By definition, Jordan multiplication commutes::
66 sage: set_random_seed()
67 sage: J = random_eja()
68 sage: x = J.random_element()
69 sage: y = J.random_element()
75 self
._natural
_basis
= natural_basis
76 self
._inner
_product
= inner_product
77 fda
= super(FiniteDimensionalEuclideanJordanAlgebra
, self
)
86 Return a string representation of ``self``.
88 fmt
= "Euclidean Jordan algebra of degree {} over {}"
89 return fmt
.format(self
.degree(), self
.base_ring())
92 def inner_product(self
, x
, y
):
94 The inner product associated with this Euclidean Jordan algebra.
96 Will default to the trace inner product if nothing else.
98 if (not x
in self
) or (not y
in self
):
99 raise ArgumentError("arguments must live in this algebra")
100 if self
._inner
_product
is None:
101 return x
.trace_inner_product(y
)
103 return self
._inner
_product
(x
,y
)
106 def natural_basis(self
):
108 Return a more-natural representation of this algebra's basis.
110 Every finite-dimensional Euclidean Jordan Algebra is a direct
111 sum of five simple algebras, four of which comprise Hermitian
112 matrices. This method returns the original "natural" basis
113 for our underlying vector space. (Typically, the natural basis
114 is used to construct the multiplication table in the first place.)
116 Note that this will always return a matrix. The standard basis
117 in `R^n` will be returned as `n`-by-`1` column matrices.
121 sage: J = RealSymmetricSimpleEJA(2)
124 sage: J.natural_basis()
132 sage: J = JordanSpinSimpleEJA(2)
135 sage: J.natural_basis()
142 if self
._natural
_basis
is None:
143 return tuple( b
.vector().column() for b
in self
.basis() )
145 return self
._natural
_basis
150 Return the rank of this EJA.
152 if self
._rank
is None:
153 raise ValueError("no rank specified at genesis")
158 class Element(FiniteDimensionalAlgebraElement
):
160 An element of a Euclidean Jordan algebra.
163 def __pow__(self
, n
):
165 Return ``self`` raised to the power ``n``.
167 Jordan algebras are always power-associative; see for
168 example Faraut and Koranyi, Proposition II.1.2 (ii).
172 We have to override this because our superclass uses row vectors
173 instead of column vectors! We, on the other hand, assume column
178 sage: set_random_seed()
179 sage: x = random_eja().random_element()
180 sage: x.operator_matrix()*x.vector() == (x^2).vector()
183 A few examples of power-associativity::
185 sage: set_random_seed()
186 sage: x = random_eja().random_element()
187 sage: x*(x*x)*(x*x) == x^5
189 sage: (x*x)*(x*x*x) == x^5
192 We also know that powers operator-commute (Koecher, Chapter
195 sage: set_random_seed()
196 sage: x = random_eja().random_element()
197 sage: m = ZZ.random_element(0,10)
198 sage: n = ZZ.random_element(0,10)
199 sage: Lxm = (x^m).operator_matrix()
200 sage: Lxn = (x^n).operator_matrix()
201 sage: Lxm*Lxn == Lxn*Lxm
211 return A( (self
.operator_matrix()**(n
-1))*self
.vector() )
214 def characteristic_polynomial(self
):
216 Return my characteristic polynomial (if I'm a regular
219 Eventually this should be implemented in terms of the parent
220 algebra's characteristic polynomial that works for ALL
223 if self
.is_regular():
224 return self
.minimal_polynomial()
226 raise NotImplementedError('irregular element')
229 def inner_product(self
, other
):
231 Return the parent algebra's inner product of myself and ``other``.
235 The inner product in the Jordan spin algebra is the usual
236 inner product on `R^n` (this example only works because the
237 basis for the Jordan algebra is the standard basis in `R^n`)::
239 sage: J = JordanSpinSimpleEJA(3)
240 sage: x = vector(QQ,[1,2,3])
241 sage: y = vector(QQ,[4,5,6])
242 sage: x.inner_product(y)
244 sage: J(x).inner_product(J(y))
249 Ensure that we can always compute an inner product, and that
250 it gives us back a real number::
252 sage: set_random_seed()
253 sage: J = random_eja()
254 sage: x = J.random_element()
255 sage: y = J.random_element()
256 sage: x.inner_product(y) in RR
262 raise ArgumentError("'other' must live in the same algebra")
264 return P
.inner_product(self
, other
)
267 def operator_commutes_with(self
, other
):
269 Return whether or not this element operator-commutes
274 The definition of a Jordan algebra says that any element
275 operator-commutes with its square::
277 sage: set_random_seed()
278 sage: x = random_eja().random_element()
279 sage: x.operator_commutes_with(x^2)
284 Test Lemma 1 from Chapter III of Koecher::
286 sage: set_random_seed()
287 sage: J = random_eja()
288 sage: u = J.random_element()
289 sage: v = J.random_element()
290 sage: lhs = u.operator_commutes_with(u*v)
291 sage: rhs = v.operator_commutes_with(u^2)
296 if not other
in self
.parent():
297 raise ArgumentError("'other' must live in the same algebra")
299 A
= self
.operator_matrix()
300 B
= other
.operator_matrix()
306 Return my determinant, the product of my eigenvalues.
310 sage: J = JordanSpinSimpleEJA(2)
311 sage: e0,e1 = J.gens()
315 sage: J = JordanSpinSimpleEJA(3)
316 sage: e0,e1,e2 = J.gens()
317 sage: x = e0 + e1 + e2
322 cs
= self
.characteristic_polynomial().coefficients(sparse
=False)
325 return cs
[0] * (-1)**r
327 raise ValueError('charpoly had no coefficients')
332 Return the Jordan-multiplicative inverse of this element.
334 We can't use the superclass method because it relies on the
335 algebra being associative.
339 The inverse in the spin factor algebra is given in Alizadeh's
342 sage: set_random_seed()
343 sage: n = ZZ.random_element(1,10)
344 sage: J = JordanSpinSimpleEJA(n)
345 sage: x = J.random_element()
346 sage: while x.is_zero():
347 ....: x = J.random_element()
348 sage: x_vec = x.vector()
350 sage: x_bar = x_vec[1:]
351 sage: coeff = 1/(x0^2 - x_bar.inner_product(x_bar))
352 sage: inv_vec = x_vec.parent()([x0] + (-x_bar).list())
353 sage: x_inverse = coeff*inv_vec
354 sage: x.inverse() == J(x_inverse)
359 The identity element is its own inverse::
361 sage: set_random_seed()
362 sage: J = random_eja()
363 sage: J.one().inverse() == J.one()
366 If an element has an inverse, it acts like one. TODO: this
367 can be a lot less ugly once ``is_invertible`` doesn't crash
368 on irregular elements::
370 sage: set_random_seed()
371 sage: J = random_eja()
372 sage: x = J.random_element()
374 ....: x.inverse()*x == J.one()
380 if self
.parent().is_associative():
381 elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
384 # TODO: we can do better once the call to is_invertible()
385 # doesn't crash on irregular elements.
386 #if not self.is_invertible():
387 # raise ArgumentError('element is not invertible')
389 # We do this a little different than the usual recursive
390 # call to a finite-dimensional algebra element, because we
391 # wind up with an inverse that lives in the subalgebra and
392 # we need information about the parent to convert it back.
393 V
= self
.span_of_powers()
394 assoc_subalg
= self
.subalgebra_generated_by()
395 # Mis-design warning: the basis used for span_of_powers()
396 # and subalgebra_generated_by() must be the same, and in
398 elt
= assoc_subalg(V
.coordinates(self
.vector()))
400 # This will be in the subalgebra's coordinates...
401 fda_elt
= FiniteDimensionalAlgebraElement(assoc_subalg
, elt
)
402 subalg_inverse
= fda_elt
.inverse()
404 # So we have to convert back...
405 basis
= [ self
.parent(v
) for v
in V
.basis() ]
406 pairs
= zip(subalg_inverse
.vector(), basis
)
407 return self
.parent().linear_combination(pairs
)
410 def is_invertible(self
):
412 Return whether or not this element is invertible.
414 We can't use the superclass method because it relies on
415 the algebra being associative.
417 return not self
.det().is_zero()
420 def is_nilpotent(self
):
422 Return whether or not some power of this element is zero.
424 The superclass method won't work unless we're in an
425 associative algebra, and we aren't. However, we generate
426 an assocoative subalgebra and we're nilpotent there if and
427 only if we're nilpotent here (probably).
431 The identity element is never nilpotent::
433 sage: set_random_seed()
434 sage: random_eja().one().is_nilpotent()
437 The additive identity is always nilpotent::
439 sage: set_random_seed()
440 sage: random_eja().zero().is_nilpotent()
444 # The element we're going to call "is_nilpotent()" on.
445 # Either myself, interpreted as an element of a finite-
446 # dimensional algebra, or an element of an associative
450 if self
.parent().is_associative():
451 elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
453 V
= self
.span_of_powers()
454 assoc_subalg
= self
.subalgebra_generated_by()
455 # Mis-design warning: the basis used for span_of_powers()
456 # and subalgebra_generated_by() must be the same, and in
458 elt
= assoc_subalg(V
.coordinates(self
.vector()))
460 # Recursive call, but should work since elt lives in an
461 # associative algebra.
462 return elt
.is_nilpotent()
465 def is_regular(self
):
467 Return whether or not this is a regular element.
471 The identity element always has degree one, but any element
472 linearly-independent from it is regular::
474 sage: J = JordanSpinSimpleEJA(5)
475 sage: J.one().is_regular()
477 sage: e0, e1, e2, e3, e4 = J.gens() # e0 is the identity
478 sage: for x in J.gens():
479 ....: (J.one() + x).is_regular()
487 return self
.degree() == self
.parent().rank()
492 Compute the degree of this element the straightforward way
493 according to the definition; by appending powers to a list
494 and figuring out its dimension (that is, whether or not
495 they're linearly dependent).
499 sage: J = JordanSpinSimpleEJA(4)
500 sage: J.one().degree()
502 sage: e0,e1,e2,e3 = J.gens()
503 sage: (e0 - e1).degree()
506 In the spin factor algebra (of rank two), all elements that
507 aren't multiples of the identity are regular::
509 sage: set_random_seed()
510 sage: n = ZZ.random_element(1,10)
511 sage: J = JordanSpinSimpleEJA(n)
512 sage: x = J.random_element()
513 sage: x == x.coefficient(0)*J.one() or x.degree() == 2
517 return self
.span_of_powers().dimension()
520 def minimal_polynomial(self
):
524 sage: set_random_seed()
525 sage: x = random_eja().random_element()
526 sage: x.degree() == x.minimal_polynomial().degree()
531 sage: set_random_seed()
532 sage: x = random_eja().random_element()
533 sage: x.degree() == x.minimal_polynomial().degree()
536 The minimal polynomial and the characteristic polynomial coincide
537 and are known (see Alizadeh, Example 11.11) for all elements of
538 the spin factor algebra that aren't scalar multiples of the
541 sage: set_random_seed()
542 sage: n = ZZ.random_element(2,10)
543 sage: J = JordanSpinSimpleEJA(n)
544 sage: y = J.random_element()
545 sage: while y == y.coefficient(0)*J.one():
546 ....: y = J.random_element()
547 sage: y0 = y.vector()[0]
548 sage: y_bar = y.vector()[1:]
549 sage: actual = y.minimal_polynomial()
550 sage: x = SR.symbol('x', domain='real')
551 sage: expected = x^2 - 2*y0*x + (y0^2 - norm(y_bar)^2)
552 sage: bool(actual == expected)
556 # The element we're going to call "minimal_polynomial()" on.
557 # Either myself, interpreted as an element of a finite-
558 # dimensional algebra, or an element of an associative
562 if self
.parent().is_associative():
563 elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
565 V
= self
.span_of_powers()
566 assoc_subalg
= self
.subalgebra_generated_by()
567 # Mis-design warning: the basis used for span_of_powers()
568 # and subalgebra_generated_by() must be the same, and in
570 elt
= assoc_subalg(V
.coordinates(self
.vector()))
572 # Recursive call, but should work since elt lives in an
573 # associative algebra.
574 return elt
.minimal_polynomial()
577 def natural_representation(self
):
579 Return a more-natural representation of this element.
581 Every finite-dimensional Euclidean Jordan Algebra is a
582 direct sum of five simple algebras, four of which comprise
583 Hermitian matrices. This method returns the original
584 "natural" representation of this element as a Hermitian
585 matrix, if it has one. If not, you get the usual representation.
589 sage: J = ComplexHermitianSimpleEJA(3)
592 sage: J.one().natural_representation()
601 B
= self
.parent().natural_basis()
602 W
= B
[0].matrix_space()
603 return W
.linear_combination(zip(self
.vector(), B
))
606 def operator_matrix(self
):
608 Return the matrix that represents left- (or right-)
609 multiplication by this element in the parent algebra.
611 We have to override this because the superclass method
612 returns a matrix that acts on row vectors (that is, on
617 Test the first polarization identity from my notes, Koecher Chapter
618 III, or from Baes (2.3)::
620 sage: set_random_seed()
621 sage: J = random_eja()
622 sage: x = J.random_element()
623 sage: y = J.random_element()
624 sage: Lx = x.operator_matrix()
625 sage: Ly = y.operator_matrix()
626 sage: Lxx = (x*x).operator_matrix()
627 sage: Lxy = (x*y).operator_matrix()
628 sage: bool(2*Lx*Lxy + Ly*Lxx == 2*Lxy*Lx + Lxx*Ly)
631 Test the second polarization identity from my notes or from
634 sage: set_random_seed()
635 sage: J = random_eja()
636 sage: x = J.random_element()
637 sage: y = J.random_element()
638 sage: z = J.random_element()
639 sage: Lx = x.operator_matrix()
640 sage: Ly = y.operator_matrix()
641 sage: Lz = z.operator_matrix()
642 sage: Lzy = (z*y).operator_matrix()
643 sage: Lxy = (x*y).operator_matrix()
644 sage: Lxz = (x*z).operator_matrix()
645 sage: bool(Lx*Lzy + Lz*Lxy + Ly*Lxz == Lzy*Lx + Lxy*Lz + Lxz*Ly)
648 Test the third polarization identity from my notes or from
651 sage: set_random_seed()
652 sage: J = random_eja()
653 sage: u = J.random_element()
654 sage: y = J.random_element()
655 sage: z = J.random_element()
656 sage: Lu = u.operator_matrix()
657 sage: Ly = y.operator_matrix()
658 sage: Lz = z.operator_matrix()
659 sage: Lzy = (z*y).operator_matrix()
660 sage: Luy = (u*y).operator_matrix()
661 sage: Luz = (u*z).operator_matrix()
662 sage: Luyz = (u*(y*z)).operator_matrix()
663 sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz
664 sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly
665 sage: bool(lhs == rhs)
669 fda_elt
= FiniteDimensionalAlgebraElement(self
.parent(), self
)
670 return fda_elt
.matrix().transpose()
673 def quadratic_representation(self
, other
=None):
675 Return the quadratic representation of this element.
679 The explicit form in the spin factor algebra is given by
680 Alizadeh's Example 11.12::
682 sage: set_random_seed()
683 sage: n = ZZ.random_element(1,10)
684 sage: J = JordanSpinSimpleEJA(n)
685 sage: x = J.random_element()
686 sage: x_vec = x.vector()
688 sage: x_bar = x_vec[1:]
689 sage: A = matrix(QQ, 1, [x_vec.inner_product(x_vec)])
690 sage: B = 2*x0*x_bar.row()
691 sage: C = 2*x0*x_bar.column()
692 sage: D = identity_matrix(QQ, n-1)
693 sage: D = (x0^2 - x_bar.inner_product(x_bar))*D
694 sage: D = D + 2*x_bar.tensor_product(x_bar)
695 sage: Q = block_matrix(2,2,[A,B,C,D])
696 sage: Q == x.quadratic_representation()
699 Test all of the properties from Theorem 11.2 in Alizadeh::
701 sage: set_random_seed()
702 sage: J = random_eja()
703 sage: x = J.random_element()
704 sage: y = J.random_element()
708 sage: actual = x.quadratic_representation(y)
709 sage: expected = ( (x+y).quadratic_representation()
710 ....: -x.quadratic_representation()
711 ....: -y.quadratic_representation() ) / 2
712 sage: actual == expected
717 sage: alpha = QQ.random_element()
718 sage: actual = (alpha*x).quadratic_representation()
719 sage: expected = (alpha^2)*x.quadratic_representation()
720 sage: actual == expected
725 sage: Qy = y.quadratic_representation()
726 sage: actual = J(Qy*x.vector()).quadratic_representation()
727 sage: expected = Qy*x.quadratic_representation()*Qy
728 sage: actual == expected
733 sage: k = ZZ.random_element(1,10)
734 sage: actual = (x^k).quadratic_representation()
735 sage: expected = (x.quadratic_representation())^k
736 sage: actual == expected
742 elif not other
in self
.parent():
743 raise ArgumentError("'other' must live in the same algebra")
745 L
= self
.operator_matrix()
746 M
= other
.operator_matrix()
747 return ( L
*M
+ M
*L
- (self
*other
).operator_matrix() )
750 def span_of_powers(self
):
752 Return the vector space spanned by successive powers of
755 # The dimension of the subalgebra can't be greater than
756 # the big algebra, so just put everything into a list
757 # and let span() get rid of the excess.
758 V
= self
.vector().parent()
759 return V
.span( (self
**d
).vector() for d
in xrange(V
.dimension()) )
762 def subalgebra_generated_by(self
):
764 Return the associative subalgebra of the parent EJA generated
769 sage: set_random_seed()
770 sage: x = random_eja().random_element()
771 sage: x.subalgebra_generated_by().is_associative()
774 Squaring in the subalgebra should be the same thing as
775 squaring in the superalgebra::
777 sage: set_random_seed()
778 sage: x = random_eja().random_element()
779 sage: u = x.subalgebra_generated_by().random_element()
780 sage: u.operator_matrix()*u.vector() == (u**2).vector()
784 # First get the subspace spanned by the powers of myself...
785 V
= self
.span_of_powers()
788 # Now figure out the entries of the right-multiplication
789 # matrix for the successive basis elements b0, b1,... of
792 for b_right
in V
.basis():
793 eja_b_right
= self
.parent()(b_right
)
795 # The first row of the right-multiplication matrix by
796 # b1 is what we get if we apply that matrix to b1. The
797 # second row of the right multiplication matrix by b1
798 # is what we get when we apply that matrix to b2...
800 # IMPORTANT: this assumes that all vectors are COLUMN
801 # vectors, unlike our superclass (which uses row vectors).
802 for b_left
in V
.basis():
803 eja_b_left
= self
.parent()(b_left
)
804 # Multiply in the original EJA, but then get the
805 # coordinates from the subalgebra in terms of its
807 this_row
= V
.coordinates((eja_b_left
*eja_b_right
).vector())
808 b_right_rows
.append(this_row
)
809 b_right_matrix
= matrix(F
, b_right_rows
)
810 mats
.append(b_right_matrix
)
812 # It's an algebra of polynomials in one element, and EJAs
813 # are power-associative.
815 # TODO: choose generator names intelligently.
816 return FiniteDimensionalEuclideanJordanAlgebra(F
, mats
, assume_associative
=True, names
='f')
819 def subalgebra_idempotent(self
):
821 Find an idempotent in the associative subalgebra I generate
822 using Proposition 2.3.5 in Baes.
826 sage: set_random_seed()
828 sage: c = J.random_element().subalgebra_idempotent()
831 sage: J = JordanSpinSimpleEJA(5)
832 sage: c = J.random_element().subalgebra_idempotent()
837 if self
.is_nilpotent():
838 raise ValueError("this only works with non-nilpotent elements!")
840 V
= self
.span_of_powers()
841 J
= self
.subalgebra_generated_by()
842 # Mis-design warning: the basis used for span_of_powers()
843 # and subalgebra_generated_by() must be the same, and in
845 u
= J(V
.coordinates(self
.vector()))
847 # The image of the matrix of left-u^m-multiplication
848 # will be minimal for some natural number s...
850 minimal_dim
= V
.dimension()
851 for i
in xrange(1, V
.dimension()):
852 this_dim
= (u
**i
).operator_matrix().image().dimension()
853 if this_dim
< minimal_dim
:
854 minimal_dim
= this_dim
857 # Now minimal_matrix should correspond to the smallest
858 # non-zero subspace in Baes's (or really, Koecher's)
861 # However, we need to restrict the matrix to work on the
862 # subspace... or do we? Can't we just solve, knowing that
863 # A(c) = u^(s+1) should have a solution in the big space,
866 # Beware, solve_right() means that we're using COLUMN vectors.
867 # Our FiniteDimensionalAlgebraElement superclass uses rows.
869 A
= u_next
.operator_matrix()
870 c_coordinates
= A
.solve_right(u_next
.vector())
872 # Now c_coordinates is the idempotent we want, but it's in
873 # the coordinate system of the subalgebra.
875 # We need the basis for J, but as elements of the parent algebra.
877 basis
= [self
.parent(v
) for v
in V
.basis()]
878 return self
.parent().linear_combination(zip(c_coordinates
, basis
))
883 Return my trace, the sum of my eigenvalues.
887 sage: J = JordanSpinSimpleEJA(3)
888 sage: e0,e1,e2 = J.gens()
889 sage: x = e0 + e1 + e2
894 cs
= self
.characteristic_polynomial().coefficients(sparse
=False)
898 raise ValueError('charpoly had fewer than 2 coefficients')
901 def trace_inner_product(self
, other
):
903 Return the trace inner product of myself and ``other``.
905 if not other
in self
.parent():
906 raise ArgumentError("'other' must live in the same algebra")
908 return (self
*other
).trace()
911 def eja_rn(dimension
, field
=QQ
):
913 Return the Euclidean Jordan Algebra corresponding to the set
914 `R^n` under the Hadamard product.
918 This multiplication table can be verified by hand::
921 sage: e0,e1,e2 = J.gens()
936 # The FiniteDimensionalAlgebra constructor takes a list of
937 # matrices, the ith representing right multiplication by the ith
938 # basis element in the vector space. So if e_1 = (1,0,0), then
939 # right (Hadamard) multiplication of x by e_1 picks out the first
940 # component of x; and likewise for the ith basis element e_i.
941 Qs
= [ matrix(field
, dimension
, dimension
, lambda k
,j
: 1*(k
== j
== i
))
942 for i
in xrange(dimension
) ]
944 # The usual inner product on R^n.
945 ip
= lambda x
, y
: x
.vector().inner_product(y
.vector())
947 return FiniteDimensionalEuclideanJordanAlgebra(field
,
956 Return a "random" finite-dimensional Euclidean Jordan Algebra.
960 For now, we choose a random natural number ``n`` (greater than zero)
961 and then give you back one of the following:
963 * The cartesian product of the rational numbers ``n`` times; this is
964 ``QQ^n`` with the Hadamard product.
966 * The Jordan spin algebra on ``QQ^n``.
968 * The ``n``-by-``n`` rational symmetric matrices with the symmetric
971 Later this might be extended to return Cartesian products of the
977 Euclidean Jordan algebra of degree...
980 n
= ZZ
.random_element(1,5)
981 constructor
= choice([eja_rn
,
983 RealSymmetricSimpleEJA
,
984 ComplexHermitianSimpleEJA
])
985 return constructor(n
, field
=QQ
)
989 def _real_symmetric_basis(n
, field
=QQ
):
991 Return a basis for the space of real symmetric n-by-n matrices.
993 # The basis of symmetric matrices, as matrices, in their R^(n-by-n)
997 for j
in xrange(i
+1):
998 Eij
= matrix(field
, n
, lambda k
,l
: k
==i
and l
==j
)
1002 # Beware, orthogonal but not normalized!
1003 Sij
= Eij
+ Eij
.transpose()
1008 def _complex_hermitian_basis(n
, field
=QQ
):
1010 Returns a basis for the space of complex Hermitian n-by-n matrices.
1014 sage: set_random_seed()
1015 sage: n = ZZ.random_element(1,5)
1016 sage: all( M.is_symmetric() for M in _complex_hermitian_basis(n) )
1020 F
= QuadraticField(-1, 'I')
1023 # This is like the symmetric case, but we need to be careful:
1025 # * We want conjugate-symmetry, not just symmetry.
1026 # * The diagonal will (as a result) be real.
1030 for j
in xrange(i
+1):
1031 Eij
= matrix(field
, n
, lambda k
,l
: k
==i
and l
==j
)
1033 Sij
= _embed_complex_matrix(Eij
)
1036 # Beware, orthogonal but not normalized! The second one
1037 # has a minus because it's conjugated.
1038 Sij_real
= _embed_complex_matrix(Eij
+ Eij
.transpose())
1040 Sij_imag
= _embed_complex_matrix(I
*Eij
- I
*Eij
.transpose())
1045 def _multiplication_table_from_matrix_basis(basis
):
1047 At least three of the five simple Euclidean Jordan algebras have the
1048 symmetric multiplication (A,B) |-> (AB + BA)/2, where the
1049 multiplication on the right is matrix multiplication. Given a basis
1050 for the underlying matrix space, this function returns a
1051 multiplication table (obtained by looping through the basis
1052 elements) for an algebra of those matrices. A reordered copy
1053 of the basis is also returned to work around the fact that
1054 the ``span()`` in this function will change the order of the basis
1055 from what we think it is, to... something else.
1057 # In S^2, for example, we nominally have four coordinates even
1058 # though the space is of dimension three only. The vector space V
1059 # is supposed to hold the entire long vector, and the subspace W
1060 # of V will be spanned by the vectors that arise from symmetric
1061 # matrices. Thus for S^2, dim(V) == 4 and dim(W) == 3.
1062 field
= basis
[0].base_ring()
1063 dimension
= basis
[0].nrows()
1066 return vector(field
, m
.list())
1069 return matrix(field
, dimension
, v
.list())
1071 V
= VectorSpace(field
, dimension
**2)
1072 W
= V
.span( mat2vec(s
) for s
in basis
)
1074 # Taking the span above reorders our basis (thanks, jerk!) so we
1075 # need to put our "matrix basis" in the same order as the
1076 # (reordered) vector basis.
1077 S
= tuple( vec2mat(b
) for b
in W
.basis() )
1081 # Brute force the multiplication-by-s matrix by looping
1082 # through all elements of the basis and doing the computation
1083 # to find out what the corresponding row should be. BEWARE:
1084 # these multiplication tables won't be symmetric! It therefore
1085 # becomes REALLY IMPORTANT that the underlying algebra
1086 # constructor uses ROW vectors and not COLUMN vectors. That's
1087 # why we're computing rows here and not columns.
1090 this_row
= mat2vec((s
*t
+ t
*s
)/2)
1091 Q_rows
.append(W
.coordinates(this_row
))
1092 Q
= matrix(field
, W
.dimension(), Q_rows
)
1098 def _embed_complex_matrix(M
):
1100 Embed the n-by-n complex matrix ``M`` into the space of real
1101 matrices of size 2n-by-2n via the map the sends each entry `z = a +
1102 bi` to the block matrix ``[[a,b],[-b,a]]``.
1106 sage: F = QuadraticField(-1,'i')
1107 sage: x1 = F(4 - 2*i)
1108 sage: x2 = F(1 + 2*i)
1111 sage: M = matrix(F,2,[x1,x2,x3,x4])
1112 sage: _embed_complex_matrix(M)
1122 raise ArgumentError("the matrix 'M' must be square")
1123 field
= M
.base_ring()
1128 blocks
.append(matrix(field
, 2, [[a
,-b
],[b
,a
]]))
1130 # We can drop the imaginaries here.
1131 return block_matrix(field
.base_ring(), n
, blocks
)
1134 def _unembed_complex_matrix(M
):
1136 The inverse of _embed_complex_matrix().
1140 sage: A = matrix(QQ,[ [ 1, 2, 3, 4],
1141 ....: [-2, 1, -4, 3],
1142 ....: [ 9, 10, 11, 12],
1143 ....: [-10, 9, -12, 11] ])
1144 sage: _unembed_complex_matrix(A)
1145 [ -2*i + 1 -4*i + 3]
1146 [ -10*i + 9 -12*i + 11]
1150 raise ArgumentError("the matrix 'M' must be square")
1151 if not n
.mod(2).is_zero():
1152 raise ArgumentError("the matrix 'M' must be a complex embedding")
1154 F
= QuadraticField(-1, 'i')
1157 # Go top-left to bottom-right (reading order), converting every
1158 # 2-by-2 block we see to a single complex element.
1160 for k
in xrange(n
/2):
1161 for j
in xrange(n
/2):
1162 submat
= M
[2*k
:2*k
+2,2*j
:2*j
+2]
1163 if submat
[0,0] != submat
[1,1]:
1164 raise ArgumentError('bad real submatrix')
1165 if submat
[0,1] != -submat
[1,0]:
1166 raise ArgumentError('bad imag submatrix')
1167 z
= submat
[0,0] + submat
[1,0]*i
1170 return matrix(F
, n
/2, elements
)
1173 def RealSymmetricSimpleEJA(n
, field
=QQ
):
1175 The rank-n simple EJA consisting of real symmetric n-by-n
1176 matrices, the usual symmetric Jordan product, and the trace inner
1177 product. It has dimension `(n^2 + n)/2` over the reals.
1181 sage: J = RealSymmetricSimpleEJA(2)
1182 sage: e0, e1, e2 = J.gens()
1192 The degree of this algebra is `(n^2 + n) / 2`::
1194 sage: set_random_seed()
1195 sage: n = ZZ.random_element(1,5)
1196 sage: J = RealSymmetricSimpleEJA(n)
1197 sage: J.degree() == (n^2 + n)/2
1201 S
= _real_symmetric_basis(n
, field
=field
)
1202 (Qs
, T
) = _multiplication_table_from_matrix_basis(S
)
1204 return FiniteDimensionalEuclideanJordanAlgebra(field
,
1210 def ComplexHermitianSimpleEJA(n
, field
=QQ
):
1212 The rank-n simple EJA consisting of complex Hermitian n-by-n
1213 matrices over the real numbers, the usual symmetric Jordan product,
1214 and the real-part-of-trace inner product. It has dimension `n^2` over
1219 The degree of this algebra is `n^2`::
1221 sage: set_random_seed()
1222 sage: n = ZZ.random_element(1,5)
1223 sage: J = ComplexHermitianSimpleEJA(n)
1224 sage: J.degree() == n^2
1228 S
= _complex_hermitian_basis(n
)
1229 (Qs
, T
) = _multiplication_table_from_matrix_basis(S
)
1230 return FiniteDimensionalEuclideanJordanAlgebra(field
,
1236 def QuaternionHermitianSimpleEJA(n
):
1238 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
1239 matrices, the usual symmetric Jordan product, and the
1240 real-part-of-trace inner product. It has dimension `2n^2 - n` over
1245 def OctonionHermitianSimpleEJA(n
):
1247 This shit be crazy. It has dimension 27 over the reals.
1252 def JordanSpinSimpleEJA(n
, field
=QQ
):
1254 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
1255 with the usual inner product and jordan product ``x*y =
1256 (<x_bar,y_bar>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
1261 This multiplication table can be verified by hand::
1263 sage: J = JordanSpinSimpleEJA(4)
1264 sage: e0,e1,e2,e3 = J.gens()
1282 id_matrix
= identity_matrix(field
, n
)
1284 ei
= id_matrix
.column(i
)
1285 Qi
= zero_matrix(field
, n
)
1287 Qi
.set_column(0, ei
)
1288 Qi
+= diagonal_matrix(n
, [ei
[0]]*n
)
1289 # The addition of the diagonal matrix adds an extra ei[0] in the
1290 # upper-left corner of the matrix.
1291 Qi
[0,0] = Qi
[0,0] * ~
field(2)
1294 # The usual inner product on R^n.
1295 ip
= lambda x
, y
: x
.vector().inner_product(y
.vector())
1297 # The rank of the spin factor algebra is two, UNLESS we're in a
1298 # one-dimensional ambient space (the rank is bounded by the
1299 # ambient dimension).
1300 return FiniteDimensionalEuclideanJordanAlgebra(field
,