1 from sage
.matrix
.constructor
import matrix
2 from sage
.modules
.free_module
import VectorSpace
3 from sage
.modules
.with_basis
.indexed_element
import IndexedFreeModuleElement
5 # TODO: make this unnecessary somehow.
6 from sage
.misc
.lazy_import
import lazy_import
7 lazy_import('mjo.eja.eja_algebra', 'FiniteDimensionalEuclideanJordanAlgebra')
8 lazy_import('mjo.eja.eja_element_subalgebra',
9 'FiniteDimensionalEuclideanJordanElementSubalgebra')
10 from mjo
.eja
.eja_operator
import FiniteDimensionalEuclideanJordanAlgebraOperator
11 from mjo
.eja
.eja_utils
import _mat2vec
13 class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement
):
15 An element of a Euclidean Jordan algebra.
20 Oh man, I should not be doing this. This hides the "disabled"
21 methods ``left_matrix`` and ``matrix`` from introspection;
22 in particular it removes them from tab-completion.
24 return filter(lambda s
: s
not in ['left_matrix', 'matrix'],
32 Return ``self`` raised to the power ``n``.
34 Jordan algebras are always power-associative; see for
35 example Faraut and Korányi, Proposition II.1.2 (ii).
37 We have to override this because our superclass uses row
38 vectors instead of column vectors! We, on the other hand,
39 assume column vectors everywhere.
43 sage: from mjo.eja.eja_algebra import random_eja
47 The definition of `x^2` is the unambiguous `x*x`::
49 sage: set_random_seed()
50 sage: x = random_eja().random_element()
54 A few examples of power-associativity::
56 sage: set_random_seed()
57 sage: x = random_eja().random_element()
58 sage: x*(x*x)*(x*x) == x^5
60 sage: (x*x)*(x*x*x) == x^5
63 We also know that powers operator-commute (Koecher, Chapter
66 sage: set_random_seed()
67 sage: x = random_eja().random_element()
68 sage: m = ZZ.random_element(0,10)
69 sage: n = ZZ.random_element(0,10)
70 sage: Lxm = (x^m).operator()
71 sage: Lxn = (x^n).operator()
72 sage: Lxm*Lxn == Lxn*Lxm
77 return self
.parent().one()
81 return (self
**(n
-1))*self
84 def apply_univariate_polynomial(self
, p
):
86 Apply the univariate polynomial ``p`` to this element.
88 A priori, SageMath won't allow us to apply a univariate
89 polynomial to an element of an EJA, because we don't know
90 that EJAs are rings (they are usually not associative). Of
91 course, we know that EJAs are power-associative, so the
92 operation is ultimately kosher. This function sidesteps
93 the CAS to get the answer we want and expect.
97 sage: from mjo.eja.eja_algebra import (HadamardEJA,
102 sage: R = PolynomialRing(QQ, 't')
104 sage: p = t^4 - t^3 + 5*t - 2
105 sage: J = HadamardEJA(5)
106 sage: J.one().apply_univariate_polynomial(p) == 3*J.one()
111 We should always get back an element of the algebra::
113 sage: set_random_seed()
114 sage: p = PolynomialRing(AA, 't').random_element()
115 sage: J = random_eja()
116 sage: x = J.random_element()
117 sage: x.apply_univariate_polynomial(p) in J
121 if len(p
.variables()) > 1:
122 raise ValueError("not a univariate polynomial")
125 # Convert the coeficcients to the parent's base ring,
126 # because a priori they might live in an (unnecessarily)
127 # larger ring for which P.sum() would fail below.
128 cs
= [ R(c
) for c
in p
.coefficients(sparse
=False) ]
129 return P
.sum( cs
[k
]*(self
**k
) for k
in range(len(cs
)) )
132 def characteristic_polynomial(self
):
134 Return the characteristic polynomial of this element.
138 sage: from mjo.eja.eja_algebra import HadamardEJA
142 The rank of `R^3` is three, and the minimal polynomial of
143 the identity element is `(t-1)` from which it follows that
144 the characteristic polynomial should be `(t-1)^3`::
146 sage: J = HadamardEJA(3)
147 sage: J.one().characteristic_polynomial()
148 t^3 - 3*t^2 + 3*t - 1
150 Likewise, the characteristic of the zero element in the
151 rank-three algebra `R^{n}` should be `t^{3}`::
153 sage: J = HadamardEJA(3)
154 sage: J.zero().characteristic_polynomial()
159 The characteristic polynomial of an element should evaluate
160 to zero on that element::
162 sage: set_random_seed()
163 sage: x = HadamardEJA(3).random_element()
164 sage: p = x.characteristic_polynomial()
165 sage: x.apply_univariate_polynomial(p)
168 The characteristic polynomials of the zero and unit elements
169 should be what we think they are in a subalgebra, too::
171 sage: J = HadamardEJA(3)
172 sage: p1 = J.one().characteristic_polynomial()
173 sage: q1 = J.zero().characteristic_polynomial()
174 sage: e0,e1,e2 = J.gens()
175 sage: A = (e0 + 2*e1 + 3*e2).subalgebra_generated_by() # dim 3
176 sage: p2 = A.one().characteristic_polynomial()
177 sage: q2 = A.zero().characteristic_polynomial()
184 p
= self
.parent().characteristic_polynomial_of()
185 return p(*self
.to_vector())
188 def inner_product(self
, other
):
190 Return the parent algebra's inner product of myself and ``other``.
194 sage: from mjo.eja.eja_algebra import (
195 ....: ComplexHermitianEJA,
197 ....: QuaternionHermitianEJA,
198 ....: RealSymmetricEJA,
203 The inner product in the Jordan spin algebra is the usual
204 inner product on `R^n` (this example only works because the
205 basis for the Jordan algebra is the standard basis in `R^n`)::
207 sage: J = JordanSpinEJA(3)
208 sage: x = vector(QQ,[1,2,3])
209 sage: y = vector(QQ,[4,5,6])
210 sage: x.inner_product(y)
212 sage: J.from_vector(x).inner_product(J.from_vector(y))
215 The inner product on `S^n` is `<X,Y> = trace(X*Y)`, where
216 multiplication is the usual matrix multiplication in `S^n`,
217 so the inner product of the identity matrix with itself
220 sage: J = RealSymmetricEJA(3)
221 sage: J.one().inner_product(J.one())
224 Likewise, the inner product on `C^n` is `<X,Y> =
225 Re(trace(X*Y))`, where we must necessarily take the real
226 part because the product of Hermitian matrices may not be
229 sage: J = ComplexHermitianEJA(3)
230 sage: J.one().inner_product(J.one())
233 Ditto for the quaternions::
235 sage: J = QuaternionHermitianEJA(3)
236 sage: J.one().inner_product(J.one())
241 Ensure that we can always compute an inner product, and that
242 it gives us back a real number::
244 sage: set_random_seed()
245 sage: J = random_eja()
246 sage: x,y = J.random_elements(2)
247 sage: x.inner_product(y) in RLF
253 raise TypeError("'other' must live in the same algebra")
255 return P
.inner_product(self
, other
)
258 def operator_commutes_with(self
, other
):
260 Return whether or not this element operator-commutes
265 sage: from mjo.eja.eja_algebra import random_eja
269 The definition of a Jordan algebra says that any element
270 operator-commutes with its square::
272 sage: set_random_seed()
273 sage: x = random_eja().random_element()
274 sage: x.operator_commutes_with(x^2)
279 Test Lemma 1 from Chapter III of Koecher::
281 sage: set_random_seed()
282 sage: u,v = random_eja().random_elements(2)
283 sage: lhs = u.operator_commutes_with(u*v)
284 sage: rhs = v.operator_commutes_with(u^2)
288 Test the first polarization identity from my notes, Koecher
289 Chapter III, or from Baes (2.3)::
291 sage: set_random_seed()
292 sage: x,y = random_eja().random_elements(2)
293 sage: Lx = x.operator()
294 sage: Ly = y.operator()
295 sage: Lxx = (x*x).operator()
296 sage: Lxy = (x*y).operator()
297 sage: bool(2*Lx*Lxy + Ly*Lxx == 2*Lxy*Lx + Lxx*Ly)
300 Test the second polarization identity from my notes or from
303 sage: set_random_seed()
304 sage: x,y,z = random_eja().random_elements(3)
305 sage: Lx = x.operator()
306 sage: Ly = y.operator()
307 sage: Lz = z.operator()
308 sage: Lzy = (z*y).operator()
309 sage: Lxy = (x*y).operator()
310 sage: Lxz = (x*z).operator()
311 sage: bool(Lx*Lzy + Lz*Lxy + Ly*Lxz == Lzy*Lx + Lxy*Lz + Lxz*Ly)
314 Test the third polarization identity from my notes or from
317 sage: set_random_seed()
318 sage: u,y,z = random_eja().random_elements(3)
319 sage: Lu = u.operator()
320 sage: Ly = y.operator()
321 sage: Lz = z.operator()
322 sage: Lzy = (z*y).operator()
323 sage: Luy = (u*y).operator()
324 sage: Luz = (u*z).operator()
325 sage: Luyz = (u*(y*z)).operator()
326 sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz
327 sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly
328 sage: bool(lhs == rhs)
332 if not other
in self
.parent():
333 raise TypeError("'other' must live in the same algebra")
342 Return my determinant, the product of my eigenvalues.
346 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
352 sage: J = JordanSpinEJA(2)
353 sage: e0,e1 = J.gens()
354 sage: x = sum( J.gens() )
360 sage: J = JordanSpinEJA(3)
361 sage: e0,e1,e2 = J.gens()
362 sage: x = sum( J.gens() )
366 The determinant of the sole element in the rank-zero trivial
367 algebra is ``1``, by three paths of reasoning. First, its
368 characteristic polynomial is a constant ``1``, so the constant
369 term in that polynomial is ``1``. Second, the characteristic
370 polynomial evaluated at zero is again ``1``. And finally, the
371 (empty) product of its eigenvalues is likewise just unity::
373 sage: J = TrivialEJA()
379 An element is invertible if and only if its determinant is
382 sage: set_random_seed()
383 sage: x = random_eja().random_element()
384 sage: x.is_invertible() == (x.det() != 0)
387 Ensure that the determinant is multiplicative on an associative
388 subalgebra as in Faraut and Korányi's Proposition II.2.2::
390 sage: set_random_seed()
391 sage: J = random_eja().random_element().subalgebra_generated_by()
392 sage: x,y = J.random_elements(2)
393 sage: (x*y).det() == x.det()*y.det()
400 # Special case, since we don't get the a0=1
401 # coefficient when the rank of the algebra
403 return P
.base_ring().one()
405 p
= P
._charpoly
_coefficients
()[0]
406 # The _charpoly_coeff function already adds the factor of -1
407 # to ensure that _charpoly_coefficients()[0] is really what
408 # appears in front of t^{0} in the charpoly. However, we want
409 # (-1)^r times THAT for the determinant.
410 return ((-1)**r
)*p(*self
.to_vector())
415 Return the Jordan-multiplicative inverse of this element.
419 We appeal to the quadratic representation as in Koecher's
420 Theorem 12 in Chapter III, Section 5.
424 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
430 The inverse in the spin factor algebra is given in Alizadeh's
433 sage: set_random_seed()
434 sage: J = JordanSpinEJA.random_instance()
435 sage: x = J.random_element()
436 sage: while not x.is_invertible():
437 ....: x = J.random_element()
438 sage: x_vec = x.to_vector()
440 sage: x_bar = x_vec[1:]
441 sage: coeff = x0.inner_product(x0) - x_bar.inner_product(x_bar)
442 sage: x_inverse = x_vec.parent()(x0.list() + (-x_bar).list())
443 sage: if not coeff.is_zero(): x_inverse = x_inverse/coeff
444 sage: x.inverse() == J.from_vector(x_inverse)
447 Trying to invert a non-invertible element throws an error:
449 sage: JordanSpinEJA(3).zero().inverse()
450 Traceback (most recent call last):
452 ValueError: element is not invertible
456 The identity element is its own inverse::
458 sage: set_random_seed()
459 sage: J = random_eja()
460 sage: J.one().inverse() == J.one()
463 If an element has an inverse, it acts like one::
465 sage: set_random_seed()
466 sage: J = random_eja()
467 sage: x = J.random_element()
468 sage: (not x.is_invertible()) or (x.inverse()*x == J.one())
471 The inverse of the inverse is what we started with::
473 sage: set_random_seed()
474 sage: J = random_eja()
475 sage: x = J.random_element()
476 sage: (not x.is_invertible()) or (x.inverse().inverse() == x)
479 Proposition II.2.3 in Faraut and Korányi says that the inverse
480 of an element is the inverse of its left-multiplication operator
481 applied to the algebra's identity, when that inverse exists::
483 sage: set_random_seed()
484 sage: J = random_eja()
485 sage: x = J.random_element()
486 sage: (not x.operator().is_invertible()) or (
487 ....: x.operator().inverse()(J.one()) == x.inverse() )
490 Proposition II.2.4 in Faraut and Korányi gives a formula for
491 the inverse based on the characteristic polynomial and the
492 Cayley-Hamilton theorem for Euclidean Jordan algebras::
494 sage: set_random_seed()
495 sage: J = ComplexHermitianEJA(3)
496 sage: x = J.random_element()
497 sage: while not x.is_invertible():
498 ....: x = J.random_element()
500 sage: a = x.characteristic_polynomial().coefficients(sparse=False)
501 sage: expected = (-1)^(r+1)/x.det()
502 sage: expected *= sum( a[i+1]*x^i for i in range(r) )
503 sage: x.inverse() == expected
507 if not self
.is_invertible():
508 raise ValueError("element is not invertible")
510 return (~self
.quadratic_representation())(self
)
513 def is_invertible(self
):
515 Return whether or not this element is invertible.
519 The usual way to do this is to check if the determinant is
520 zero, but we need the characteristic polynomial for the
521 determinant. The minimal polynomial is a lot easier to get,
522 so we use Corollary 2 in Chapter V of Koecher to check
523 whether or not the paren't algebra's zero element is a root
524 of this element's minimal polynomial.
526 That is... unless the coefficients of our algebra's
527 "characteristic polynomial of" function are already cached!
528 In that case, we just use the determinant (which will be fast
531 Beware that we can't use the superclass method, because it
532 relies on the algebra being associative.
536 sage: from mjo.eja.eja_algebra import random_eja
540 The identity element is always invertible::
542 sage: set_random_seed()
543 sage: J = random_eja()
544 sage: J.one().is_invertible()
547 The zero element is never invertible in a non-trivial algebra::
549 sage: set_random_seed()
550 sage: J = random_eja()
551 sage: (not J.is_trivial()) and J.zero().is_invertible()
556 if self
.parent().is_trivial():
561 if self
.parent()._charpoly
_coefficients
.is_in_cache():
562 # The determinant will be quicker than computing the minimal
563 # polynomial from scratch, most likely.
564 return (not self
.det().is_zero())
566 # In fact, we only need to know if the constant term is non-zero,
567 # so we can pass in the field's zero element instead.
568 zero
= self
.base_ring().zero()
569 p
= self
.minimal_polynomial()
570 return not (p(zero
) == zero
)
573 def is_primitive_idempotent(self
):
575 Return whether or not this element is a primitive (or minimal)
578 A primitive idempotent is a non-zero idempotent that is not
579 the sum of two other non-zero idempotents. Remark 2.7.15 in
580 Baes shows that this is what he refers to as a "minimal
583 An element of a Euclidean Jordan algebra is a minimal idempotent
584 if it :meth:`is_idempotent` and if its Peirce subalgebra
585 corresponding to the eigenvalue ``1`` has dimension ``1`` (Baes,
590 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
591 ....: RealSymmetricEJA,
597 This method is sloooooow.
601 The spectral decomposition of a non-regular element should always
602 contain at least one non-minimal idempotent::
604 sage: J = RealSymmetricEJA(3)
605 sage: x = sum(J.gens())
608 sage: [ c.is_primitive_idempotent()
609 ....: for (l,c) in x.spectral_decomposition() ]
612 On the other hand, the spectral decomposition of a regular
613 element should always be in terms of minimal idempotents::
615 sage: J = JordanSpinEJA(4)
616 sage: x = sum( i*J.gens()[i] for i in range(len(J.gens())) )
619 sage: [ c.is_primitive_idempotent()
620 ....: for (l,c) in x.spectral_decomposition() ]
625 The identity element is minimal only in an EJA of rank one::
627 sage: set_random_seed()
628 sage: J = random_eja()
629 sage: J.rank() == 1 or not J.one().is_primitive_idempotent()
632 A non-idempotent cannot be a minimal idempotent::
634 sage: set_random_seed()
635 sage: J = JordanSpinEJA(4)
636 sage: x = J.random_element()
637 sage: (not x.is_idempotent()) and x.is_primitive_idempotent()
640 Proposition 2.7.19 in Baes says that an element is a minimal
641 idempotent if and only if it's idempotent with trace equal to
644 sage: set_random_seed()
645 sage: J = JordanSpinEJA(4)
646 sage: x = J.random_element()
647 sage: expected = (x.is_idempotent() and x.trace() == 1)
648 sage: actual = x.is_primitive_idempotent()
649 sage: actual == expected
652 Primitive idempotents must be non-zero::
654 sage: set_random_seed()
655 sage: J = random_eja()
656 sage: J.zero().is_idempotent()
658 sage: J.zero().is_primitive_idempotent()
661 As a consequence of the fact that primitive idempotents must
662 be non-zero, there are no primitive idempotents in a trivial
663 Euclidean Jordan algebra::
665 sage: J = TrivialEJA()
666 sage: J.one().is_idempotent()
668 sage: J.one().is_primitive_idempotent()
672 if not self
.is_idempotent():
678 (_
,_
,J1
) = self
.parent().peirce_decomposition(self
)
679 return (J1
.dimension() == 1)
682 def is_nilpotent(self
):
684 Return whether or not some power of this element is zero.
688 We use Theorem 5 in Chapter III of Koecher, which says that
689 an element ``x`` is nilpotent if and only if ``x.operator()``
690 is nilpotent. And it is a basic fact of linear algebra that
691 an operator on an `n`-dimensional space is nilpotent if and
692 only if, when raised to the `n`th power, it equals the zero
693 operator (for example, see Axler Corollary 8.8).
697 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
702 sage: J = JordanSpinEJA(3)
703 sage: x = sum(J.gens())
704 sage: x.is_nilpotent()
709 The identity element is never nilpotent, except in a trivial EJA::
711 sage: set_random_seed()
712 sage: J = random_eja()
713 sage: J.one().is_nilpotent() and not J.is_trivial()
716 The additive identity is always nilpotent::
718 sage: set_random_seed()
719 sage: random_eja().zero().is_nilpotent()
724 zero_operator
= P
.zero().operator()
725 return self
.operator()**P
.dimension() == zero_operator
728 def is_regular(self
):
730 Return whether or not this is a regular element.
734 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
739 The identity element always has degree one, but any element
740 linearly-independent from it is regular::
742 sage: J = JordanSpinEJA(5)
743 sage: J.one().is_regular()
745 sage: e0, e1, e2, e3, e4 = J.gens() # e0 is the identity
746 sage: for x in J.gens():
747 ....: (J.one() + x).is_regular()
756 The zero element should never be regular, unless the parent
757 algebra has dimension less than or equal to one::
759 sage: set_random_seed()
760 sage: J = random_eja()
761 sage: J.dimension() <= 1 or not J.zero().is_regular()
764 The unit element isn't regular unless the algebra happens to
765 consist of only its scalar multiples::
767 sage: set_random_seed()
768 sage: J = random_eja()
769 sage: J.dimension() <= 1 or not J.one().is_regular()
773 return self
.degree() == self
.parent().rank()
778 Return the degree of this element, which is defined to be
779 the degree of its minimal polynomial.
783 For now, we skip the messy minimal polynomial computation
784 and instead return the dimension of the vector space spanned
785 by the powers of this element. The latter is a bit more
786 straightforward to compute.
790 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
795 sage: J = JordanSpinEJA(4)
796 sage: J.one().degree()
798 sage: e0,e1,e2,e3 = J.gens()
799 sage: (e0 - e1).degree()
802 In the spin factor algebra (of rank two), all elements that
803 aren't multiples of the identity are regular::
805 sage: set_random_seed()
806 sage: J = JordanSpinEJA.random_instance()
807 sage: n = J.dimension()
808 sage: x = J.random_element()
809 sage: x.degree() == min(n,2) or (x == x.coefficient(0)*J.one())
814 The zero and unit elements are both of degree one in nontrivial
817 sage: set_random_seed()
818 sage: J = random_eja()
819 sage: d = J.zero().degree()
820 sage: (J.is_trivial() and d == 0) or d == 1
822 sage: d = J.one().degree()
823 sage: (J.is_trivial() and d == 0) or d == 1
826 Our implementation agrees with the definition::
828 sage: set_random_seed()
829 sage: x = random_eja().random_element()
830 sage: x.degree() == x.minimal_polynomial().degree()
834 if self
.is_zero() and not self
.parent().is_trivial():
835 # The minimal polynomial of zero in a nontrivial algebra
836 # is "t"; in a trivial algebra it's "1" by convention
837 # (it's an empty product).
839 return self
.subalgebra_generated_by().dimension()
842 def left_matrix(self
):
844 Our parent class defines ``left_matrix`` and ``matrix``
845 methods whose names are misleading. We don't want them.
847 raise NotImplementedError("use operator().matrix() instead")
852 def minimal_polynomial(self
):
854 Return the minimal polynomial of this element,
855 as a function of the variable `t`.
859 We restrict ourselves to the associative subalgebra
860 generated by this element, and then return the minimal
861 polynomial of this element's operator matrix (in that
862 subalgebra). This works by Baes Proposition 2.3.16.
866 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
867 ....: RealSymmetricEJA,
873 Keeping in mind that the polynomial ``1`` evaluates the identity
874 element (also the zero element) of the trivial algebra, it is clear
875 that the polynomial ``1`` is the minimal polynomial of the only
876 element in a trivial algebra::
878 sage: J = TrivialEJA()
879 sage: J.one().minimal_polynomial()
881 sage: J.zero().minimal_polynomial()
886 The minimal polynomial of the identity and zero elements are
887 always the same, except in trivial algebras where the minimal
888 polynomial of the unit/zero element is ``1``::
890 sage: set_random_seed()
891 sage: J = random_eja()
892 sage: mu = J.one().minimal_polynomial()
893 sage: t = mu.parent().gen()
894 sage: mu + int(J.is_trivial())*(t-2)
896 sage: mu = J.zero().minimal_polynomial()
897 sage: t = mu.parent().gen()
898 sage: mu + int(J.is_trivial())*(t-1)
901 The degree of an element is (by one definition) the degree
902 of its minimal polynomial::
904 sage: set_random_seed()
905 sage: x = random_eja().random_element()
906 sage: x.degree() == x.minimal_polynomial().degree()
909 The minimal polynomial and the characteristic polynomial coincide
910 and are known (see Alizadeh, Example 11.11) for all elements of
911 the spin factor algebra that aren't scalar multiples of the
912 identity. We require the dimension of the algebra to be at least
913 two here so that said elements actually exist::
915 sage: set_random_seed()
916 sage: n_max = max(2, JordanSpinEJA._max_random_instance_size())
917 sage: n = ZZ.random_element(2, n_max)
918 sage: J = JordanSpinEJA(n)
919 sage: y = J.random_element()
920 sage: while y == y.coefficient(0)*J.one():
921 ....: y = J.random_element()
922 sage: y0 = y.to_vector()[0]
923 sage: y_bar = y.to_vector()[1:]
924 sage: actual = y.minimal_polynomial()
925 sage: t = PolynomialRing(J.base_ring(),'t').gen(0)
926 sage: expected = t^2 - 2*y0*t + (y0^2 - norm(y_bar)^2)
927 sage: bool(actual == expected)
930 The minimal polynomial should always kill its element::
932 sage: set_random_seed()
933 sage: x = random_eja().random_element()
934 sage: p = x.minimal_polynomial()
935 sage: x.apply_univariate_polynomial(p)
938 The minimal polynomial is invariant under a change of basis,
939 and in particular, a re-scaling of the basis::
941 sage: set_random_seed()
942 sage: n_max = RealSymmetricEJA._max_random_instance_size()
943 sage: n = ZZ.random_element(1, n_max)
944 sage: J1 = RealSymmetricEJA(n)
945 sage: J2 = RealSymmetricEJA(n,normalize_basis=False)
946 sage: X = random_matrix(AA,n)
947 sage: X = X*X.transpose()
950 sage: x1.minimal_polynomial() == x2.minimal_polynomial()
955 # We would generate a zero-dimensional subalgebra
956 # where the minimal polynomial would be constant.
957 # That might be correct, but only if *this* algebra
959 if not self
.parent().is_trivial():
960 # Pretty sure we know what the minimal polynomial of
961 # the zero operator is going to be. This ensures
962 # consistency of e.g. the polynomial variable returned
963 # in the "normal" case without us having to think about it.
964 return self
.operator().minimal_polynomial()
966 A
= self
.subalgebra_generated_by(orthonormalize_basis
=False)
967 return A(self
).operator().minimal_polynomial()
971 def natural_representation(self
):
973 Return a more-natural representation of this element.
975 Every finite-dimensional Euclidean Jordan Algebra is a
976 direct sum of five simple algebras, four of which comprise
977 Hermitian matrices. This method returns the original
978 "natural" representation of this element as a Hermitian
979 matrix, if it has one. If not, you get the usual representation.
983 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
984 ....: QuaternionHermitianEJA)
988 sage: J = ComplexHermitianEJA(3)
991 sage: J.one().natural_representation()
1001 sage: J = QuaternionHermitianEJA(3)
1004 sage: J.one().natural_representation()
1005 [1 0 0 0 0 0 0 0 0 0 0 0]
1006 [0 1 0 0 0 0 0 0 0 0 0 0]
1007 [0 0 1 0 0 0 0 0 0 0 0 0]
1008 [0 0 0 1 0 0 0 0 0 0 0 0]
1009 [0 0 0 0 1 0 0 0 0 0 0 0]
1010 [0 0 0 0 0 1 0 0 0 0 0 0]
1011 [0 0 0 0 0 0 1 0 0 0 0 0]
1012 [0 0 0 0 0 0 0 1 0 0 0 0]
1013 [0 0 0 0 0 0 0 0 1 0 0 0]
1014 [0 0 0 0 0 0 0 0 0 1 0 0]
1015 [0 0 0 0 0 0 0 0 0 0 1 0]
1016 [0 0 0 0 0 0 0 0 0 0 0 1]
1019 B
= self
.parent().natural_basis()
1020 W
= self
.parent().natural_basis_space()
1022 # This is just a manual "from_vector()", but of course
1023 # matrix spaces aren't vector spaces in sage, so they
1024 # don't have a from_vector() method.
1025 return W
.linear_combination(zip(B
,self
.to_vector()))
1030 The norm of this element with respect to :meth:`inner_product`.
1034 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1039 sage: J = HadamardEJA(2)
1040 sage: x = sum(J.gens())
1046 sage: J = JordanSpinEJA(4)
1047 sage: x = sum(J.gens())
1052 return self
.inner_product(self
).sqrt()
1057 Return the left-multiplication-by-this-element
1058 operator on the ambient algebra.
1062 sage: from mjo.eja.eja_algebra import random_eja
1066 sage: set_random_seed()
1067 sage: J = random_eja()
1068 sage: x,y = J.random_elements(2)
1069 sage: x.operator()(y) == x*y
1071 sage: y.operator()(x) == x*y
1076 left_mult_by_self
= lambda y
: self
*y
1077 L
= P
.module_morphism(function
=left_mult_by_self
, codomain
=P
)
1078 return FiniteDimensionalEuclideanJordanAlgebraOperator(
1084 def quadratic_representation(self
, other
=None):
1086 Return the quadratic representation of this element.
1090 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1095 The explicit form in the spin factor algebra is given by
1096 Alizadeh's Example 11.12::
1098 sage: set_random_seed()
1099 sage: x = JordanSpinEJA.random_instance().random_element()
1100 sage: x_vec = x.to_vector()
1101 sage: Q = matrix.identity(x.base_ring(), 0)
1102 sage: n = x_vec.degree()
1105 ....: x_bar = x_vec[1:]
1106 ....: A = matrix(x.base_ring(), 1, [x_vec.inner_product(x_vec)])
1107 ....: B = 2*x0*x_bar.row()
1108 ....: C = 2*x0*x_bar.column()
1109 ....: D = matrix.identity(x.base_ring(), n-1)
1110 ....: D = (x0^2 - x_bar.inner_product(x_bar))*D
1111 ....: D = D + 2*x_bar.tensor_product(x_bar)
1112 ....: Q = matrix.block(2,2,[A,B,C,D])
1113 sage: Q == x.quadratic_representation().matrix()
1116 Test all of the properties from Theorem 11.2 in Alizadeh::
1118 sage: set_random_seed()
1119 sage: J = random_eja()
1120 sage: x,y = J.random_elements(2)
1121 sage: Lx = x.operator()
1122 sage: Lxx = (x*x).operator()
1123 sage: Qx = x.quadratic_representation()
1124 sage: Qy = y.quadratic_representation()
1125 sage: Qxy = x.quadratic_representation(y)
1126 sage: Qex = J.one().quadratic_representation(x)
1127 sage: n = ZZ.random_element(10)
1128 sage: Qxn = (x^n).quadratic_representation()
1132 sage: 2*Qxy == (x+y).quadratic_representation() - Qx - Qy
1135 Property 2 (multiply on the right for :trac:`28272`):
1137 sage: alpha = J.base_ring().random_element()
1138 sage: (alpha*x).quadratic_representation() == Qx*(alpha^2)
1143 sage: not x.is_invertible() or ( Qx(x.inverse()) == x )
1146 sage: not x.is_invertible() or (
1149 ....: x.inverse().quadratic_representation() )
1152 sage: Qxy(J.one()) == x*y
1157 sage: not x.is_invertible() or (
1158 ....: x.quadratic_representation(x.inverse())*Qx
1159 ....: == Qx*x.quadratic_representation(x.inverse()) )
1162 sage: not x.is_invertible() or (
1163 ....: x.quadratic_representation(x.inverse())*Qx
1165 ....: 2*Lx*Qex - Qx )
1168 sage: 2*Lx*Qex - Qx == Lxx
1173 sage: Qy(x).quadratic_representation() == Qy*Qx*Qy
1183 sage: not x.is_invertible() or (
1184 ....: Qx*x.inverse().operator() == Lx )
1189 sage: not x.operator_commutes_with(y) or (
1190 ....: Qx(y)^n == Qxn(y^n) )
1196 elif not other
in self
.parent():
1197 raise TypeError("'other' must live in the same algebra")
1200 M
= other
.operator()
1201 return ( L
*M
+ M
*L
- (self
*other
).operator() )
1205 def spectral_decomposition(self
):
1207 Return the unique spectral decomposition of this element.
1211 Following Faraut and Korányi's Theorem III.1.1, we restrict this
1212 element's left-multiplication-by operator to the subalgebra it
1213 generates. We then compute the spectral decomposition of that
1214 operator, and the spectral projectors we get back must be the
1215 left-multiplication-by operators for the idempotents we
1216 seek. Thus applying them to the identity element gives us those
1219 Since the eigenvalues are required to be distinct, we take
1220 the spectral decomposition of the zero element to be zero
1221 times the identity element of the algebra (which is idempotent,
1226 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1230 The spectral decomposition of the identity is ``1`` times itself,
1231 and the spectral decomposition of zero is ``0`` times the identity::
1233 sage: J = RealSymmetricEJA(3)
1236 sage: J.one().spectral_decomposition()
1238 sage: J.zero().spectral_decomposition()
1243 sage: J = RealSymmetricEJA(4)
1244 sage: x = sum(J.gens())
1245 sage: sd = x.spectral_decomposition()
1250 sage: c0.inner_product(c1) == 0
1252 sage: c0.is_idempotent()
1254 sage: c1.is_idempotent()
1256 sage: c0 + c1 == J.one()
1258 sage: l0*c0 + l1*c1 == x
1261 The spectral decomposition should work in subalgebras, too::
1263 sage: J = RealSymmetricEJA(4)
1264 sage: (e0, e1, e2, e3, e4, e5, e6, e7, e8, e9) = J.gens()
1265 sage: A = 2*e5 - 2*e8
1266 sage: (lambda1, c1) = A.spectral_decomposition()[1]
1267 sage: (J0, J5, J1) = J.peirce_decomposition(c1)
1268 sage: (f0, f1, f2) = J1.gens()
1269 sage: f0.spectral_decomposition()
1273 A
= self
.subalgebra_generated_by(orthonormalize_basis
=True)
1275 for (evalue
, proj
) in A(self
).operator().spectral_decomposition():
1276 result
.append( (evalue
, proj(A
.one()).superalgebra_element()) )
1279 def subalgebra_generated_by(self
, orthonormalize_basis
=False):
1281 Return the associative subalgebra of the parent EJA generated
1284 Since our parent algebra is unital, we want "subalgebra" to mean
1285 "unital subalgebra" as well; thus the subalgebra that an element
1286 generates will itself be a Euclidean Jordan algebra after
1287 restricting the algebra operations appropriately. This is the
1288 subalgebra that Faraut and Korányi work with in section II.2, for
1293 sage: from mjo.eja.eja_algebra import random_eja
1297 This subalgebra, being composed of only powers, is associative::
1299 sage: set_random_seed()
1300 sage: x0 = random_eja().random_element()
1301 sage: A = x0.subalgebra_generated_by()
1302 sage: x,y,z = A.random_elements(3)
1303 sage: (x*y)*z == x*(y*z)
1306 Squaring in the subalgebra should work the same as in
1309 sage: set_random_seed()
1310 sage: x = random_eja().random_element()
1311 sage: A = x.subalgebra_generated_by()
1312 sage: A(x^2) == A(x)*A(x)
1315 By definition, the subalgebra generated by the zero element is
1316 the one-dimensional algebra generated by the identity
1317 element... unless the original algebra was trivial, in which
1318 case the subalgebra is trivial too::
1320 sage: set_random_seed()
1321 sage: A = random_eja().zero().subalgebra_generated_by()
1322 sage: (A.is_trivial() and A.dimension() == 0) or A.dimension() == 1
1326 return FiniteDimensionalEuclideanJordanElementSubalgebra(self
, orthonormalize_basis
)
1329 def subalgebra_idempotent(self
):
1331 Find an idempotent in the associative subalgebra I generate
1332 using Proposition 2.3.5 in Baes.
1336 sage: from mjo.eja.eja_algebra import random_eja
1340 Ensure that we can find an idempotent in a non-trivial algebra
1341 where there are non-nilpotent elements, or that we get the dumb
1342 solution in the trivial algebra::
1344 sage: set_random_seed()
1345 sage: J = random_eja()
1346 sage: x = J.random_element()
1347 sage: while x.is_nilpotent() and not J.is_trivial():
1348 ....: x = J.random_element()
1349 sage: c = x.subalgebra_idempotent()
1354 if self
.parent().is_trivial():
1357 if self
.is_nilpotent():
1358 raise ValueError("this only works with non-nilpotent elements!")
1360 J
= self
.subalgebra_generated_by()
1363 # The image of the matrix of left-u^m-multiplication
1364 # will be minimal for some natural number s...
1366 minimal_dim
= J
.dimension()
1367 for i
in range(1, minimal_dim
):
1368 this_dim
= (u
**i
).operator().matrix().image().dimension()
1369 if this_dim
< minimal_dim
:
1370 minimal_dim
= this_dim
1373 # Now minimal_matrix should correspond to the smallest
1374 # non-zero subspace in Baes's (or really, Koecher's)
1377 # However, we need to restrict the matrix to work on the
1378 # subspace... or do we? Can't we just solve, knowing that
1379 # A(c) = u^(s+1) should have a solution in the big space,
1382 # Beware, solve_right() means that we're using COLUMN vectors.
1383 # Our FiniteDimensionalAlgebraElement superclass uses rows.
1385 A
= u_next
.operator().matrix()
1386 c
= J
.from_vector(A
.solve_right(u_next
.to_vector()))
1388 # Now c is the idempotent we want, but it still lives in the subalgebra.
1389 return c
.superalgebra_element()
1394 Return my trace, the sum of my eigenvalues.
1396 In a trivial algebra, however you want to look at it, the trace is
1397 an empty sum for which we declare the result to be zero.
1401 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1408 sage: J = TrivialEJA()
1409 sage: J.zero().trace()
1413 sage: J = JordanSpinEJA(3)
1414 sage: x = sum(J.gens())
1420 sage: J = HadamardEJA(5)
1421 sage: J.one().trace()
1426 The trace of an element is a real number::
1428 sage: set_random_seed()
1429 sage: J = random_eja()
1430 sage: J.random_element().trace() in RLF
1438 # Special case for the trivial algebra where
1439 # the trace is an empty sum.
1440 return P
.base_ring().zero()
1442 p
= P
._charpoly
_coefficients
()[r
-1]
1443 # The _charpoly_coeff function already adds the factor of
1444 # -1 to ensure that _charpoly_coeff(r-1) is really what
1445 # appears in front of t^{r-1} in the charpoly. However,
1446 # we want the negative of THAT for the trace.
1447 return -p(*self
.to_vector())
1450 def trace_inner_product(self
, other
):
1452 Return the trace inner product of myself and ``other``.
1456 sage: from mjo.eja.eja_algebra import random_eja
1460 The trace inner product is commutative, bilinear, and associative::
1462 sage: set_random_seed()
1463 sage: J = random_eja()
1464 sage: x,y,z = J.random_elements(3)
1466 sage: x.trace_inner_product(y) == y.trace_inner_product(x)
1469 sage: a = J.base_ring().random_element();
1470 sage: actual = (a*(x+z)).trace_inner_product(y)
1471 sage: expected = ( a*x.trace_inner_product(y) +
1472 ....: a*z.trace_inner_product(y) )
1473 sage: actual == expected
1475 sage: actual = x.trace_inner_product(a*(y+z))
1476 sage: expected = ( a*x.trace_inner_product(y) +
1477 ....: a*x.trace_inner_product(z) )
1478 sage: actual == expected
1481 sage: (x*y).trace_inner_product(z) == y.trace_inner_product(x*z)
1485 if not other
in self
.parent():
1486 raise TypeError("'other' must live in the same algebra")
1488 return (self
*other
).trace()
1491 def trace_norm(self
):
1493 The norm of this element with respect to :meth:`trace_inner_product`.
1497 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1502 sage: J = HadamardEJA(2)
1503 sage: x = sum(J.gens())
1504 sage: x.trace_norm()
1509 sage: J = JordanSpinEJA(4)
1510 sage: x = sum(J.gens())
1511 sage: x.trace_norm()
1515 return self
.trace_inner_product(self
).sqrt()