1 # -*- coding: utf-8 -*-
3 from itertools
import izip
5 from sage
.matrix
.constructor
import matrix
6 from sage
.modules
.free_module
import VectorSpace
7 from sage
.modules
.with_basis
.indexed_element
import IndexedFreeModuleElement
9 # TODO: make this unnecessary somehow.
10 from sage
.misc
.lazy_import
import lazy_import
11 lazy_import('mjo.eja.eja_algebra', 'FiniteDimensionalEuclideanJordanAlgebra')
12 lazy_import('mjo.eja.eja_subalgebra',
13 'FiniteDimensionalEuclideanJordanElementSubalgebra')
14 from mjo
.eja
.eja_operator
import FiniteDimensionalEuclideanJordanAlgebraOperator
15 from mjo
.eja
.eja_utils
import _mat2vec
17 class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement
):
19 An element of a Euclidean Jordan algebra.
24 Oh man, I should not be doing this. This hides the "disabled"
25 methods ``left_matrix`` and ``matrix`` from introspection;
26 in particular it removes them from tab-completion.
28 return filter(lambda s
: s
not in ['left_matrix', 'matrix'],
36 Return ``self`` raised to the power ``n``.
38 Jordan algebras are always power-associative; see for
39 example Faraut and Korányi, Proposition II.1.2 (ii).
41 We have to override this because our superclass uses row
42 vectors instead of column vectors! We, on the other hand,
43 assume column vectors everywhere.
47 sage: from mjo.eja.eja_algebra import random_eja
51 The definition of `x^2` is the unambiguous `x*x`::
53 sage: set_random_seed()
54 sage: x = random_eja().random_element()
58 A few examples of power-associativity::
60 sage: set_random_seed()
61 sage: x = random_eja().random_element()
62 sage: x*(x*x)*(x*x) == x^5
64 sage: (x*x)*(x*x*x) == x^5
67 We also know that powers operator-commute (Koecher, Chapter
70 sage: set_random_seed()
71 sage: x = random_eja().random_element()
72 sage: m = ZZ.random_element(0,10)
73 sage: n = ZZ.random_element(0,10)
74 sage: Lxm = (x^m).operator()
75 sage: Lxn = (x^n).operator()
76 sage: Lxm*Lxn == Lxn*Lxm
81 return self
.parent().one()
85 return (self
**(n
-1))*self
88 def apply_univariate_polynomial(self
, p
):
90 Apply the univariate polynomial ``p`` to this element.
92 A priori, SageMath won't allow us to apply a univariate
93 polynomial to an element of an EJA, because we don't know
94 that EJAs are rings (they are usually not associative). Of
95 course, we know that EJAs are power-associative, so the
96 operation is ultimately kosher. This function sidesteps
97 the CAS to get the answer we want and expect.
101 sage: from mjo.eja.eja_algebra import (RealCartesianProductEJA,
106 sage: R = PolynomialRing(QQ, 't')
108 sage: p = t^4 - t^3 + 5*t - 2
109 sage: J = RealCartesianProductEJA(5)
110 sage: J.one().apply_univariate_polynomial(p) == 3*J.one()
115 We should always get back an element of the algebra::
117 sage: set_random_seed()
118 sage: p = PolynomialRing(QQ, 't').random_element()
119 sage: J = random_eja()
120 sage: x = J.random_element()
121 sage: x.apply_univariate_polynomial(p) in J
125 if len(p
.variables()) > 1:
126 raise ValueError("not a univariate polynomial")
129 # Convert the coeficcients to the parent's base ring,
130 # because a priori they might live in an (unnecessarily)
131 # larger ring for which P.sum() would fail below.
132 cs
= [ R(c
) for c
in p
.coefficients(sparse
=False) ]
133 return P
.sum( cs
[k
]*(self
**k
) for k
in range(len(cs
)) )
136 def characteristic_polynomial(self
):
138 Return the characteristic polynomial of this element.
142 sage: from mjo.eja.eja_algebra import RealCartesianProductEJA
146 The rank of `R^3` is three, and the minimal polynomial of
147 the identity element is `(t-1)` from which it follows that
148 the characteristic polynomial should be `(t-1)^3`::
150 sage: J = RealCartesianProductEJA(3)
151 sage: J.one().characteristic_polynomial()
152 t^3 - 3*t^2 + 3*t - 1
154 Likewise, the characteristic of the zero element in the
155 rank-three algebra `R^{n}` should be `t^{3}`::
157 sage: J = RealCartesianProductEJA(3)
158 sage: J.zero().characteristic_polynomial()
163 The characteristic polynomial of an element should evaluate
164 to zero on that element::
166 sage: set_random_seed()
167 sage: x = RealCartesianProductEJA(3).random_element()
168 sage: p = x.characteristic_polynomial()
169 sage: x.apply_univariate_polynomial(p)
172 The characteristic polynomials of the zero and unit elements
173 should be what we think they are in a subalgebra, too::
175 sage: J = RealCartesianProductEJA(3)
176 sage: p1 = J.one().characteristic_polynomial()
177 sage: q1 = J.zero().characteristic_polynomial()
178 sage: e0,e1,e2 = J.gens()
179 sage: A = (e0 + 2*e1 + 3*e2).subalgebra_generated_by() # dim 3
180 sage: p2 = A.one().characteristic_polynomial()
181 sage: q2 = A.zero().characteristic_polynomial()
188 p
= self
.parent().characteristic_polynomial()
189 return p(*self
.to_vector())
192 def inner_product(self
, other
):
194 Return the parent algebra's inner product of myself and ``other``.
198 sage: from mjo.eja.eja_algebra import (
199 ....: ComplexHermitianEJA,
201 ....: QuaternionHermitianEJA,
202 ....: RealSymmetricEJA,
207 The inner product in the Jordan spin algebra is the usual
208 inner product on `R^n` (this example only works because the
209 basis for the Jordan algebra is the standard basis in `R^n`)::
211 sage: J = JordanSpinEJA(3)
212 sage: x = vector(QQ,[1,2,3])
213 sage: y = vector(QQ,[4,5,6])
214 sage: x.inner_product(y)
216 sage: J.from_vector(x).inner_product(J.from_vector(y))
219 The inner product on `S^n` is `<X,Y> = trace(X*Y)`, where
220 multiplication is the usual matrix multiplication in `S^n`,
221 so the inner product of the identity matrix with itself
224 sage: J = RealSymmetricEJA(3)
225 sage: J.one().inner_product(J.one())
228 Likewise, the inner product on `C^n` is `<X,Y> =
229 Re(trace(X*Y))`, where we must necessarily take the real
230 part because the product of Hermitian matrices may not be
233 sage: J = ComplexHermitianEJA(3)
234 sage: J.one().inner_product(J.one())
237 Ditto for the quaternions::
239 sage: J = QuaternionHermitianEJA(3)
240 sage: J.one().inner_product(J.one())
245 Ensure that we can always compute an inner product, and that
246 it gives us back a real number::
248 sage: set_random_seed()
249 sage: J = random_eja()
250 sage: x,y = J.random_elements(2)
251 sage: x.inner_product(y) in RLF
257 raise TypeError("'other' must live in the same algebra")
259 return P
.inner_product(self
, other
)
262 def operator_commutes_with(self
, other
):
264 Return whether or not this element operator-commutes
269 sage: from mjo.eja.eja_algebra import random_eja
273 The definition of a Jordan algebra says that any element
274 operator-commutes with its square::
276 sage: set_random_seed()
277 sage: x = random_eja().random_element()
278 sage: x.operator_commutes_with(x^2)
283 Test Lemma 1 from Chapter III of Koecher::
285 sage: set_random_seed()
286 sage: u,v = random_eja().random_elements(2)
287 sage: lhs = u.operator_commutes_with(u*v)
288 sage: rhs = v.operator_commutes_with(u^2)
292 Test the first polarization identity from my notes, Koecher
293 Chapter III, or from Baes (2.3)::
295 sage: set_random_seed()
296 sage: x,y = random_eja().random_elements(2)
297 sage: Lx = x.operator()
298 sage: Ly = y.operator()
299 sage: Lxx = (x*x).operator()
300 sage: Lxy = (x*y).operator()
301 sage: bool(2*Lx*Lxy + Ly*Lxx == 2*Lxy*Lx + Lxx*Ly)
304 Test the second polarization identity from my notes or from
307 sage: set_random_seed()
308 sage: x,y,z = random_eja().random_elements(3)
309 sage: Lx = x.operator()
310 sage: Ly = y.operator()
311 sage: Lz = z.operator()
312 sage: Lzy = (z*y).operator()
313 sage: Lxy = (x*y).operator()
314 sage: Lxz = (x*z).operator()
315 sage: bool(Lx*Lzy + Lz*Lxy + Ly*Lxz == Lzy*Lx + Lxy*Lz + Lxz*Ly)
318 Test the third polarization identity from my notes or from
321 sage: set_random_seed()
322 sage: u,y,z = random_eja().random_elements(3)
323 sage: Lu = u.operator()
324 sage: Ly = y.operator()
325 sage: Lz = z.operator()
326 sage: Lzy = (z*y).operator()
327 sage: Luy = (u*y).operator()
328 sage: Luz = (u*z).operator()
329 sage: Luyz = (u*(y*z)).operator()
330 sage: lhs = Lu*Lzy + Lz*Luy + Ly*Luz
331 sage: rhs = Luyz + Ly*Lu*Lz + Lz*Lu*Ly
332 sage: bool(lhs == rhs)
336 if not other
in self
.parent():
337 raise TypeError("'other' must live in the same algebra")
346 Return my determinant, the product of my eigenvalues.
350 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
355 sage: J = JordanSpinEJA(2)
356 sage: e0,e1 = J.gens()
357 sage: x = sum( J.gens() )
363 sage: J = JordanSpinEJA(3)
364 sage: e0,e1,e2 = J.gens()
365 sage: x = sum( J.gens() )
371 An element is invertible if and only if its determinant is
374 sage: set_random_seed()
375 sage: x = random_eja().random_element()
376 sage: x.is_invertible() == (x.det() != 0)
379 Ensure that the determinant is multiplicative on an associative
380 subalgebra as in Faraut and Korányi's Proposition II.2.2::
382 sage: set_random_seed()
383 sage: J = random_eja().random_element().subalgebra_generated_by()
384 sage: x,y = J.random_elements(2)
385 sage: (x*y).det() == x.det()*y.det()
391 p
= P
._charpoly
_coeff
(0)
392 # The _charpoly_coeff function already adds the factor of
393 # -1 to ensure that _charpoly_coeff(0) is really what
394 # appears in front of t^{0} in the charpoly. However,
395 # we want (-1)^r times THAT for the determinant.
396 return ((-1)**r
)*p(*self
.to_vector())
401 Return the Jordan-multiplicative inverse of this element.
405 We appeal to the quadratic representation as in Koecher's
406 Theorem 12 in Chapter III, Section 5.
410 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
416 The inverse in the spin factor algebra is given in Alizadeh's
419 sage: set_random_seed()
420 sage: J = JordanSpinEJA.random_instance()
421 sage: x = J.random_element()
422 sage: while not x.is_invertible():
423 ....: x = J.random_element()
424 sage: x_vec = x.to_vector()
426 sage: x_bar = x_vec[1:]
427 sage: coeff = ~(x0^2 - x_bar.inner_product(x_bar))
428 sage: inv_vec = x_vec.parent()([x0] + (-x_bar).list())
429 sage: x_inverse = coeff*inv_vec
430 sage: x.inverse() == J.from_vector(x_inverse)
435 The identity element is its own inverse::
437 sage: set_random_seed()
438 sage: J = random_eja()
439 sage: J.one().inverse() == J.one()
442 If an element has an inverse, it acts like one::
444 sage: set_random_seed()
445 sage: J = random_eja()
446 sage: x = J.random_element()
447 sage: (not x.is_invertible()) or (x.inverse()*x == J.one())
450 The inverse of the inverse is what we started with::
452 sage: set_random_seed()
453 sage: J = random_eja()
454 sage: x = J.random_element()
455 sage: (not x.is_invertible()) or (x.inverse().inverse() == x)
458 The zero element is never invertible::
460 sage: set_random_seed()
461 sage: J = random_eja().zero().inverse()
462 Traceback (most recent call last):
464 ValueError: element is not invertible
466 Proposition II.2.3 in Faraut and Korányi says that the inverse
467 of an element is the inverse of its left-multiplication operator
468 applied to the algebra's identity, when that inverse exists::
470 sage: set_random_seed()
471 sage: J = random_eja()
472 sage: x = J.random_element()
473 sage: (not x.operator().is_invertible()) or (
474 ....: x.operator().inverse()(J.one()) == x.inverse() )
477 Proposition II.2.4 in Faraut and Korányi gives a formula for
478 the inverse based on the characteristic polynomial and the
479 Cayley-Hamilton theorem for Euclidean Jordan algebras::
481 sage: set_random_seed()
482 sage: J = ComplexHermitianEJA(3)
483 sage: x = J.random_element()
484 sage: while not x.is_invertible():
485 ....: x = J.random_element()
487 sage: a = x.characteristic_polynomial().coefficients(sparse=False)
488 sage: expected = (-1)^(r+1)/x.det()
489 sage: expected *= sum( a[i+1]*x^i for i in range(r) )
490 sage: x.inverse() == expected
494 if not self
.is_invertible():
495 raise ValueError("element is not invertible")
497 return (~self
.quadratic_representation())(self
)
500 def is_invertible(self
):
502 Return whether or not this element is invertible.
506 The usual way to do this is to check if the determinant is
507 zero, but we need the characteristic polynomial for the
508 determinant. The minimal polynomial is a lot easier to get,
509 so we use Corollary 2 in Chapter V of Koecher to check
510 whether or not the paren't algebra's zero element is a root
511 of this element's minimal polynomial.
513 Beware that we can't use the superclass method, because it
514 relies on the algebra being associative.
518 sage: from mjo.eja.eja_algebra import random_eja
522 The identity element is always invertible::
524 sage: set_random_seed()
525 sage: J = random_eja()
526 sage: J.one().is_invertible()
529 The zero element is never invertible in a non-trivial algebra::
531 sage: set_random_seed()
532 sage: J = random_eja()
533 sage: (not J.is_trivial()) and J.zero().is_invertible()
538 if self
.parent().is_trivial():
543 # In fact, we only need to know if the constant term is non-zero,
544 # so we can pass in the field's zero element instead.
545 zero
= self
.base_ring().zero()
546 p
= self
.minimal_polynomial()
547 return not (p(zero
) == zero
)
550 def is_nilpotent(self
):
552 Return whether or not some power of this element is zero.
556 We use Theorem 5 in Chapter III of Koecher, which says that
557 an element ``x`` is nilpotent if and only if ``x.operator()``
558 is nilpotent. And it is a basic fact of linear algebra that
559 an operator on an `n`-dimensional space is nilpotent if and
560 only if, when raised to the `n`th power, it equals the zero
561 operator (for example, see Axler Corollary 8.8).
565 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
570 sage: J = JordanSpinEJA(3)
571 sage: x = sum(J.gens())
572 sage: x.is_nilpotent()
577 The identity element is never nilpotent::
579 sage: set_random_seed()
580 sage: random_eja().one().is_nilpotent()
583 The additive identity is always nilpotent::
585 sage: set_random_seed()
586 sage: random_eja().zero().is_nilpotent()
591 zero_operator
= P
.zero().operator()
592 return self
.operator()**P
.dimension() == zero_operator
595 def is_regular(self
):
597 Return whether or not this is a regular element.
601 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
606 The identity element always has degree one, but any element
607 linearly-independent from it is regular::
609 sage: J = JordanSpinEJA(5)
610 sage: J.one().is_regular()
612 sage: e0, e1, e2, e3, e4 = J.gens() # e0 is the identity
613 sage: for x in J.gens():
614 ....: (J.one() + x).is_regular()
623 The zero element should never be regular, unless the parent
624 algebra has dimension one::
626 sage: set_random_seed()
627 sage: J = random_eja()
628 sage: J.dimension() == 1 or not J.zero().is_regular()
631 The unit element isn't regular unless the algebra happens to
632 consist of only its scalar multiples::
634 sage: set_random_seed()
635 sage: J = random_eja()
636 sage: J.dimension() == 1 or not J.one().is_regular()
640 return self
.degree() == self
.parent().rank()
645 Return the degree of this element, which is defined to be
646 the degree of its minimal polynomial.
650 For now, we skip the messy minimal polynomial computation
651 and instead return the dimension of the vector space spanned
652 by the powers of this element. The latter is a bit more
653 straightforward to compute.
657 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
662 sage: J = JordanSpinEJA(4)
663 sage: J.one().degree()
665 sage: e0,e1,e2,e3 = J.gens()
666 sage: (e0 - e1).degree()
669 In the spin factor algebra (of rank two), all elements that
670 aren't multiples of the identity are regular::
672 sage: set_random_seed()
673 sage: J = JordanSpinEJA.random_instance()
674 sage: x = J.random_element()
675 sage: x == x.coefficient(0)*J.one() or x.degree() == 2
680 The zero and unit elements are both of degree one::
682 sage: set_random_seed()
683 sage: J = random_eja()
684 sage: J.zero().degree()
686 sage: J.one().degree()
689 Our implementation agrees with the definition::
691 sage: set_random_seed()
692 sage: x = random_eja().random_element()
693 sage: x.degree() == x.minimal_polynomial().degree()
697 if self
.is_zero() and not self
.parent().is_trivial():
698 # The minimal polynomial of zero in a nontrivial algebra
699 # is "t"; in a trivial algebra it's "1" by convention
700 # (it's an empty product).
702 return self
.subalgebra_generated_by().dimension()
705 def left_matrix(self
):
707 Our parent class defines ``left_matrix`` and ``matrix``
708 methods whose names are misleading. We don't want them.
710 raise NotImplementedError("use operator().matrix() instead")
715 def minimal_polynomial(self
):
717 Return the minimal polynomial of this element,
718 as a function of the variable `t`.
722 We restrict ourselves to the associative subalgebra
723 generated by this element, and then return the minimal
724 polynomial of this element's operator matrix (in that
725 subalgebra). This works by Baes Proposition 2.3.16.
729 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
730 ....: RealSymmetricEJA,
735 The minimal polynomial of the identity and zero elements are
738 sage: set_random_seed()
739 sage: J = random_eja()
740 sage: J.one().minimal_polynomial()
742 sage: J.zero().minimal_polynomial()
745 The degree of an element is (by one definition) the degree
746 of its minimal polynomial::
748 sage: set_random_seed()
749 sage: x = random_eja().random_element()
750 sage: x.degree() == x.minimal_polynomial().degree()
753 The minimal polynomial and the characteristic polynomial coincide
754 and are known (see Alizadeh, Example 11.11) for all elements of
755 the spin factor algebra that aren't scalar multiples of the
756 identity. We require the dimension of the algebra to be at least
757 two here so that said elements actually exist::
759 sage: set_random_seed()
760 sage: n_max = max(2, JordanSpinEJA._max_test_case_size())
761 sage: n = ZZ.random_element(2, n_max)
762 sage: J = JordanSpinEJA(n)
763 sage: y = J.random_element()
764 sage: while y == y.coefficient(0)*J.one():
765 ....: y = J.random_element()
766 sage: y0 = y.to_vector()[0]
767 sage: y_bar = y.to_vector()[1:]
768 sage: actual = y.minimal_polynomial()
769 sage: t = PolynomialRing(J.base_ring(),'t').gen(0)
770 sage: expected = t^2 - 2*y0*t + (y0^2 - norm(y_bar)^2)
771 sage: bool(actual == expected)
774 The minimal polynomial should always kill its element::
776 sage: set_random_seed()
777 sage: x = random_eja().random_element()
778 sage: p = x.minimal_polynomial()
779 sage: x.apply_univariate_polynomial(p)
782 The minimal polynomial is invariant under a change of basis,
783 and in particular, a re-scaling of the basis::
785 sage: set_random_seed()
786 sage: n_max = RealSymmetricEJA._max_test_case_size()
787 sage: n = ZZ.random_element(1, n_max)
788 sage: J1 = RealSymmetricEJA(n,QQ)
789 sage: J2 = RealSymmetricEJA(n,QQ,normalize_basis=False)
790 sage: X = random_matrix(QQ,n)
791 sage: X = X*X.transpose()
794 sage: x1.minimal_polynomial() == x2.minimal_polynomial()
799 # We would generate a zero-dimensional subalgebra
800 # where the minimal polynomial would be constant.
801 # That might be correct, but only if *this* algebra
803 if not self
.parent().is_trivial():
804 # Pretty sure we know what the minimal polynomial of
805 # the zero operator is going to be. This ensures
806 # consistency of e.g. the polynomial variable returned
807 # in the "normal" case without us having to think about it.
808 return self
.operator().minimal_polynomial()
810 A
= self
.subalgebra_generated_by()
811 return A(self
).operator().minimal_polynomial()
815 def natural_representation(self
):
817 Return a more-natural representation of this element.
819 Every finite-dimensional Euclidean Jordan Algebra is a
820 direct sum of five simple algebras, four of which comprise
821 Hermitian matrices. This method returns the original
822 "natural" representation of this element as a Hermitian
823 matrix, if it has one. If not, you get the usual representation.
827 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
828 ....: QuaternionHermitianEJA)
832 sage: J = ComplexHermitianEJA(3)
835 sage: J.one().natural_representation()
845 sage: J = QuaternionHermitianEJA(3)
848 sage: J.one().natural_representation()
849 [1 0 0 0 0 0 0 0 0 0 0 0]
850 [0 1 0 0 0 0 0 0 0 0 0 0]
851 [0 0 1 0 0 0 0 0 0 0 0 0]
852 [0 0 0 1 0 0 0 0 0 0 0 0]
853 [0 0 0 0 1 0 0 0 0 0 0 0]
854 [0 0 0 0 0 1 0 0 0 0 0 0]
855 [0 0 0 0 0 0 1 0 0 0 0 0]
856 [0 0 0 0 0 0 0 1 0 0 0 0]
857 [0 0 0 0 0 0 0 0 1 0 0 0]
858 [0 0 0 0 0 0 0 0 0 1 0 0]
859 [0 0 0 0 0 0 0 0 0 0 1 0]
860 [0 0 0 0 0 0 0 0 0 0 0 1]
863 B
= self
.parent().natural_basis()
864 W
= self
.parent().natural_basis_space()
865 return W
.linear_combination(izip(B
,self
.to_vector()))
870 The norm of this element with respect to :meth:`inner_product`.
874 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
875 ....: RealCartesianProductEJA)
879 sage: J = RealCartesianProductEJA(2)
880 sage: x = sum(J.gens())
886 sage: J = JordanSpinEJA(4)
887 sage: x = sum(J.gens())
892 return self
.inner_product(self
).sqrt()
897 Return the left-multiplication-by-this-element
898 operator on the ambient algebra.
902 sage: from mjo.eja.eja_algebra import random_eja
906 sage: set_random_seed()
907 sage: J = random_eja()
908 sage: x,y = J.random_elements(2)
909 sage: x.operator()(y) == x*y
911 sage: y.operator()(x) == x*y
916 left_mult_by_self
= lambda y
: self
*y
917 L
= P
.module_morphism(function
=left_mult_by_self
, codomain
=P
)
918 return FiniteDimensionalEuclideanJordanAlgebraOperator(
924 def quadratic_representation(self
, other
=None):
926 Return the quadratic representation of this element.
930 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
935 The explicit form in the spin factor algebra is given by
936 Alizadeh's Example 11.12::
938 sage: set_random_seed()
939 sage: x = JordanSpinEJA.random_instance().random_element()
940 sage: x_vec = x.to_vector()
941 sage: n = x_vec.degree()
943 sage: x_bar = x_vec[1:]
944 sage: A = matrix(QQ, 1, [x_vec.inner_product(x_vec)])
945 sage: B = 2*x0*x_bar.row()
946 sage: C = 2*x0*x_bar.column()
947 sage: D = matrix.identity(QQ, n-1)
948 sage: D = (x0^2 - x_bar.inner_product(x_bar))*D
949 sage: D = D + 2*x_bar.tensor_product(x_bar)
950 sage: Q = matrix.block(2,2,[A,B,C,D])
951 sage: Q == x.quadratic_representation().matrix()
954 Test all of the properties from Theorem 11.2 in Alizadeh::
956 sage: set_random_seed()
957 sage: J = random_eja()
958 sage: x,y = J.random_elements(2)
959 sage: Lx = x.operator()
960 sage: Lxx = (x*x).operator()
961 sage: Qx = x.quadratic_representation()
962 sage: Qy = y.quadratic_representation()
963 sage: Qxy = x.quadratic_representation(y)
964 sage: Qex = J.one().quadratic_representation(x)
965 sage: n = ZZ.random_element(10)
966 sage: Qxn = (x^n).quadratic_representation()
970 sage: 2*Qxy == (x+y).quadratic_representation() - Qx - Qy
973 Property 2 (multiply on the right for :trac:`28272`):
975 sage: alpha = J.base_ring().random_element()
976 sage: (alpha*x).quadratic_representation() == Qx*(alpha^2)
981 sage: not x.is_invertible() or ( Qx(x.inverse()) == x )
984 sage: not x.is_invertible() or (
987 ....: x.inverse().quadratic_representation() )
990 sage: Qxy(J.one()) == x*y
995 sage: not x.is_invertible() or (
996 ....: x.quadratic_representation(x.inverse())*Qx
997 ....: == Qx*x.quadratic_representation(x.inverse()) )
1000 sage: not x.is_invertible() or (
1001 ....: x.quadratic_representation(x.inverse())*Qx
1003 ....: 2*Lx*Qex - Qx )
1006 sage: 2*Lx*Qex - Qx == Lxx
1011 sage: Qy(x).quadratic_representation() == Qy*Qx*Qy
1021 sage: not x.is_invertible() or (
1022 ....: Qx*x.inverse().operator() == Lx )
1027 sage: not x.operator_commutes_with(y) or (
1028 ....: Qx(y)^n == Qxn(y^n) )
1034 elif not other
in self
.parent():
1035 raise TypeError("'other' must live in the same algebra")
1038 M
= other
.operator()
1039 return ( L
*M
+ M
*L
- (self
*other
).operator() )
1043 def spectral_decomposition(self
):
1045 Return the unique spectral decomposition of this element.
1049 Following Faraut and Korányi's Theorem III.1.1, we restrict this
1050 element's left-multiplication-by operator to the subalgebra it
1051 generates. We then compute the spectral decomposition of that
1052 operator, and the spectral projectors we get back must be the
1053 left-multiplication-by operators for the idempotents we
1054 seek. Thus applying them to the identity element gives us those
1057 Since the eigenvalues are required to be distinct, we take
1058 the spectral decomposition of the zero element to be zero
1059 times the identity element of the algebra (which is idempotent,
1064 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1068 The spectral decomposition of the identity is ``1`` times itself,
1069 and the spectral decomposition of zero is ``0`` times the identity::
1071 sage: J = RealSymmetricEJA(3,AA)
1074 sage: J.one().spectral_decomposition()
1076 sage: J.zero().spectral_decomposition()
1081 sage: J = RealSymmetricEJA(4,AA)
1082 sage: x = sum(J.gens())
1083 sage: sd = x.spectral_decomposition()
1088 sage: c0.inner_product(c1) == 0
1090 sage: c0.is_idempotent()
1092 sage: c1.is_idempotent()
1094 sage: c0 + c1 == J.one()
1096 sage: l0*c0 + l1*c1 == x
1101 A
= self
.subalgebra_generated_by(orthonormalize_basis
=True)
1103 for (evalue
, proj
) in A(self
).operator().spectral_decomposition():
1104 result
.append( (evalue
, proj(A
.one()).superalgebra_element()) )
1107 def subalgebra_generated_by(self
, orthonormalize_basis
=False):
1109 Return the associative subalgebra of the parent EJA generated
1112 Since our parent algebra is unital, we want "subalgebra" to mean
1113 "unital subalgebra" as well; thus the subalgebra that an element
1114 generates will itself be a Euclidean Jordan algebra after
1115 restricting the algebra operations appropriately. This is the
1116 subalgebra that Faraut and Korányi work with in section II.2, for
1121 sage: from mjo.eja.eja_algebra import random_eja
1125 This subalgebra, being composed of only powers, is associative::
1127 sage: set_random_seed()
1128 sage: x0 = random_eja().random_element()
1129 sage: A = x0.subalgebra_generated_by()
1130 sage: x,y,z = A.random_elements(3)
1131 sage: (x*y)*z == x*(y*z)
1134 Squaring in the subalgebra should work the same as in
1137 sage: set_random_seed()
1138 sage: x = random_eja().random_element()
1139 sage: A = x.subalgebra_generated_by()
1140 sage: A(x^2) == A(x)*A(x)
1143 By definition, the subalgebra generated by the zero element is the
1144 one-dimensional algebra generated by the identity element::
1146 sage: set_random_seed()
1147 sage: A = random_eja().zero().subalgebra_generated_by()
1152 return FiniteDimensionalEuclideanJordanElementSubalgebra(self
, orthonormalize_basis
)
1155 def subalgebra_idempotent(self
):
1157 Find an idempotent in the associative subalgebra I generate
1158 using Proposition 2.3.5 in Baes.
1162 sage: from mjo.eja.eja_algebra import random_eja
1166 sage: set_random_seed()
1167 sage: J = random_eja()
1168 sage: x = J.random_element()
1169 sage: while x.is_nilpotent():
1170 ....: x = J.random_element()
1171 sage: c = x.subalgebra_idempotent()
1176 if self
.is_nilpotent():
1177 raise ValueError("this only works with non-nilpotent elements!")
1179 J
= self
.subalgebra_generated_by()
1182 # The image of the matrix of left-u^m-multiplication
1183 # will be minimal for some natural number s...
1185 minimal_dim
= J
.dimension()
1186 for i
in xrange(1, minimal_dim
):
1187 this_dim
= (u
**i
).operator().matrix().image().dimension()
1188 if this_dim
< minimal_dim
:
1189 minimal_dim
= this_dim
1192 # Now minimal_matrix should correspond to the smallest
1193 # non-zero subspace in Baes's (or really, Koecher's)
1196 # However, we need to restrict the matrix to work on the
1197 # subspace... or do we? Can't we just solve, knowing that
1198 # A(c) = u^(s+1) should have a solution in the big space,
1201 # Beware, solve_right() means that we're using COLUMN vectors.
1202 # Our FiniteDimensionalAlgebraElement superclass uses rows.
1204 A
= u_next
.operator().matrix()
1205 c
= J
.from_vector(A
.solve_right(u_next
.to_vector()))
1207 # Now c is the idempotent we want, but it still lives in the subalgebra.
1208 return c
.superalgebra_element()
1213 Return my trace, the sum of my eigenvalues.
1217 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1218 ....: RealCartesianProductEJA,
1223 sage: J = JordanSpinEJA(3)
1224 sage: x = sum(J.gens())
1230 sage: J = RealCartesianProductEJA(5)
1231 sage: J.one().trace()
1236 The trace of an element is a real number::
1238 sage: set_random_seed()
1239 sage: J = random_eja()
1240 sage: J.random_element().trace() in RLF
1246 p
= P
._charpoly
_coeff
(r
-1)
1247 # The _charpoly_coeff function already adds the factor of
1248 # -1 to ensure that _charpoly_coeff(r-1) is really what
1249 # appears in front of t^{r-1} in the charpoly. However,
1250 # we want the negative of THAT for the trace.
1251 return -p(*self
.to_vector())
1254 def trace_inner_product(self
, other
):
1256 Return the trace inner product of myself and ``other``.
1260 sage: from mjo.eja.eja_algebra import random_eja
1264 The trace inner product is commutative, bilinear, and associative::
1266 sage: set_random_seed()
1267 sage: J = random_eja()
1268 sage: x,y,z = J.random_elements(3)
1270 sage: x.trace_inner_product(y) == y.trace_inner_product(x)
1273 sage: a = J.base_ring().random_element();
1274 sage: actual = (a*(x+z)).trace_inner_product(y)
1275 sage: expected = ( a*x.trace_inner_product(y) +
1276 ....: a*z.trace_inner_product(y) )
1277 sage: actual == expected
1279 sage: actual = x.trace_inner_product(a*(y+z))
1280 sage: expected = ( a*x.trace_inner_product(y) +
1281 ....: a*x.trace_inner_product(z) )
1282 sage: actual == expected
1285 sage: (x*y).trace_inner_product(z) == y.trace_inner_product(x*z)
1289 if not other
in self
.parent():
1290 raise TypeError("'other' must live in the same algebra")
1292 return (self
*other
).trace()
1295 def trace_norm(self
):
1297 The norm of this element with respect to :meth:`trace_inner_product`.
1301 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1302 ....: RealCartesianProductEJA)
1306 sage: J = RealCartesianProductEJA(2)
1307 sage: x = sum(J.gens())
1308 sage: x.trace_norm()
1313 sage: J = JordanSpinEJA(4)
1314 sage: x = sum(J.gens())
1315 sage: x.trace_norm()
1319 return self
.trace_inner_product(self
).sqrt()