2 Representations and constructions for Euclidean Jordan algebras.
4 A Euclidean Jordan algebra is a Jordan algebra that has some
7 1. It is finite-dimensional.
8 2. Its scalar field is the real numbers.
9 3a. An inner product is defined on it, and...
10 3b. That inner product is compatible with the Jordan product
11 in the sense that `<x*y,z> = <y,x*z>` for all elements
12 `x,y,z` in the algebra.
14 Every Euclidean Jordan algebra is formally-real: for any two elements
15 `x` and `y` in the algebra, `x^{2} + y^{2} = 0` implies that `x = y =
16 0`. Conversely, every finite-dimensional formally-real Jordan algebra
17 can be made into a Euclidean Jordan algebra with an appropriate choice
20 Formally-real Jordan algebras were originally studied as a framework
21 for quantum mechanics. Today, Euclidean Jordan algebras are crucial in
22 symmetric cone optimization, since every symmetric cone arises as the
23 cone of squares in some Euclidean Jordan algebra.
25 It is known that every Euclidean Jordan algebra decomposes into an
26 orthogonal direct sum (essentially, a Cartesian product) of simple
27 algebras, and that moreover, up to Jordan-algebra isomorphism, there
28 are only five families of simple algebras. We provide constructions
29 for these simple algebras:
31 * :class:`BilinearFormEJA`
32 * :class:`RealSymmetricEJA`
33 * :class:`ComplexHermitianEJA`
34 * :class:`QuaternionHermitianEJA`
36 Missing from this list is the algebra of three-by-three octononion
37 Hermitian matrices, as there is (as of yet) no implementation of the
38 octonions in SageMath. In addition to these, we provide two other
39 example constructions,
41 * :class:`HadamardEJA`
44 The Jordan spin algebra is a bilinear form algebra where the bilinear
45 form is the identity. The Hadamard EJA is simply a Cartesian product
46 of one-dimensional spin algebras. And last but not least, the trivial
47 EJA is exactly what you think. Cartesian products of these are also
48 supported using the usual ``cartesian_product()`` function; as a
49 result, we support (up to isomorphism) all Euclidean Jordan algebras
50 that don't involve octonions.
54 sage: from mjo.eja.eja_algebra import random_eja
59 Euclidean Jordan algebra of dimension...
62 from itertools
import repeat
64 from sage
.algebras
.quatalg
.quaternion_algebra
import QuaternionAlgebra
65 from sage
.categories
.magmatic_algebras
import MagmaticAlgebras
66 from sage
.categories
.sets_cat
import cartesian_product
67 from sage
.combinat
.free_module
import (CombinatorialFreeModule
,
68 CombinatorialFreeModule_CartesianProduct
)
69 from sage
.matrix
.constructor
import matrix
70 from sage
.matrix
.matrix_space
import MatrixSpace
71 from sage
.misc
.cachefunc
import cached_method
72 from sage
.misc
.table
import table
73 from sage
.modules
.free_module
import FreeModule
, VectorSpace
74 from sage
.rings
.all
import (ZZ
, QQ
, AA
, QQbar
, RR
, RLF
, CLF
,
77 from mjo
.eja
.eja_element
import FiniteDimensionalEJAElement
78 from mjo
.eja
.eja_operator
import FiniteDimensionalEJAOperator
79 from mjo
.eja
.eja_utils
import _all2list
, _mat2vec
81 class FiniteDimensionalEJA(CombinatorialFreeModule
):
83 A finite-dimensional Euclidean Jordan algebra.
87 - ``basis`` -- a tuple; a tuple of basis elements in "matrix
88 form," which must be the same form as the arguments to
89 ``jordan_product`` and ``inner_product``. In reality, "matrix
90 form" can be either vectors, matrices, or a Cartesian product
91 (ordered tuple) of vectors or matrices. All of these would
92 ideally be vector spaces in sage with no special-casing
93 needed; but in reality we turn vectors into column-matrices
94 and Cartesian products `(a,b)` into column matrices
95 `(a,b)^{T}` after converting `a` and `b` themselves.
97 - ``jordan_product`` -- a function; afunction of two ``basis``
98 elements (in matrix form) that returns their jordan product,
99 also in matrix form; this will be applied to ``basis`` to
100 compute a multiplication table for the algebra.
102 - ``inner_product`` -- a function; a function of two ``basis``
103 elements (in matrix form) that returns their inner
104 product. This will be applied to ``basis`` to compute an
105 inner-product table (basically a matrix) for this algebra.
107 - ``field`` -- a subfield of the reals (default: ``AA``); the scalar
108 field for the algebra.
110 - ``orthonormalize`` -- boolean (default: ``True``); whether or
111 not to orthonormalize the basis. Doing so is expensive and
112 generally rules out using the rationals as your ``field``, but
113 is required for spectral decompositions.
116 Element
= FiniteDimensionalEJAElement
125 cartesian_product
=False,
130 # Keep track of whether or not the matrix basis consists of
131 # tuples, since we need special cases for them damned near
132 # everywhere. This is INDEPENDENT of whether or not the
133 # algebra is a cartesian product, since a subalgebra of a
134 # cartesian product will have a basis of tuples, but will not
135 # in general itself be a cartesian product algebra.
136 self
._matrix
_basis
_is
_cartesian
= False
139 if hasattr(basis
[0], 'cartesian_factors'):
140 self
._matrix
_basis
_is
_cartesian
= True
143 if not field
.is_subring(RR
):
144 # Note: this does return true for the real algebraic
145 # field, the rationals, and any quadratic field where
146 # we've specified a real embedding.
147 raise ValueError("scalar field is not real")
149 # If the basis given to us wasn't over the field that it's
150 # supposed to be over, fix that. Or, you know, crash.
151 if not cartesian_product
:
152 # The field for a cartesian product algebra comes from one
153 # of its factors and is the same for all factors, so
154 # there's no need to "reapply" it on product algebras.
155 if self
._matrix
_basis
_is
_cartesian
:
156 # OK since if n == 0, the basis does not consist of tuples.
157 P
= basis
[0].parent()
158 basis
= tuple( P(tuple(b_i
.change_ring(field
) for b_i
in b
))
161 basis
= tuple( b
.change_ring(field
) for b
in basis
)
165 # Check commutativity of the Jordan and inner-products.
166 # This has to be done before we build the multiplication
167 # and inner-product tables/matrices, because we take
168 # advantage of symmetry in the process.
169 if not all( jordan_product(bi
,bj
) == jordan_product(bj
,bi
)
172 raise ValueError("Jordan product is not commutative")
174 if not all( inner_product(bi
,bj
) == inner_product(bj
,bi
)
177 raise ValueError("inner-product is not commutative")
180 category
= MagmaticAlgebras(field
).FiniteDimensional()
181 category
= category
.WithBasis().Unital().Commutative()
184 # Element subalgebras can take advantage of this.
185 category
= category
.Associative()
186 if cartesian_product
:
187 category
= category
.CartesianProducts()
189 # Call the superclass constructor so that we can use its from_vector()
190 # method to build our multiplication table.
191 CombinatorialFreeModule
.__init
__(self
,
198 # Now comes all of the hard work. We'll be constructing an
199 # ambient vector space V that our (vectorized) basis lives in,
200 # as well as a subspace W of V spanned by those (vectorized)
201 # basis elements. The W-coordinates are the coefficients that
202 # we see in things like x = 1*e1 + 2*e2.
207 degree
= len(_all2list(basis
[0]))
209 # Build an ambient space that fits our matrix basis when
210 # written out as "long vectors."
211 V
= VectorSpace(field
, degree
)
213 # The matrix that will hole the orthonormal -> unorthonormal
214 # coordinate transformation.
215 self
._deortho
_matrix
= None
218 # Save a copy of the un-orthonormalized basis for later.
219 # Convert it to ambient V (vector) coordinates while we're
220 # at it, because we'd have to do it later anyway.
221 deortho_vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
223 from mjo
.eja
.eja_utils
import gram_schmidt
224 basis
= tuple(gram_schmidt(basis
, inner_product
))
226 # Save the (possibly orthonormalized) matrix basis for
228 self
._matrix
_basis
= basis
230 # Now create the vector space for the algebra, which will have
231 # its own set of non-ambient coordinates (in terms of the
233 vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
234 W
= V
.span_of_basis( vector_basis
, check
=check_axioms
)
237 # Now "W" is the vector space of our algebra coordinates. The
238 # variables "X1", "X2",... refer to the entries of vectors in
239 # W. Thus to convert back and forth between the orthonormal
240 # coordinates and the given ones, we need to stick the original
242 U
= V
.span_of_basis( deortho_vector_basis
, check
=check_axioms
)
243 self
._deortho
_matrix
= matrix( U
.coordinate_vector(q
)
244 for q
in vector_basis
)
247 # Now we actually compute the multiplication and inner-product
248 # tables/matrices using the possibly-orthonormalized basis.
249 self
._inner
_product
_matrix
= matrix
.identity(field
, n
)
250 self
._multiplication
_table
= [ [0 for j
in range(i
+1)]
253 # Note: the Jordan and inner-products are defined in terms
254 # of the ambient basis. It's important that their arguments
255 # are in ambient coordinates as well.
258 # ortho basis w.r.t. ambient coords
262 # The jordan product returns a matrixy answer, so we
263 # have to convert it to the algebra coordinates.
264 elt
= jordan_product(q_i
, q_j
)
265 elt
= W
.coordinate_vector(V(_all2list(elt
)))
266 self
._multiplication
_table
[i
][j
] = self
.from_vector(elt
)
268 if not orthonormalize
:
269 # If we're orthonormalizing the basis with respect
270 # to an inner-product, then the inner-product
271 # matrix with respect to the resulting basis is
272 # just going to be the identity.
273 ip
= inner_product(q_i
, q_j
)
274 self
._inner
_product
_matrix
[i
,j
] = ip
275 self
._inner
_product
_matrix
[j
,i
] = ip
277 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
278 self
._inner
_product
_matrix
.set_immutable()
281 if not self
._is
_jordanian
():
282 raise ValueError("Jordan identity does not hold")
283 if not self
._inner
_product
_is
_associative
():
284 raise ValueError("inner product is not associative")
287 def _coerce_map_from_base_ring(self
):
289 Disable the map from the base ring into the algebra.
291 Performing a nonsense conversion like this automatically
292 is counterpedagogical. The fallback is to try the usual
293 element constructor, which should also fail.
297 sage: from mjo.eja.eja_algebra import random_eja
301 sage: set_random_seed()
302 sage: J = random_eja()
304 Traceback (most recent call last):
306 ValueError: not an element of this algebra
312 def product_on_basis(self
, i
, j
):
314 Returns the Jordan product of the `i` and `j`th basis elements.
316 This completely defines the Jordan product on the algebra, and
317 is used direclty by our superclass machinery to implement
322 sage: from mjo.eja.eja_algebra import random_eja
326 sage: set_random_seed()
327 sage: J = random_eja()
328 sage: n = J.dimension()
331 sage: ei_ej = J.zero()*J.zero()
333 ....: i = ZZ.random_element(n)
334 ....: j = ZZ.random_element(n)
335 ....: ei = J.gens()[i]
336 ....: ej = J.gens()[j]
337 ....: ei_ej = J.product_on_basis(i,j)
342 # We only stored the lower-triangular portion of the
343 # multiplication table.
345 return self
._multiplication
_table
[i
][j
]
347 return self
._multiplication
_table
[j
][i
]
349 def inner_product(self
, x
, y
):
351 The inner product associated with this Euclidean Jordan algebra.
353 Defaults to the trace inner product, but can be overridden by
354 subclasses if they are sure that the necessary properties are
359 sage: from mjo.eja.eja_algebra import (random_eja,
361 ....: BilinearFormEJA)
365 Our inner product is "associative," which means the following for
366 a symmetric bilinear form::
368 sage: set_random_seed()
369 sage: J = random_eja()
370 sage: x,y,z = J.random_elements(3)
371 sage: (x*y).inner_product(z) == y.inner_product(x*z)
376 Ensure that this is the usual inner product for the algebras
379 sage: set_random_seed()
380 sage: J = HadamardEJA.random_instance()
381 sage: x,y = J.random_elements(2)
382 sage: actual = x.inner_product(y)
383 sage: expected = x.to_vector().inner_product(y.to_vector())
384 sage: actual == expected
387 Ensure that this is one-half of the trace inner-product in a
388 BilinearFormEJA that isn't just the reals (when ``n`` isn't
389 one). This is in Faraut and Koranyi, and also my "On the
392 sage: set_random_seed()
393 sage: J = BilinearFormEJA.random_instance()
394 sage: n = J.dimension()
395 sage: x = J.random_element()
396 sage: y = J.random_element()
397 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
401 B
= self
._inner
_product
_matrix
402 return (B
*x
.to_vector()).inner_product(y
.to_vector())
405 def is_associative(self
):
407 Return whether or not this algebra's Jordan product is associative.
411 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
415 sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False)
416 sage: J.is_associative()
418 sage: x = sum(J.gens())
419 sage: A = x.subalgebra_generated_by(orthonormalize=False)
420 sage: A.is_associative()
424 return "Associative" in self
.category().axioms()
426 def _is_commutative(self
):
428 Whether or not this algebra's multiplication table is commutative.
430 This method should of course always return ``True``, unless
431 this algebra was constructed with ``check_axioms=False`` and
432 passed an invalid multiplication table.
434 return all( self
.product_on_basis(i
,j
) == self
.product_on_basis(i
,j
)
435 for i
in range(self
.dimension())
436 for j
in range(self
.dimension()) )
438 def _is_jordanian(self
):
440 Whether or not this algebra's multiplication table respects the
441 Jordan identity `(x^{2})(xy) = x(x^{2}y)`.
443 We only check one arrangement of `x` and `y`, so for a
444 ``True`` result to be truly true, you should also check
445 :meth:`_is_commutative`. This method should of course always
446 return ``True``, unless this algebra was constructed with
447 ``check_axioms=False`` and passed an invalid multiplication table.
449 return all( (self
.gens()[i
]**2)*(self
.gens()[i
]*self
.gens()[j
])
451 (self
.gens()[i
])*((self
.gens()[i
]**2)*self
.gens()[j
])
452 for i
in range(self
.dimension())
453 for j
in range(self
.dimension()) )
455 def _jordan_product_is_associative(self
):
457 Return whether or not this algebra's Jordan product is
458 associative; that is, whether or not `x*(y*z) = (x*y)*z`
461 This method should agree with :meth:`is_associative` unless
462 you lied about the value of the ``associative`` parameter
463 when you constructed the algebra.
467 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
468 ....: ComplexHermitianEJA,
469 ....: QuaternionHermitianEJA)
473 sage: J = RealSymmetricEJA(4, orthonormalize=False)
474 sage: J._jordan_product_is_associative()
476 sage: x = sum(J.gens())
477 sage: A = x.subalgebra_generated_by()
478 sage: A._jordan_product_is_associative()
483 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
484 sage: J._jordan_product_is_associative()
486 sage: x = sum(J.gens())
487 sage: A = x.subalgebra_generated_by(orthonormalize=False)
488 sage: A._jordan_product_is_associative()
493 sage: J = QuaternionHermitianEJA(2)
494 sage: J._jordan_product_is_associative()
496 sage: x = sum(J.gens())
497 sage: A = x.subalgebra_generated_by()
498 sage: A._jordan_product_is_associative()
504 # Used to check whether or not something is zero.
507 # I don't know of any examples that make this magnitude
508 # necessary because I don't know how to make an
509 # associative algebra when the element subalgebra
510 # construction is unreliable (as it is over RDF; we can't
511 # find the degree of an element because we can't compute
512 # the rank of a matrix). But even multiplication of floats
513 # is non-associative, so *some* epsilon is needed... let's
514 # just take the one from _inner_product_is_associative?
517 for i
in range(self
.dimension()):
518 for j
in range(self
.dimension()):
519 for k
in range(self
.dimension()):
523 diff
= (x
*y
)*z
- x
*(y
*z
)
525 if diff
.norm() > epsilon
:
530 def _inner_product_is_associative(self
):
532 Return whether or not this algebra's inner product `B` is
533 associative; that is, whether or not `B(xy,z) = B(x,yz)`.
535 This method should of course always return ``True``, unless
536 this algebra was constructed with ``check_axioms=False`` and
537 passed an invalid Jordan or inner-product.
541 # Used to check whether or not something is zero.
544 # This choice is sufficient to allow the construction of
545 # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
548 for i
in range(self
.dimension()):
549 for j
in range(self
.dimension()):
550 for k
in range(self
.dimension()):
554 diff
= (x
*y
).inner_product(z
) - x
.inner_product(y
*z
)
556 if diff
.abs() > epsilon
:
561 def _element_constructor_(self
, elt
):
563 Construct an element of this algebra from its vector or matrix
566 This gets called only after the parent element _call_ method
567 fails to find a coercion for the argument.
571 sage: from mjo.eja.eja_algebra import (random_eja,
574 ....: RealSymmetricEJA)
578 The identity in `S^n` is converted to the identity in the EJA::
580 sage: J = RealSymmetricEJA(3)
581 sage: I = matrix.identity(QQ,3)
582 sage: J(I) == J.one()
585 This skew-symmetric matrix can't be represented in the EJA::
587 sage: J = RealSymmetricEJA(3)
588 sage: A = matrix(QQ,3, lambda i,j: i-j)
590 Traceback (most recent call last):
592 ValueError: not an element of this algebra
594 Tuples work as well, provided that the matrix basis for the
595 algebra consists of them::
597 sage: J1 = HadamardEJA(3)
598 sage: J2 = RealSymmetricEJA(2)
599 sage: J = cartesian_product([J1,J2])
600 sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
605 Ensure that we can convert any element back and forth
606 faithfully between its matrix and algebra representations::
608 sage: set_random_seed()
609 sage: J = random_eja()
610 sage: x = J.random_element()
611 sage: J(x.to_matrix()) == x
614 We cannot coerce elements between algebras just because their
615 matrix representations are compatible::
617 sage: J1 = HadamardEJA(3)
618 sage: J2 = JordanSpinEJA(3)
620 Traceback (most recent call last):
622 ValueError: not an element of this algebra
624 Traceback (most recent call last):
626 ValueError: not an element of this algebra
628 msg
= "not an element of this algebra"
629 if elt
in self
.base_ring():
630 # Ensure that no base ring -> algebra coercion is performed
631 # by this method. There's some stupidity in sage that would
632 # otherwise propagate to this method; for example, sage thinks
633 # that the integer 3 belongs to the space of 2-by-2 matrices.
634 raise ValueError(msg
)
637 # Try to convert a vector into a column-matrix...
639 except (AttributeError, TypeError):
640 # and ignore failure, because we weren't really expecting
641 # a vector as an argument anyway.
644 if elt
not in self
.matrix_space():
645 raise ValueError(msg
)
647 # Thanks for nothing! Matrix spaces aren't vector spaces in
648 # Sage, so we have to figure out its matrix-basis coordinates
649 # ourselves. We use the basis space's ring instead of the
650 # element's ring because the basis space might be an algebraic
651 # closure whereas the base ring of the 3-by-3 identity matrix
652 # could be QQ instead of QQbar.
654 # And, we also have to handle Cartesian product bases (when
655 # the matrix basis consists of tuples) here. The "good news"
656 # is that we're already converting everything to long vectors,
657 # and that strategy works for tuples as well.
659 # We pass check=False because the matrix basis is "guaranteed"
660 # to be linearly independent... right? Ha ha.
662 V
= VectorSpace(self
.base_ring(), len(elt
))
663 W
= V
.span_of_basis( (V(_all2list(s
)) for s
in self
.matrix_basis()),
667 coords
= W
.coordinate_vector(V(elt
))
668 except ArithmeticError: # vector is not in free module
669 raise ValueError(msg
)
671 return self
.from_vector(coords
)
675 Return a string representation of ``self``.
679 sage: from mjo.eja.eja_algebra import JordanSpinEJA
683 Ensure that it says what we think it says::
685 sage: JordanSpinEJA(2, field=AA)
686 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
687 sage: JordanSpinEJA(3, field=RDF)
688 Euclidean Jordan algebra of dimension 3 over Real Double Field
691 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
692 return fmt
.format(self
.dimension(), self
.base_ring())
696 def characteristic_polynomial_of(self
):
698 Return the algebra's "characteristic polynomial of" function,
699 which is itself a multivariate polynomial that, when evaluated
700 at the coordinates of some algebra element, returns that
701 element's characteristic polynomial.
703 The resulting polynomial has `n+1` variables, where `n` is the
704 dimension of this algebra. The first `n` variables correspond to
705 the coordinates of an algebra element: when evaluated at the
706 coordinates of an algebra element with respect to a certain
707 basis, the result is a univariate polynomial (in the one
708 remaining variable ``t``), namely the characteristic polynomial
713 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
717 The characteristic polynomial in the spin algebra is given in
718 Alizadeh, Example 11.11::
720 sage: J = JordanSpinEJA(3)
721 sage: p = J.characteristic_polynomial_of(); p
722 X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
723 sage: xvec = J.one().to_vector()
727 By definition, the characteristic polynomial is a monic
728 degree-zero polynomial in a rank-zero algebra. Note that
729 Cayley-Hamilton is indeed satisfied since the polynomial
730 ``1`` evaluates to the identity element of the algebra on
733 sage: J = TrivialEJA()
734 sage: J.characteristic_polynomial_of()
741 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
742 a
= self
._charpoly
_coefficients
()
744 # We go to a bit of trouble here to reorder the
745 # indeterminates, so that it's easier to evaluate the
746 # characteristic polynomial at x's coordinates and get back
747 # something in terms of t, which is what we want.
748 S
= PolynomialRing(self
.base_ring(),'t')
752 S
= PolynomialRing(S
, R
.variable_names())
755 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
757 def coordinate_polynomial_ring(self
):
759 The multivariate polynomial ring in which this algebra's
760 :meth:`characteristic_polynomial_of` lives.
764 sage: from mjo.eja.eja_algebra import (HadamardEJA,
765 ....: RealSymmetricEJA)
769 sage: J = HadamardEJA(2)
770 sage: J.coordinate_polynomial_ring()
771 Multivariate Polynomial Ring in X1, X2...
772 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
773 sage: J.coordinate_polynomial_ring()
774 Multivariate Polynomial Ring in X1, X2, X3, X4, X5, X6...
777 var_names
= tuple( "X%d" % z
for z
in range(1, self
.dimension()+1) )
778 return PolynomialRing(self
.base_ring(), var_names
)
780 def inner_product(self
, x
, y
):
782 The inner product associated with this Euclidean Jordan algebra.
784 Defaults to the trace inner product, but can be overridden by
785 subclasses if they are sure that the necessary properties are
790 sage: from mjo.eja.eja_algebra import (random_eja,
792 ....: BilinearFormEJA)
796 Our inner product is "associative," which means the following for
797 a symmetric bilinear form::
799 sage: set_random_seed()
800 sage: J = random_eja()
801 sage: x,y,z = J.random_elements(3)
802 sage: (x*y).inner_product(z) == y.inner_product(x*z)
807 Ensure that this is the usual inner product for the algebras
810 sage: set_random_seed()
811 sage: J = HadamardEJA.random_instance()
812 sage: x,y = J.random_elements(2)
813 sage: actual = x.inner_product(y)
814 sage: expected = x.to_vector().inner_product(y.to_vector())
815 sage: actual == expected
818 Ensure that this is one-half of the trace inner-product in a
819 BilinearFormEJA that isn't just the reals (when ``n`` isn't
820 one). This is in Faraut and Koranyi, and also my "On the
823 sage: set_random_seed()
824 sage: J = BilinearFormEJA.random_instance()
825 sage: n = J.dimension()
826 sage: x = J.random_element()
827 sage: y = J.random_element()
828 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
831 B
= self
._inner
_product
_matrix
832 return (B
*x
.to_vector()).inner_product(y
.to_vector())
835 def is_trivial(self
):
837 Return whether or not this algebra is trivial.
839 A trivial algebra contains only the zero element.
843 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
848 sage: J = ComplexHermitianEJA(3)
854 sage: J = TrivialEJA()
859 return self
.dimension() == 0
862 def multiplication_table(self
):
864 Return a visual representation of this algebra's multiplication
865 table (on basis elements).
869 sage: from mjo.eja.eja_algebra import JordanSpinEJA
873 sage: J = JordanSpinEJA(4)
874 sage: J.multiplication_table()
875 +----++----+----+----+----+
876 | * || e0 | e1 | e2 | e3 |
877 +====++====+====+====+====+
878 | e0 || e0 | e1 | e2 | e3 |
879 +----++----+----+----+----+
880 | e1 || e1 | e0 | 0 | 0 |
881 +----++----+----+----+----+
882 | e2 || e2 | 0 | e0 | 0 |
883 +----++----+----+----+----+
884 | e3 || e3 | 0 | 0 | e0 |
885 +----++----+----+----+----+
889 # Prepend the header row.
890 M
= [["*"] + list(self
.gens())]
892 # And to each subsequent row, prepend an entry that belongs to
893 # the left-side "header column."
894 M
+= [ [self
.gens()[i
]] + [ self
.product_on_basis(i
,j
)
898 return table(M
, header_row
=True, header_column
=True, frame
=True)
901 def matrix_basis(self
):
903 Return an (often more natural) representation of this algebras
904 basis as an ordered tuple of matrices.
906 Every finite-dimensional Euclidean Jordan Algebra is a, up to
907 Jordan isomorphism, a direct sum of five simple
908 algebras---four of which comprise Hermitian matrices. And the
909 last type of algebra can of course be thought of as `n`-by-`1`
910 column matrices (ambiguusly called column vectors) to avoid
911 special cases. As a result, matrices (and column vectors) are
912 a natural representation format for Euclidean Jordan algebra
915 But, when we construct an algebra from a basis of matrices,
916 those matrix representations are lost in favor of coordinate
917 vectors *with respect to* that basis. We could eventually
918 convert back if we tried hard enough, but having the original
919 representations handy is valuable enough that we simply store
920 them and return them from this method.
922 Why implement this for non-matrix algebras? Avoiding special
923 cases for the :class:`BilinearFormEJA` pays with simplicity in
924 its own right. But mainly, we would like to be able to assume
925 that elements of a :class:`CartesianProductEJA` can be displayed
926 nicely, without having to have special classes for direct sums
927 one of whose components was a matrix algebra.
931 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
932 ....: RealSymmetricEJA)
936 sage: J = RealSymmetricEJA(2)
938 Finite family {0: e0, 1: e1, 2: e2}
939 sage: J.matrix_basis()
941 [1 0] [ 0 0.7071067811865475?] [0 0]
942 [0 0], [0.7071067811865475? 0], [0 1]
947 sage: J = JordanSpinEJA(2)
949 Finite family {0: e0, 1: e1}
950 sage: J.matrix_basis()
956 return self
._matrix
_basis
959 def matrix_space(self
):
961 Return the matrix space in which this algebra's elements live, if
962 we think of them as matrices (including column vectors of the
965 "By default" this will be an `n`-by-`1` column-matrix space,
966 except when the algebra is trivial. There it's `n`-by-`n`
967 (where `n` is zero), to ensure that two elements of the matrix
968 space (empty matrices) can be multiplied. For algebras of
969 matrices, this returns the space in which their
970 real embeddings live.
974 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
976 ....: QuaternionHermitianEJA,
981 By default, the matrix representation is just a column-matrix
982 equivalent to the vector representation::
984 sage: J = JordanSpinEJA(3)
985 sage: J.matrix_space()
986 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
989 The matrix representation in the trivial algebra is
990 zero-by-zero instead of the usual `n`-by-one::
992 sage: J = TrivialEJA()
993 sage: J.matrix_space()
994 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
997 The matrix space for complex/quaternion Hermitian matrix EJA
998 is the space in which their real-embeddings live, not the
999 original complex/quaternion matrix space::
1001 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1002 sage: J.matrix_space()
1003 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1004 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1005 sage: J.matrix_space()
1006 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1009 if self
.is_trivial():
1010 return MatrixSpace(self
.base_ring(), 0)
1012 return self
.matrix_basis()[0].parent()
1018 Return the unit element of this algebra.
1022 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1027 We can compute unit element in the Hadamard EJA::
1029 sage: J = HadamardEJA(5)
1031 e0 + e1 + e2 + e3 + e4
1033 The unit element in the Hadamard EJA is inherited in the
1034 subalgebras generated by its elements::
1036 sage: J = HadamardEJA(5)
1038 e0 + e1 + e2 + e3 + e4
1039 sage: x = sum(J.gens())
1040 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1043 sage: A.one().superalgebra_element()
1044 e0 + e1 + e2 + e3 + e4
1048 The identity element acts like the identity, regardless of
1049 whether or not we orthonormalize::
1051 sage: set_random_seed()
1052 sage: J = random_eja()
1053 sage: x = J.random_element()
1054 sage: J.one()*x == x and x*J.one() == x
1056 sage: A = x.subalgebra_generated_by()
1057 sage: y = A.random_element()
1058 sage: A.one()*y == y and y*A.one() == y
1063 sage: set_random_seed()
1064 sage: J = random_eja(field=QQ, orthonormalize=False)
1065 sage: x = J.random_element()
1066 sage: J.one()*x == x and x*J.one() == x
1068 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1069 sage: y = A.random_element()
1070 sage: A.one()*y == y and y*A.one() == y
1073 The matrix of the unit element's operator is the identity,
1074 regardless of the base field and whether or not we
1077 sage: set_random_seed()
1078 sage: J = random_eja()
1079 sage: actual = J.one().operator().matrix()
1080 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1081 sage: actual == expected
1083 sage: x = J.random_element()
1084 sage: A = x.subalgebra_generated_by()
1085 sage: actual = A.one().operator().matrix()
1086 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1087 sage: actual == expected
1092 sage: set_random_seed()
1093 sage: J = random_eja(field=QQ, orthonormalize=False)
1094 sage: actual = J.one().operator().matrix()
1095 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1096 sage: actual == expected
1098 sage: x = J.random_element()
1099 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1100 sage: actual = A.one().operator().matrix()
1101 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1102 sage: actual == expected
1105 Ensure that the cached unit element (often precomputed by
1106 hand) agrees with the computed one::
1108 sage: set_random_seed()
1109 sage: J = random_eja()
1110 sage: cached = J.one()
1111 sage: J.one.clear_cache()
1112 sage: J.one() == cached
1117 sage: set_random_seed()
1118 sage: J = random_eja(field=QQ, orthonormalize=False)
1119 sage: cached = J.one()
1120 sage: J.one.clear_cache()
1121 sage: J.one() == cached
1125 # We can brute-force compute the matrices of the operators
1126 # that correspond to the basis elements of this algebra.
1127 # If some linear combination of those basis elements is the
1128 # algebra identity, then the same linear combination of
1129 # their matrices has to be the identity matrix.
1131 # Of course, matrices aren't vectors in sage, so we have to
1132 # appeal to the "long vectors" isometry.
1133 oper_vecs
= [ _mat2vec(g
.operator().matrix()) for g
in self
.gens() ]
1135 # Now we use basic linear algebra to find the coefficients,
1136 # of the matrices-as-vectors-linear-combination, which should
1137 # work for the original algebra basis too.
1138 A
= matrix(self
.base_ring(), oper_vecs
)
1140 # We used the isometry on the left-hand side already, but we
1141 # still need to do it for the right-hand side. Recall that we
1142 # wanted something that summed to the identity matrix.
1143 b
= _mat2vec( matrix
.identity(self
.base_ring(), self
.dimension()) )
1145 # Now if there's an identity element in the algebra, this
1146 # should work. We solve on the left to avoid having to
1147 # transpose the matrix "A".
1148 return self
.from_vector(A
.solve_left(b
))
1151 def peirce_decomposition(self
, c
):
1153 The Peirce decomposition of this algebra relative to the
1156 In the future, this can be extended to a complete system of
1157 orthogonal idempotents.
1161 - ``c`` -- an idempotent of this algebra.
1165 A triple (J0, J5, J1) containing two subalgebras and one subspace
1168 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1169 corresponding to the eigenvalue zero.
1171 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1172 corresponding to the eigenvalue one-half.
1174 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1175 corresponding to the eigenvalue one.
1177 These are the only possible eigenspaces for that operator, and this
1178 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1179 orthogonal, and are subalgebras of this algebra with the appropriate
1184 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1188 The canonical example comes from the symmetric matrices, which
1189 decompose into diagonal and off-diagonal parts::
1191 sage: J = RealSymmetricEJA(3)
1192 sage: C = matrix(QQ, [ [1,0,0],
1196 sage: J0,J5,J1 = J.peirce_decomposition(c)
1198 Euclidean Jordan algebra of dimension 1...
1200 Vector space of degree 6 and dimension 2...
1202 Euclidean Jordan algebra of dimension 3...
1203 sage: J0.one().to_matrix()
1207 sage: orig_df = AA.options.display_format
1208 sage: AA.options.display_format = 'radical'
1209 sage: J.from_vector(J5.basis()[0]).to_matrix()
1213 sage: J.from_vector(J5.basis()[1]).to_matrix()
1217 sage: AA.options.display_format = orig_df
1218 sage: J1.one().to_matrix()
1225 Every algebra decomposes trivially with respect to its identity
1228 sage: set_random_seed()
1229 sage: J = random_eja()
1230 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1231 sage: J0.dimension() == 0 and J5.dimension() == 0
1233 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1236 The decomposition is into eigenspaces, and its components are
1237 therefore necessarily orthogonal. Moreover, the identity
1238 elements in the two subalgebras are the projections onto their
1239 respective subspaces of the superalgebra's identity element::
1241 sage: set_random_seed()
1242 sage: J = random_eja()
1243 sage: x = J.random_element()
1244 sage: if not J.is_trivial():
1245 ....: while x.is_nilpotent():
1246 ....: x = J.random_element()
1247 sage: c = x.subalgebra_idempotent()
1248 sage: J0,J5,J1 = J.peirce_decomposition(c)
1250 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1251 ....: w = w.superalgebra_element()
1252 ....: y = J.from_vector(y)
1253 ....: z = z.superalgebra_element()
1254 ....: ipsum += w.inner_product(y).abs()
1255 ....: ipsum += w.inner_product(z).abs()
1256 ....: ipsum += y.inner_product(z).abs()
1259 sage: J1(c) == J1.one()
1261 sage: J0(J.one() - c) == J0.one()
1265 if not c
.is_idempotent():
1266 raise ValueError("element is not idempotent: %s" % c
)
1268 # Default these to what they should be if they turn out to be
1269 # trivial, because eigenspaces_left() won't return eigenvalues
1270 # corresponding to trivial spaces (e.g. it returns only the
1271 # eigenspace corresponding to lambda=1 if you take the
1272 # decomposition relative to the identity element).
1273 trivial
= self
.subalgebra(())
1274 J0
= trivial
# eigenvalue zero
1275 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1276 J1
= trivial
# eigenvalue one
1278 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1279 if eigval
== ~
(self
.base_ring()(2)):
1282 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1283 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1289 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1294 def random_element(self
, thorough
=False):
1296 Return a random element of this algebra.
1298 Our algebra superclass method only returns a linear
1299 combination of at most two basis elements. We instead
1300 want the vector space "random element" method that
1301 returns a more diverse selection.
1305 - ``thorough`` -- (boolean; default False) whether or not we
1306 should generate irrational coefficients for the random
1307 element when our base ring is irrational; this slows the
1308 algebra operations to a crawl, but any truly random method
1312 # For a general base ring... maybe we can trust this to do the
1313 # right thing? Unlikely, but.
1314 V
= self
.vector_space()
1315 v
= V
.random_element()
1317 if self
.base_ring() is AA
:
1318 # The "random element" method of the algebraic reals is
1319 # stupid at the moment, and only returns integers between
1320 # -2 and 2, inclusive:
1322 # https://trac.sagemath.org/ticket/30875
1324 # Instead, we implement our own "random vector" method,
1325 # and then coerce that into the algebra. We use the vector
1326 # space degree here instead of the dimension because a
1327 # subalgebra could (for example) be spanned by only two
1328 # vectors, each with five coordinates. We need to
1329 # generate all five coordinates.
1331 v
*= QQbar
.random_element().real()
1333 v
*= QQ
.random_element()
1335 return self
.from_vector(V
.coordinate_vector(v
))
1337 def random_elements(self
, count
, thorough
=False):
1339 Return ``count`` random elements as a tuple.
1343 - ``thorough`` -- (boolean; default False) whether or not we
1344 should generate irrational coefficients for the random
1345 elements when our base ring is irrational; this slows the
1346 algebra operations to a crawl, but any truly random method
1351 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1355 sage: J = JordanSpinEJA(3)
1356 sage: x,y,z = J.random_elements(3)
1357 sage: all( [ x in J, y in J, z in J ])
1359 sage: len( J.random_elements(10) ) == 10
1363 return tuple( self
.random_element(thorough
)
1364 for idx
in range(count
) )
1368 def _charpoly_coefficients(self
):
1370 The `r` polynomial coefficients of the "characteristic polynomial
1375 sage: from mjo.eja.eja_algebra import random_eja
1379 The theory shows that these are all homogeneous polynomials of
1382 sage: set_random_seed()
1383 sage: J = random_eja()
1384 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1388 n
= self
.dimension()
1389 R
= self
.coordinate_polynomial_ring()
1391 F
= R
.fraction_field()
1394 # From a result in my book, these are the entries of the
1395 # basis representation of L_x.
1396 return sum( vars[k
]*self
.gens()[k
].operator().matrix()[i
,j
]
1399 L_x
= matrix(F
, n
, n
, L_x_i_j
)
1402 if self
.rank
.is_in_cache():
1404 # There's no need to pad the system with redundant
1405 # columns if we *know* they'll be redundant.
1408 # Compute an extra power in case the rank is equal to
1409 # the dimension (otherwise, we would stop at x^(r-1)).
1410 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1411 for k
in range(n
+1) ]
1412 A
= matrix
.column(F
, x_powers
[:n
])
1413 AE
= A
.extended_echelon_form()
1420 # The theory says that only the first "r" coefficients are
1421 # nonzero, and they actually live in the original polynomial
1422 # ring and not the fraction field. We negate them because in
1423 # the actual characteristic polynomial, they get moved to the
1424 # other side where x^r lives. We don't bother to trim A_rref
1425 # down to a square matrix and solve the resulting system,
1426 # because the upper-left r-by-r portion of A_rref is
1427 # guaranteed to be the identity matrix, so e.g.
1429 # A_rref.solve_right(Y)
1431 # would just be returning Y.
1432 return (-E
*b
)[:r
].change_ring(R
)
1437 Return the rank of this EJA.
1439 This is a cached method because we know the rank a priori for
1440 all of the algebras we can construct. Thus we can avoid the
1441 expensive ``_charpoly_coefficients()`` call unless we truly
1442 need to compute the whole characteristic polynomial.
1446 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1447 ....: JordanSpinEJA,
1448 ....: RealSymmetricEJA,
1449 ....: ComplexHermitianEJA,
1450 ....: QuaternionHermitianEJA,
1455 The rank of the Jordan spin algebra is always two::
1457 sage: JordanSpinEJA(2).rank()
1459 sage: JordanSpinEJA(3).rank()
1461 sage: JordanSpinEJA(4).rank()
1464 The rank of the `n`-by-`n` Hermitian real, complex, or
1465 quaternion matrices is `n`::
1467 sage: RealSymmetricEJA(4).rank()
1469 sage: ComplexHermitianEJA(3).rank()
1471 sage: QuaternionHermitianEJA(2).rank()
1476 Ensure that every EJA that we know how to construct has a
1477 positive integer rank, unless the algebra is trivial in
1478 which case its rank will be zero::
1480 sage: set_random_seed()
1481 sage: J = random_eja()
1485 sage: r > 0 or (r == 0 and J.is_trivial())
1488 Ensure that computing the rank actually works, since the ranks
1489 of all simple algebras are known and will be cached by default::
1491 sage: set_random_seed() # long time
1492 sage: J = random_eja() # long time
1493 sage: cached = J.rank() # long time
1494 sage: J.rank.clear_cache() # long time
1495 sage: J.rank() == cached # long time
1499 return len(self
._charpoly
_coefficients
())
1502 def subalgebra(self
, basis
, **kwargs
):
1504 Create a subalgebra of this algebra from the given basis.
1506 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1507 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1510 def vector_space(self
):
1512 Return the vector space that underlies this algebra.
1516 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1520 sage: J = RealSymmetricEJA(2)
1521 sage: J.vector_space()
1522 Vector space of dimension 3 over...
1525 return self
.zero().to_vector().parent().ambient_vector_space()
1529 class RationalBasisEJA(FiniteDimensionalEJA
):
1531 New class for algebras whose supplied basis elements have all rational entries.
1535 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1539 The supplied basis is orthonormalized by default::
1541 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1542 sage: J = BilinearFormEJA(B)
1543 sage: J.matrix_basis()
1560 # Abuse the check_field parameter to check that the entries of
1561 # out basis (in ambient coordinates) are in the field QQ.
1562 if not all( all(b_i
in QQ
for b_i
in b
.list()) for b
in basis
):
1563 raise TypeError("basis not rational")
1565 self
._rational
_algebra
= None
1567 # There's no point in constructing the extra algebra if this
1568 # one is already rational.
1570 # Note: the same Jordan and inner-products work here,
1571 # because they are necessarily defined with respect to
1572 # ambient coordinates and not any particular basis.
1573 self
._rational
_algebra
= FiniteDimensionalEJA(
1578 orthonormalize
=False,
1582 super().__init
__(basis
,
1586 check_field
=check_field
,
1590 def _charpoly_coefficients(self
):
1594 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1595 ....: JordanSpinEJA)
1599 The base ring of the resulting polynomial coefficients is what
1600 it should be, and not the rationals (unless the algebra was
1601 already over the rationals)::
1603 sage: J = JordanSpinEJA(3)
1604 sage: J._charpoly_coefficients()
1605 (X1^2 - X2^2 - X3^2, -2*X1)
1606 sage: a0 = J._charpoly_coefficients()[0]
1608 Algebraic Real Field
1609 sage: a0.base_ring()
1610 Algebraic Real Field
1613 if self
._rational
_algebra
is None:
1614 # There's no need to construct *another* algebra over the
1615 # rationals if this one is already over the
1616 # rationals. Likewise, if we never orthonormalized our
1617 # basis, we might as well just use the given one.
1618 return super()._charpoly
_coefficients
()
1620 # Do the computation over the rationals. The answer will be
1621 # the same, because all we've done is a change of basis.
1622 # Then, change back from QQ to our real base ring
1623 a
= ( a_i
.change_ring(self
.base_ring())
1624 for a_i
in self
._rational
_algebra
._charpoly
_coefficients
() )
1626 if self
._deortho
_matrix
is None:
1627 # This can happen if our base ring was, say, AA and we
1628 # chose not to (or didn't need to) orthonormalize. It's
1629 # still faster to do the computations over QQ even if
1630 # the numbers in the boxes stay the same.
1633 # Otherwise, convert the coordinate variables back to the
1634 # deorthonormalized ones.
1635 R
= self
.coordinate_polynomial_ring()
1636 from sage
.modules
.free_module_element
import vector
1637 X
= vector(R
, R
.gens())
1638 BX
= self
._deortho
_matrix
*X
1640 subs_dict
= { X[i]: BX[i] for i in range(len(X)) }
1641 return tuple( a_i
.subs(subs_dict
) for a_i
in a
)
1643 class ConcreteEJA(RationalBasisEJA
):
1645 A class for the Euclidean Jordan algebras that we know by name.
1647 These are the Jordan algebras whose basis, multiplication table,
1648 rank, and so on are known a priori. More to the point, they are
1649 the Euclidean Jordan algebras for which we are able to conjure up
1650 a "random instance."
1654 sage: from mjo.eja.eja_algebra import ConcreteEJA
1658 Our basis is normalized with respect to the algebra's inner
1659 product, unless we specify otherwise::
1661 sage: set_random_seed()
1662 sage: J = ConcreteEJA.random_instance()
1663 sage: all( b.norm() == 1 for b in J.gens() )
1666 Since our basis is orthonormal with respect to the algebra's inner
1667 product, and since we know that this algebra is an EJA, any
1668 left-multiplication operator's matrix will be symmetric because
1669 natural->EJA basis representation is an isometry and within the
1670 EJA the operator is self-adjoint by the Jordan axiom::
1672 sage: set_random_seed()
1673 sage: J = ConcreteEJA.random_instance()
1674 sage: x = J.random_element()
1675 sage: x.operator().is_self_adjoint()
1680 def _max_random_instance_size():
1682 Return an integer "size" that is an upper bound on the size of
1683 this algebra when it is used in a random test
1684 case. Unfortunately, the term "size" is ambiguous -- when
1685 dealing with `R^n` under either the Hadamard or Jordan spin
1686 product, the "size" refers to the dimension `n`. When dealing
1687 with a matrix algebra (real symmetric or complex/quaternion
1688 Hermitian), it refers to the size of the matrix, which is far
1689 less than the dimension of the underlying vector space.
1691 This method must be implemented in each subclass.
1693 raise NotImplementedError
1696 def random_instance(cls
, *args
, **kwargs
):
1698 Return a random instance of this type of algebra.
1700 This method should be implemented in each subclass.
1702 from sage
.misc
.prandom
import choice
1703 eja_class
= choice(cls
.__subclasses
__())
1705 # These all bubble up to the RationalBasisEJA superclass
1706 # constructor, so any (kw)args valid there are also valid
1708 return eja_class
.random_instance(*args
, **kwargs
)
1713 def dimension_over_reals():
1715 The dimension of this matrix's base ring over the reals.
1717 The reals are dimension one over themselves, obviously; that's
1718 just `\mathbb{R}^{1}`. Likewise, the complex numbers `a + bi`
1719 have dimension two. Finally, the quaternions have dimension
1720 four over the reals.
1722 This is used to determine the size of the matrix returned from
1723 :meth:`real_embed`, among other things.
1725 raise NotImplementedError
1728 def real_embed(cls
,M
):
1730 Embed the matrix ``M`` into a space of real matrices.
1732 The matrix ``M`` can have entries in any field at the moment:
1733 the real numbers, complex numbers, or quaternions. And although
1734 they are not a field, we can probably support octonions at some
1735 point, too. This function returns a real matrix that "acts like"
1736 the original with respect to matrix multiplication; i.e.
1738 real_embed(M*N) = real_embed(M)*real_embed(N)
1741 if M
.ncols() != M
.nrows():
1742 raise ValueError("the matrix 'M' must be square")
1747 def real_unembed(cls
,M
):
1749 The inverse of :meth:`real_embed`.
1751 if M
.ncols() != M
.nrows():
1752 raise ValueError("the matrix 'M' must be square")
1753 if not ZZ(M
.nrows()).mod(cls
.dimension_over_reals()).is_zero():
1754 raise ValueError("the matrix 'M' must be a real embedding")
1758 def jordan_product(X
,Y
):
1759 return (X
*Y
+ Y
*X
)/2
1762 def trace_inner_product(cls
,X
,Y
):
1764 Compute the trace inner-product of two real-embeddings.
1768 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
1769 ....: ComplexHermitianEJA,
1770 ....: QuaternionHermitianEJA)
1774 This gives the same answer as it would if we computed the trace
1775 from the unembedded (original) matrices::
1777 sage: set_random_seed()
1778 sage: J = RealSymmetricEJA.random_instance()
1779 sage: x,y = J.random_elements(2)
1780 sage: Xe = x.to_matrix()
1781 sage: Ye = y.to_matrix()
1782 sage: X = J.real_unembed(Xe)
1783 sage: Y = J.real_unembed(Ye)
1784 sage: expected = (X*Y).trace()
1785 sage: actual = J.trace_inner_product(Xe,Ye)
1786 sage: actual == expected
1791 sage: set_random_seed()
1792 sage: J = ComplexHermitianEJA.random_instance()
1793 sage: x,y = J.random_elements(2)
1794 sage: Xe = x.to_matrix()
1795 sage: Ye = y.to_matrix()
1796 sage: X = J.real_unembed(Xe)
1797 sage: Y = J.real_unembed(Ye)
1798 sage: expected = (X*Y).trace().real()
1799 sage: actual = J.trace_inner_product(Xe,Ye)
1800 sage: actual == expected
1805 sage: set_random_seed()
1806 sage: J = QuaternionHermitianEJA.random_instance()
1807 sage: x,y = J.random_elements(2)
1808 sage: Xe = x.to_matrix()
1809 sage: Ye = y.to_matrix()
1810 sage: X = J.real_unembed(Xe)
1811 sage: Y = J.real_unembed(Ye)
1812 sage: expected = (X*Y).trace().coefficient_tuple()[0]
1813 sage: actual = J.trace_inner_product(Xe,Ye)
1814 sage: actual == expected
1818 Xu
= cls
.real_unembed(X
)
1819 Yu
= cls
.real_unembed(Y
)
1820 tr
= (Xu
*Yu
).trace()
1823 # Works in QQ, AA, RDF, et cetera.
1825 except AttributeError:
1826 # A quaternion doesn't have a real() method, but does
1827 # have coefficient_tuple() method that returns the
1828 # coefficients of 1, i, j, and k -- in that order.
1829 return tr
.coefficient_tuple()[0]
1832 class RealMatrixEJA(MatrixEJA
):
1834 def dimension_over_reals():
1838 class RealSymmetricEJA(ConcreteEJA
, RealMatrixEJA
):
1840 The rank-n simple EJA consisting of real symmetric n-by-n
1841 matrices, the usual symmetric Jordan product, and the trace inner
1842 product. It has dimension `(n^2 + n)/2` over the reals.
1846 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1850 sage: J = RealSymmetricEJA(2)
1851 sage: e0, e1, e2 = J.gens()
1859 In theory, our "field" can be any subfield of the reals::
1861 sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
1862 Euclidean Jordan algebra of dimension 3 over Real Double Field
1863 sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
1864 Euclidean Jordan algebra of dimension 3 over Real Field with
1865 53 bits of precision
1869 The dimension of this algebra is `(n^2 + n) / 2`::
1871 sage: set_random_seed()
1872 sage: n_max = RealSymmetricEJA._max_random_instance_size()
1873 sage: n = ZZ.random_element(1, n_max)
1874 sage: J = RealSymmetricEJA(n)
1875 sage: J.dimension() == (n^2 + n)/2
1878 The Jordan multiplication is what we think it is::
1880 sage: set_random_seed()
1881 sage: J = RealSymmetricEJA.random_instance()
1882 sage: x,y = J.random_elements(2)
1883 sage: actual = (x*y).to_matrix()
1884 sage: X = x.to_matrix()
1885 sage: Y = y.to_matrix()
1886 sage: expected = (X*Y + Y*X)/2
1887 sage: actual == expected
1889 sage: J(expected) == x*y
1892 We can change the generator prefix::
1894 sage: RealSymmetricEJA(3, prefix='q').gens()
1895 (q0, q1, q2, q3, q4, q5)
1897 We can construct the (trivial) algebra of rank zero::
1899 sage: RealSymmetricEJA(0)
1900 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
1904 def _denormalized_basis(cls
, n
):
1906 Return a basis for the space of real symmetric n-by-n matrices.
1910 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1914 sage: set_random_seed()
1915 sage: n = ZZ.random_element(1,5)
1916 sage: B = RealSymmetricEJA._denormalized_basis(n)
1917 sage: all( M.is_symmetric() for M in B)
1921 # The basis of symmetric matrices, as matrices, in their R^(n-by-n)
1925 for j
in range(i
+1):
1926 Eij
= matrix(ZZ
, n
, lambda k
,l
: k
==i
and l
==j
)
1930 Sij
= Eij
+ Eij
.transpose()
1936 def _max_random_instance_size():
1937 return 4 # Dimension 10
1940 def random_instance(cls
, **kwargs
):
1942 Return a random instance of this type of algebra.
1944 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
1945 return cls(n
, **kwargs
)
1947 def __init__(self
, n
, **kwargs
):
1948 # We know this is a valid EJA, but will double-check
1949 # if the user passes check_axioms=True.
1950 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
1952 super(RealSymmetricEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
1953 self
.jordan_product
,
1954 self
.trace_inner_product
,
1957 # TODO: this could be factored out somehow, but is left here
1958 # because the MatrixEJA is not presently a subclass of the
1959 # FDEJA class that defines rank() and one().
1960 self
.rank
.set_cache(n
)
1961 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
1962 self
.one
.set_cache(self(idV
))
1966 class ComplexMatrixEJA(MatrixEJA
):
1967 # A manual dictionary-cache for the complex_extension() method,
1968 # since apparently @classmethods can't also be @cached_methods.
1969 _complex_extension
= {}
1972 def complex_extension(cls
,field
):
1974 The complex field that we embed/unembed, as an extension
1975 of the given ``field``.
1977 if field
in cls
._complex
_extension
:
1978 return cls
._complex
_extension
[field
]
1980 # Sage doesn't know how to adjoin the complex "i" (the root of
1981 # x^2 + 1) to a field in a general way. Here, we just enumerate
1982 # all of the cases that I have cared to support so far.
1984 # Sage doesn't know how to embed AA into QQbar, i.e. how
1985 # to adjoin sqrt(-1) to AA.
1987 elif not field
.is_exact():
1989 F
= field
.complex_field()
1991 # Works for QQ and... maybe some other fields.
1992 R
= PolynomialRing(field
, 'z')
1994 F
= field
.extension(z
**2 + 1, 'I', embedding
=CLF(-1).sqrt())
1996 cls
._complex
_extension
[field
] = F
2000 def dimension_over_reals():
2004 def real_embed(cls
,M
):
2006 Embed the n-by-n complex matrix ``M`` into the space of real
2007 matrices of size 2n-by-2n via the map the sends each entry `z = a +
2008 bi` to the block matrix ``[[a,b],[-b,a]]``.
2012 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2016 sage: F = QuadraticField(-1, 'I')
2017 sage: x1 = F(4 - 2*i)
2018 sage: x2 = F(1 + 2*i)
2021 sage: M = matrix(F,2,[[x1,x2],[x3,x4]])
2022 sage: ComplexMatrixEJA.real_embed(M)
2031 Embedding is a homomorphism (isomorphism, in fact)::
2033 sage: set_random_seed()
2034 sage: n = ZZ.random_element(3)
2035 sage: F = QuadraticField(-1, 'I')
2036 sage: X = random_matrix(F, n)
2037 sage: Y = random_matrix(F, n)
2038 sage: Xe = ComplexMatrixEJA.real_embed(X)
2039 sage: Ye = ComplexMatrixEJA.real_embed(Y)
2040 sage: XYe = ComplexMatrixEJA.real_embed(X*Y)
2045 super(ComplexMatrixEJA
,cls
).real_embed(M
)
2048 # We don't need any adjoined elements...
2049 field
= M
.base_ring().base_ring()
2055 blocks
.append(matrix(field
, 2, [ [ a
, b
],
2058 return matrix
.block(field
, n
, blocks
)
2062 def real_unembed(cls
,M
):
2064 The inverse of _embed_complex_matrix().
2068 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2072 sage: A = matrix(QQ,[ [ 1, 2, 3, 4],
2073 ....: [-2, 1, -4, 3],
2074 ....: [ 9, 10, 11, 12],
2075 ....: [-10, 9, -12, 11] ])
2076 sage: ComplexMatrixEJA.real_unembed(A)
2078 [ 10*I + 9 12*I + 11]
2082 Unembedding is the inverse of embedding::
2084 sage: set_random_seed()
2085 sage: F = QuadraticField(-1, 'I')
2086 sage: M = random_matrix(F, 3)
2087 sage: Me = ComplexMatrixEJA.real_embed(M)
2088 sage: ComplexMatrixEJA.real_unembed(Me) == M
2092 super(ComplexMatrixEJA
,cls
).real_unembed(M
)
2094 d
= cls
.dimension_over_reals()
2095 F
= cls
.complex_extension(M
.base_ring())
2098 # Go top-left to bottom-right (reading order), converting every
2099 # 2-by-2 block we see to a single complex element.
2101 for k
in range(n
/d
):
2102 for j
in range(n
/d
):
2103 submat
= M
[d
*k
:d
*k
+d
,d
*j
:d
*j
+d
]
2104 if submat
[0,0] != submat
[1,1]:
2105 raise ValueError('bad on-diagonal submatrix')
2106 if submat
[0,1] != -submat
[1,0]:
2107 raise ValueError('bad off-diagonal submatrix')
2108 z
= submat
[0,0] + submat
[0,1]*i
2111 return matrix(F
, n
/d
, elements
)
2114 class ComplexHermitianEJA(ConcreteEJA
, ComplexMatrixEJA
):
2116 The rank-n simple EJA consisting of complex Hermitian n-by-n
2117 matrices over the real numbers, the usual symmetric Jordan product,
2118 and the real-part-of-trace inner product. It has dimension `n^2` over
2123 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2127 In theory, our "field" can be any subfield of the reals::
2129 sage: ComplexHermitianEJA(2, field=RDF, check_axioms=True)
2130 Euclidean Jordan algebra of dimension 4 over Real Double Field
2131 sage: ComplexHermitianEJA(2, field=RR, check_axioms=True)
2132 Euclidean Jordan algebra of dimension 4 over Real Field with
2133 53 bits of precision
2137 The dimension of this algebra is `n^2`::
2139 sage: set_random_seed()
2140 sage: n_max = ComplexHermitianEJA._max_random_instance_size()
2141 sage: n = ZZ.random_element(1, n_max)
2142 sage: J = ComplexHermitianEJA(n)
2143 sage: J.dimension() == n^2
2146 The Jordan multiplication is what we think it is::
2148 sage: set_random_seed()
2149 sage: J = ComplexHermitianEJA.random_instance()
2150 sage: x,y = J.random_elements(2)
2151 sage: actual = (x*y).to_matrix()
2152 sage: X = x.to_matrix()
2153 sage: Y = y.to_matrix()
2154 sage: expected = (X*Y + Y*X)/2
2155 sage: actual == expected
2157 sage: J(expected) == x*y
2160 We can change the generator prefix::
2162 sage: ComplexHermitianEJA(2, prefix='z').gens()
2165 We can construct the (trivial) algebra of rank zero::
2167 sage: ComplexHermitianEJA(0)
2168 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2173 def _denormalized_basis(cls
, n
):
2175 Returns a basis for the space of complex Hermitian n-by-n matrices.
2177 Why do we embed these? Basically, because all of numerical linear
2178 algebra assumes that you're working with vectors consisting of `n`
2179 entries from a field and scalars from the same field. There's no way
2180 to tell SageMath that (for example) the vectors contain complex
2181 numbers, while the scalar field is real.
2185 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2189 sage: set_random_seed()
2190 sage: n = ZZ.random_element(1,5)
2191 sage: B = ComplexHermitianEJA._denormalized_basis(n)
2192 sage: all( M.is_symmetric() for M in B)
2197 R
= PolynomialRing(field
, 'z')
2199 F
= field
.extension(z
**2 + 1, 'I')
2202 # This is like the symmetric case, but we need to be careful:
2204 # * We want conjugate-symmetry, not just symmetry.
2205 # * The diagonal will (as a result) be real.
2208 Eij
= matrix
.zero(F
,n
)
2210 for j
in range(i
+1):
2214 Sij
= cls
.real_embed(Eij
)
2217 # The second one has a minus because it's conjugated.
2218 Eij
[j
,i
] = 1 # Eij = Eij + Eij.transpose()
2219 Sij_real
= cls
.real_embed(Eij
)
2221 # Eij = I*Eij - I*Eij.transpose()
2224 Sij_imag
= cls
.real_embed(Eij
)
2230 # Since we embedded these, we can drop back to the "field" that we
2231 # started with instead of the complex extension "F".
2232 return tuple( s
.change_ring(field
) for s
in S
)
2235 def __init__(self
, n
, **kwargs
):
2236 # We know this is a valid EJA, but will double-check
2237 # if the user passes check_axioms=True.
2238 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2240 super(ComplexHermitianEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
2241 self
.jordan_product
,
2242 self
.trace_inner_product
,
2244 # TODO: this could be factored out somehow, but is left here
2245 # because the MatrixEJA is not presently a subclass of the
2246 # FDEJA class that defines rank() and one().
2247 self
.rank
.set_cache(n
)
2248 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2249 self
.one
.set_cache(self(idV
))
2252 def _max_random_instance_size():
2253 return 3 # Dimension 9
2256 def random_instance(cls
, **kwargs
):
2258 Return a random instance of this type of algebra.
2260 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2261 return cls(n
, **kwargs
)
2263 class QuaternionMatrixEJA(MatrixEJA
):
2265 # A manual dictionary-cache for the quaternion_extension() method,
2266 # since apparently @classmethods can't also be @cached_methods.
2267 _quaternion_extension
= {}
2270 def quaternion_extension(cls
,field
):
2272 The quaternion field that we embed/unembed, as an extension
2273 of the given ``field``.
2275 if field
in cls
._quaternion
_extension
:
2276 return cls
._quaternion
_extension
[field
]
2278 Q
= QuaternionAlgebra(field
,-1,-1)
2280 cls
._quaternion
_extension
[field
] = Q
2284 def dimension_over_reals():
2288 def real_embed(cls
,M
):
2290 Embed the n-by-n quaternion matrix ``M`` into the space of real
2291 matrices of size 4n-by-4n by first sending each quaternion entry `z
2292 = a + bi + cj + dk` to the block-complex matrix ``[[a + bi,
2293 c+di],[-c + di, a-bi]]`, and then embedding those into a real
2298 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2302 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2303 sage: i,j,k = Q.gens()
2304 sage: x = 1 + 2*i + 3*j + 4*k
2305 sage: M = matrix(Q, 1, [[x]])
2306 sage: QuaternionMatrixEJA.real_embed(M)
2312 Embedding is a homomorphism (isomorphism, in fact)::
2314 sage: set_random_seed()
2315 sage: n = ZZ.random_element(2)
2316 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2317 sage: X = random_matrix(Q, n)
2318 sage: Y = random_matrix(Q, n)
2319 sage: Xe = QuaternionMatrixEJA.real_embed(X)
2320 sage: Ye = QuaternionMatrixEJA.real_embed(Y)
2321 sage: XYe = QuaternionMatrixEJA.real_embed(X*Y)
2326 super(QuaternionMatrixEJA
,cls
).real_embed(M
)
2327 quaternions
= M
.base_ring()
2330 F
= QuadraticField(-1, 'I')
2335 t
= z
.coefficient_tuple()
2340 cplxM
= matrix(F
, 2, [[ a
+ b
*i
, c
+ d
*i
],
2341 [-c
+ d
*i
, a
- b
*i
]])
2342 realM
= ComplexMatrixEJA
.real_embed(cplxM
)
2343 blocks
.append(realM
)
2345 # We should have real entries by now, so use the realest field
2346 # we've got for the return value.
2347 return matrix
.block(quaternions
.base_ring(), n
, blocks
)
2352 def real_unembed(cls
,M
):
2354 The inverse of _embed_quaternion_matrix().
2358 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2362 sage: M = matrix(QQ, [[ 1, 2, 3, 4],
2363 ....: [-2, 1, -4, 3],
2364 ....: [-3, 4, 1, -2],
2365 ....: [-4, -3, 2, 1]])
2366 sage: QuaternionMatrixEJA.real_unembed(M)
2367 [1 + 2*i + 3*j + 4*k]
2371 Unembedding is the inverse of embedding::
2373 sage: set_random_seed()
2374 sage: Q = QuaternionAlgebra(QQ, -1, -1)
2375 sage: M = random_matrix(Q, 3)
2376 sage: Me = QuaternionMatrixEJA.real_embed(M)
2377 sage: QuaternionMatrixEJA.real_unembed(Me) == M
2381 super(QuaternionMatrixEJA
,cls
).real_unembed(M
)
2383 d
= cls
.dimension_over_reals()
2385 # Use the base ring of the matrix to ensure that its entries can be
2386 # multiplied by elements of the quaternion algebra.
2387 Q
= cls
.quaternion_extension(M
.base_ring())
2390 # Go top-left to bottom-right (reading order), converting every
2391 # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
2394 for l
in range(n
/d
):
2395 for m
in range(n
/d
):
2396 submat
= ComplexMatrixEJA
.real_unembed(
2397 M
[d
*l
:d
*l
+d
,d
*m
:d
*m
+d
] )
2398 if submat
[0,0] != submat
[1,1].conjugate():
2399 raise ValueError('bad on-diagonal submatrix')
2400 if submat
[0,1] != -submat
[1,0].conjugate():
2401 raise ValueError('bad off-diagonal submatrix')
2402 z
= submat
[0,0].real()
2403 z
+= submat
[0,0].imag()*i
2404 z
+= submat
[0,1].real()*j
2405 z
+= submat
[0,1].imag()*k
2408 return matrix(Q
, n
/d
, elements
)
2411 class QuaternionHermitianEJA(ConcreteEJA
, QuaternionMatrixEJA
):
2413 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2414 matrices, the usual symmetric Jordan product, and the
2415 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2420 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2424 In theory, our "field" can be any subfield of the reals::
2426 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2427 Euclidean Jordan algebra of dimension 6 over Real Double Field
2428 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2429 Euclidean Jordan algebra of dimension 6 over Real Field with
2430 53 bits of precision
2434 The dimension of this algebra is `2*n^2 - n`::
2436 sage: set_random_seed()
2437 sage: n_max = QuaternionHermitianEJA._max_random_instance_size()
2438 sage: n = ZZ.random_element(1, n_max)
2439 sage: J = QuaternionHermitianEJA(n)
2440 sage: J.dimension() == 2*(n^2) - n
2443 The Jordan multiplication is what we think it is::
2445 sage: set_random_seed()
2446 sage: J = QuaternionHermitianEJA.random_instance()
2447 sage: x,y = J.random_elements(2)
2448 sage: actual = (x*y).to_matrix()
2449 sage: X = x.to_matrix()
2450 sage: Y = y.to_matrix()
2451 sage: expected = (X*Y + Y*X)/2
2452 sage: actual == expected
2454 sage: J(expected) == x*y
2457 We can change the generator prefix::
2459 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2460 (a0, a1, a2, a3, a4, a5)
2462 We can construct the (trivial) algebra of rank zero::
2464 sage: QuaternionHermitianEJA(0)
2465 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2469 def _denormalized_basis(cls
, n
):
2471 Returns a basis for the space of quaternion Hermitian n-by-n matrices.
2473 Why do we embed these? Basically, because all of numerical
2474 linear algebra assumes that you're working with vectors consisting
2475 of `n` entries from a field and scalars from the same field. There's
2476 no way to tell SageMath that (for example) the vectors contain
2477 complex numbers, while the scalar field is real.
2481 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2485 sage: set_random_seed()
2486 sage: n = ZZ.random_element(1,5)
2487 sage: B = QuaternionHermitianEJA._denormalized_basis(n)
2488 sage: all( M.is_symmetric() for M in B )
2493 Q
= QuaternionAlgebra(QQ
,-1,-1)
2496 # This is like the symmetric case, but we need to be careful:
2498 # * We want conjugate-symmetry, not just symmetry.
2499 # * The diagonal will (as a result) be real.
2502 Eij
= matrix
.zero(Q
,n
)
2504 for j
in range(i
+1):
2508 Sij
= cls
.real_embed(Eij
)
2511 # The second, third, and fourth ones have a minus
2512 # because they're conjugated.
2513 # Eij = Eij + Eij.transpose()
2515 Sij_real
= cls
.real_embed(Eij
)
2517 # Eij = I*(Eij - Eij.transpose())
2520 Sij_I
= cls
.real_embed(Eij
)
2522 # Eij = J*(Eij - Eij.transpose())
2525 Sij_J
= cls
.real_embed(Eij
)
2527 # Eij = K*(Eij - Eij.transpose())
2530 Sij_K
= cls
.real_embed(Eij
)
2536 # Since we embedded these, we can drop back to the "field" that we
2537 # started with instead of the quaternion algebra "Q".
2538 return tuple( s
.change_ring(field
) for s
in S
)
2541 def __init__(self
, n
, **kwargs
):
2542 # We know this is a valid EJA, but will double-check
2543 # if the user passes check_axioms=True.
2544 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2546 super(QuaternionHermitianEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
2547 self
.jordan_product
,
2548 self
.trace_inner_product
,
2550 # TODO: this could be factored out somehow, but is left here
2551 # because the MatrixEJA is not presently a subclass of the
2552 # FDEJA class that defines rank() and one().
2553 self
.rank
.set_cache(n
)
2554 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2555 self
.one
.set_cache(self(idV
))
2559 def _max_random_instance_size():
2561 The maximum rank of a random QuaternionHermitianEJA.
2563 return 2 # Dimension 6
2566 def random_instance(cls
, **kwargs
):
2568 Return a random instance of this type of algebra.
2570 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2571 return cls(n
, **kwargs
)
2574 class HadamardEJA(ConcreteEJA
):
2576 Return the Euclidean Jordan Algebra corresponding to the set
2577 `R^n` under the Hadamard product.
2579 Note: this is nothing more than the Cartesian product of ``n``
2580 copies of the spin algebra. Once Cartesian product algebras
2581 are implemented, this can go.
2585 sage: from mjo.eja.eja_algebra import HadamardEJA
2589 This multiplication table can be verified by hand::
2591 sage: J = HadamardEJA(3)
2592 sage: e0,e1,e2 = J.gens()
2608 We can change the generator prefix::
2610 sage: HadamardEJA(3, prefix='r').gens()
2614 def __init__(self
, n
, **kwargs
):
2616 jordan_product
= lambda x
,y
: x
2617 inner_product
= lambda x
,y
: x
2619 def jordan_product(x
,y
):
2621 return P( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2623 def inner_product(x
,y
):
2626 # New defaults for keyword arguments. Don't orthonormalize
2627 # because our basis is already orthonormal with respect to our
2628 # inner-product. Don't check the axioms, because we know this
2629 # is a valid EJA... but do double-check if the user passes
2630 # check_axioms=True. Note: we DON'T override the "check_field"
2631 # default here, because the user can pass in a field!
2632 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2633 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2635 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2636 super().__init
__(column_basis
,
2641 self
.rank
.set_cache(n
)
2644 self
.one
.set_cache( self
.zero() )
2646 self
.one
.set_cache( sum(self
.gens()) )
2649 def _max_random_instance_size():
2651 The maximum dimension of a random HadamardEJA.
2656 def random_instance(cls
, **kwargs
):
2658 Return a random instance of this type of algebra.
2660 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2661 return cls(n
, **kwargs
)
2664 class BilinearFormEJA(ConcreteEJA
):
2666 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2667 with the half-trace inner product and jordan product ``x*y =
2668 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2669 a symmetric positive-definite "bilinear form" matrix. Its
2670 dimension is the size of `B`, and it has rank two in dimensions
2671 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2672 the identity matrix of order ``n``.
2674 We insist that the one-by-one upper-left identity block of `B` be
2675 passed in as well so that we can be passed a matrix of size zero
2676 to construct a trivial algebra.
2680 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2681 ....: JordanSpinEJA)
2685 When no bilinear form is specified, the identity matrix is used,
2686 and the resulting algebra is the Jordan spin algebra::
2688 sage: B = matrix.identity(AA,3)
2689 sage: J0 = BilinearFormEJA(B)
2690 sage: J1 = JordanSpinEJA(3)
2691 sage: J0.multiplication_table() == J0.multiplication_table()
2694 An error is raised if the matrix `B` does not correspond to a
2695 positive-definite bilinear form::
2697 sage: B = matrix.random(QQ,2,3)
2698 sage: J = BilinearFormEJA(B)
2699 Traceback (most recent call last):
2701 ValueError: bilinear form is not positive-definite
2702 sage: B = matrix.zero(QQ,3)
2703 sage: J = BilinearFormEJA(B)
2704 Traceback (most recent call last):
2706 ValueError: bilinear form is not positive-definite
2710 We can create a zero-dimensional algebra::
2712 sage: B = matrix.identity(AA,0)
2713 sage: J = BilinearFormEJA(B)
2717 We can check the multiplication condition given in the Jordan, von
2718 Neumann, and Wigner paper (and also discussed on my "On the
2719 symmetry..." paper). Note that this relies heavily on the standard
2720 choice of basis, as does anything utilizing the bilinear form
2721 matrix. We opt not to orthonormalize the basis, because if we
2722 did, we would have to normalize the `s_{i}` in a similar manner::
2724 sage: set_random_seed()
2725 sage: n = ZZ.random_element(5)
2726 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2727 sage: B11 = matrix.identity(QQ,1)
2728 sage: B22 = M.transpose()*M
2729 sage: B = block_matrix(2,2,[ [B11,0 ],
2731 sage: J = BilinearFormEJA(B, orthonormalize=False)
2732 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2733 sage: V = J.vector_space()
2734 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2735 ....: for ei in eis ]
2736 sage: actual = [ sis[i]*sis[j]
2737 ....: for i in range(n-1)
2738 ....: for j in range(n-1) ]
2739 sage: expected = [ J.one() if i == j else J.zero()
2740 ....: for i in range(n-1)
2741 ....: for j in range(n-1) ]
2742 sage: actual == expected
2746 def __init__(self
, B
, **kwargs
):
2747 # The matrix "B" is supplied by the user in most cases,
2748 # so it makes sense to check whether or not its positive-
2749 # definite unless we are specifically asked not to...
2750 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2751 if not B
.is_positive_definite():
2752 raise ValueError("bilinear form is not positive-definite")
2754 # However, all of the other data for this EJA is computed
2755 # by us in manner that guarantees the axioms are
2756 # satisfied. So, again, unless we are specifically asked to
2757 # verify things, we'll skip the rest of the checks.
2758 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2760 def inner_product(x
,y
):
2761 return (y
.T
*B
*x
)[0,0]
2763 def jordan_product(x
,y
):
2769 z0
= inner_product(y
,x
)
2770 zbar
= y0
*xbar
+ x0
*ybar
2771 return P([z0
] + zbar
.list())
2774 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2775 super(BilinearFormEJA
, self
).__init
__(column_basis
,
2780 # The rank of this algebra is two, unless we're in a
2781 # one-dimensional ambient space (because the rank is bounded
2782 # by the ambient dimension).
2783 self
.rank
.set_cache(min(n
,2))
2786 self
.one
.set_cache( self
.zero() )
2788 self
.one
.set_cache( self
.monomial(0) )
2791 def _max_random_instance_size():
2793 The maximum dimension of a random BilinearFormEJA.
2798 def random_instance(cls
, **kwargs
):
2800 Return a random instance of this algebra.
2802 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2804 B
= matrix
.identity(ZZ
, n
)
2805 return cls(B
, **kwargs
)
2807 B11
= matrix
.identity(ZZ
, 1)
2808 M
= matrix
.random(ZZ
, n
-1)
2809 I
= matrix
.identity(ZZ
, n
-1)
2811 while alpha
.is_zero():
2812 alpha
= ZZ
.random_element().abs()
2813 B22
= M
.transpose()*M
+ alpha
*I
2815 from sage
.matrix
.special
import block_matrix
2816 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2819 return cls(B
, **kwargs
)
2822 class JordanSpinEJA(BilinearFormEJA
):
2824 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2825 with the usual inner product and jordan product ``x*y =
2826 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2831 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2835 This multiplication table can be verified by hand::
2837 sage: J = JordanSpinEJA(4)
2838 sage: e0,e1,e2,e3 = J.gens()
2854 We can change the generator prefix::
2856 sage: JordanSpinEJA(2, prefix='B').gens()
2861 Ensure that we have the usual inner product on `R^n`::
2863 sage: set_random_seed()
2864 sage: J = JordanSpinEJA.random_instance()
2865 sage: x,y = J.random_elements(2)
2866 sage: actual = x.inner_product(y)
2867 sage: expected = x.to_vector().inner_product(y.to_vector())
2868 sage: actual == expected
2872 def __init__(self
, n
, **kwargs
):
2873 # This is a special case of the BilinearFormEJA with the
2874 # identity matrix as its bilinear form.
2875 B
= matrix
.identity(ZZ
, n
)
2877 # Don't orthonormalize because our basis is already
2878 # orthonormal with respect to our inner-product.
2879 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2881 # But also don't pass check_field=False here, because the user
2882 # can pass in a field!
2883 super(JordanSpinEJA
, self
).__init
__(B
, **kwargs
)
2886 def _max_random_instance_size():
2888 The maximum dimension of a random JordanSpinEJA.
2893 def random_instance(cls
, **kwargs
):
2895 Return a random instance of this type of algebra.
2897 Needed here to override the implementation for ``BilinearFormEJA``.
2899 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2900 return cls(n
, **kwargs
)
2903 class TrivialEJA(ConcreteEJA
):
2905 The trivial Euclidean Jordan algebra consisting of only a zero element.
2909 sage: from mjo.eja.eja_algebra import TrivialEJA
2913 sage: J = TrivialEJA()
2920 sage: 7*J.one()*12*J.one()
2922 sage: J.one().inner_product(J.one())
2924 sage: J.one().norm()
2926 sage: J.one().subalgebra_generated_by()
2927 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2932 def __init__(self
, **kwargs
):
2933 jordan_product
= lambda x
,y
: x
2934 inner_product
= lambda x
,y
: 0
2937 # New defaults for keyword arguments
2938 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2939 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2941 super(TrivialEJA
, self
).__init
__(basis
,
2945 # The rank is zero using my definition, namely the dimension of the
2946 # largest subalgebra generated by any element.
2947 self
.rank
.set_cache(0)
2948 self
.one
.set_cache( self
.zero() )
2951 def random_instance(cls
, **kwargs
):
2952 # We don't take a "size" argument so the superclass method is
2953 # inappropriate for us.
2954 return cls(**kwargs
)
2957 class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct
,
2958 FiniteDimensionalEJA
):
2960 The external (orthogonal) direct sum of two or more Euclidean
2961 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
2962 orthogonal direct sum of simple Euclidean Jordan algebras which is
2963 then isometric to a Cartesian product, so no generality is lost by
2964 providing only this construction.
2968 sage: from mjo.eja.eja_algebra import (random_eja,
2969 ....: CartesianProductEJA,
2971 ....: JordanSpinEJA,
2972 ....: RealSymmetricEJA)
2976 The Jordan product is inherited from our factors and implemented by
2977 our CombinatorialFreeModule Cartesian product superclass::
2979 sage: set_random_seed()
2980 sage: J1 = HadamardEJA(2)
2981 sage: J2 = RealSymmetricEJA(2)
2982 sage: J = cartesian_product([J1,J2])
2983 sage: x,y = J.random_elements(2)
2987 The ability to retrieve the original factors is implemented by our
2988 CombinatorialFreeModule Cartesian product superclass::
2990 sage: J1 = HadamardEJA(2, field=QQ)
2991 sage: J2 = JordanSpinEJA(3, field=QQ)
2992 sage: J = cartesian_product([J1,J2])
2993 sage: J.cartesian_factors()
2994 (Euclidean Jordan algebra of dimension 2 over Rational Field,
2995 Euclidean Jordan algebra of dimension 3 over Rational Field)
2997 You can provide more than two factors::
2999 sage: J1 = HadamardEJA(2)
3000 sage: J2 = JordanSpinEJA(3)
3001 sage: J3 = RealSymmetricEJA(3)
3002 sage: cartesian_product([J1,J2,J3])
3003 Euclidean Jordan algebra of dimension 2 over Algebraic Real
3004 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
3005 Real Field (+) Euclidean Jordan algebra of dimension 6 over
3006 Algebraic Real Field
3008 Rank is additive on a Cartesian product::
3010 sage: J1 = HadamardEJA(1)
3011 sage: J2 = RealSymmetricEJA(2)
3012 sage: J = cartesian_product([J1,J2])
3013 sage: J1.rank.clear_cache()
3014 sage: J2.rank.clear_cache()
3015 sage: J.rank.clear_cache()
3018 sage: J.rank() == J1.rank() + J2.rank()
3021 The same rank computation works over the rationals, with whatever
3024 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
3025 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
3026 sage: J = cartesian_product([J1,J2])
3027 sage: J1.rank.clear_cache()
3028 sage: J2.rank.clear_cache()
3029 sage: J.rank.clear_cache()
3032 sage: J.rank() == J1.rank() + J2.rank()
3035 The product algebra will be associative if and only if all of its
3036 components are associative::
3038 sage: J1 = HadamardEJA(2)
3039 sage: J1.is_associative()
3041 sage: J2 = HadamardEJA(3)
3042 sage: J2.is_associative()
3044 sage: J3 = RealSymmetricEJA(3)
3045 sage: J3.is_associative()
3047 sage: CP1 = cartesian_product([J1,J2])
3048 sage: CP1.is_associative()
3050 sage: CP2 = cartesian_product([J1,J3])
3051 sage: CP2.is_associative()
3056 All factors must share the same base field::
3058 sage: J1 = HadamardEJA(2, field=QQ)
3059 sage: J2 = RealSymmetricEJA(2)
3060 sage: CartesianProductEJA((J1,J2))
3061 Traceback (most recent call last):
3063 ValueError: all factors must share the same base field
3065 The cached unit element is the same one that would be computed::
3067 sage: set_random_seed() # long time
3068 sage: J1 = random_eja() # long time
3069 sage: J2 = random_eja() # long time
3070 sage: J = cartesian_product([J1,J2]) # long time
3071 sage: actual = J.one() # long time
3072 sage: J.one.clear_cache() # long time
3073 sage: expected = J.one() # long time
3074 sage: actual == expected # long time
3078 Element
= FiniteDimensionalEJAElement
3081 def __init__(self
, algebras
, **kwargs
):
3082 CombinatorialFreeModule_CartesianProduct
.__init
__(self
,
3085 field
= algebras
[0].base_ring()
3086 if not all( J
.base_ring() == field
for J
in algebras
):
3087 raise ValueError("all factors must share the same base field")
3089 associative
= all( m
.is_associative() for m
in algebras
)
3091 # The definition of matrix_space() and self.basis() relies
3092 # only on the stuff in the CFM_CartesianProduct class, which
3093 # we've already initialized.
3094 Js
= self
.cartesian_factors()
3096 MS
= self
.matrix_space()
3098 MS(tuple( self
.cartesian_projection(i
)(b
).to_matrix()
3099 for i
in range(m
) ))
3100 for b
in self
.basis()
3103 # Define jordan/inner products that operate on that matrix_basis.
3104 def jordan_product(x
,y
):
3106 (Js
[i
](x
[i
])*Js
[i
](y
[i
])).to_matrix() for i
in range(m
)
3109 def inner_product(x
, y
):
3111 Js
[i
](x
[i
]).inner_product(Js
[i
](y
[i
])) for i
in range(m
)
3114 # There's no need to check the field since it already came
3115 # from an EJA. Likewise the axioms are guaranteed to be
3116 # satisfied, unless the guy writing this class sucks.
3118 # If you want the basis to be orthonormalized, orthonormalize
3120 FiniteDimensionalEJA
.__init
__(self
,
3125 orthonormalize
=False,
3126 associative
=associative
,
3127 cartesian_product
=True,
3131 ones
= tuple(J
.one() for J
in algebras
)
3132 self
.one
.set_cache(self
._cartesian
_product
_of
_elements
(ones
))
3133 self
.rank
.set_cache(sum(J
.rank() for J
in algebras
))
3135 def matrix_space(self
):
3137 Return the space that our matrix basis lives in as a Cartesian
3142 sage: from mjo.eja.eja_algebra import (HadamardEJA,
3143 ....: RealSymmetricEJA)
3147 sage: J1 = HadamardEJA(1)
3148 sage: J2 = RealSymmetricEJA(2)
3149 sage: J = cartesian_product([J1,J2])
3150 sage: J.matrix_space()
3151 The Cartesian product of (Full MatrixSpace of 1 by 1 dense
3152 matrices over Algebraic Real Field, Full MatrixSpace of 2
3153 by 2 dense matrices over Algebraic Real Field)
3156 from sage
.categories
.cartesian_product
import cartesian_product
3157 return cartesian_product( [J
.matrix_space()
3158 for J
in self
.cartesian_factors()] )
3161 def cartesian_projection(self
, i
):
3165 sage: from mjo.eja.eja_algebra import (random_eja,
3166 ....: JordanSpinEJA,
3168 ....: RealSymmetricEJA,
3169 ....: ComplexHermitianEJA)
3173 The projection morphisms are Euclidean Jordan algebra
3176 sage: J1 = HadamardEJA(2)
3177 sage: J2 = RealSymmetricEJA(2)
3178 sage: J = cartesian_product([J1,J2])
3179 sage: J.cartesian_projection(0)
3180 Linear operator between finite-dimensional Euclidean Jordan
3181 algebras represented by the matrix:
3184 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3185 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3186 Algebraic Real Field
3187 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3189 sage: J.cartesian_projection(1)
3190 Linear operator between finite-dimensional Euclidean Jordan
3191 algebras represented by the matrix:
3195 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3196 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3197 Algebraic Real Field
3198 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3201 The projections work the way you'd expect on the vector
3202 representation of an element::
3204 sage: J1 = JordanSpinEJA(2)
3205 sage: J2 = ComplexHermitianEJA(2)
3206 sage: J = cartesian_product([J1,J2])
3207 sage: pi_left = J.cartesian_projection(0)
3208 sage: pi_right = J.cartesian_projection(1)
3209 sage: pi_left(J.one()).to_vector()
3211 sage: pi_right(J.one()).to_vector()
3213 sage: J.one().to_vector()
3218 The answer never changes::
3220 sage: set_random_seed()
3221 sage: J1 = random_eja()
3222 sage: J2 = random_eja()
3223 sage: J = cartesian_product([J1,J2])
3224 sage: P0 = J.cartesian_projection(0)
3225 sage: P1 = J.cartesian_projection(0)
3230 Ji
= self
.cartesian_factors()[i
]
3231 # Requires the fix on Trac 31421/31422 to work!
3232 Pi
= super().cartesian_projection(i
)
3233 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3236 def cartesian_embedding(self
, i
):
3240 sage: from mjo.eja.eja_algebra import (random_eja,
3241 ....: JordanSpinEJA,
3243 ....: RealSymmetricEJA)
3247 The embedding morphisms are Euclidean Jordan algebra
3250 sage: J1 = HadamardEJA(2)
3251 sage: J2 = RealSymmetricEJA(2)
3252 sage: J = cartesian_product([J1,J2])
3253 sage: J.cartesian_embedding(0)
3254 Linear operator between finite-dimensional Euclidean Jordan
3255 algebras represented by the matrix:
3261 Domain: Euclidean Jordan algebra of dimension 2 over
3262 Algebraic Real Field
3263 Codomain: Euclidean Jordan algebra of dimension 2 over
3264 Algebraic Real Field (+) Euclidean Jordan algebra of
3265 dimension 3 over Algebraic Real Field
3266 sage: J.cartesian_embedding(1)
3267 Linear operator between finite-dimensional Euclidean Jordan
3268 algebras represented by the matrix:
3274 Domain: Euclidean Jordan algebra of dimension 3 over
3275 Algebraic Real Field
3276 Codomain: Euclidean Jordan algebra of dimension 2 over
3277 Algebraic Real Field (+) Euclidean Jordan algebra of
3278 dimension 3 over Algebraic Real Field
3280 The embeddings work the way you'd expect on the vector
3281 representation of an element::
3283 sage: J1 = JordanSpinEJA(3)
3284 sage: J2 = RealSymmetricEJA(2)
3285 sage: J = cartesian_product([J1,J2])
3286 sage: iota_left = J.cartesian_embedding(0)
3287 sage: iota_right = J.cartesian_embedding(1)
3288 sage: iota_left(J1.zero()) == J.zero()
3290 sage: iota_right(J2.zero()) == J.zero()
3292 sage: J1.one().to_vector()
3294 sage: iota_left(J1.one()).to_vector()
3296 sage: J2.one().to_vector()
3298 sage: iota_right(J2.one()).to_vector()
3300 sage: J.one().to_vector()
3305 The answer never changes::
3307 sage: set_random_seed()
3308 sage: J1 = random_eja()
3309 sage: J2 = random_eja()
3310 sage: J = cartesian_product([J1,J2])
3311 sage: E0 = J.cartesian_embedding(0)
3312 sage: E1 = J.cartesian_embedding(0)
3316 Composing a projection with the corresponding inclusion should
3317 produce the identity map, and mismatching them should produce
3320 sage: set_random_seed()
3321 sage: J1 = random_eja()
3322 sage: J2 = random_eja()
3323 sage: J = cartesian_product([J1,J2])
3324 sage: iota_left = J.cartesian_embedding(0)
3325 sage: iota_right = J.cartesian_embedding(1)
3326 sage: pi_left = J.cartesian_projection(0)
3327 sage: pi_right = J.cartesian_projection(1)
3328 sage: pi_left*iota_left == J1.one().operator()
3330 sage: pi_right*iota_right == J2.one().operator()
3332 sage: (pi_left*iota_right).is_zero()
3334 sage: (pi_right*iota_left).is_zero()
3338 Ji
= self
.cartesian_factors()[i
]
3339 # Requires the fix on Trac 31421/31422 to work!
3340 Ei
= super().cartesian_embedding(i
)
3341 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3345 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3347 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3350 A separate class for products of algebras for which we know a
3355 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
3356 ....: RealSymmetricEJA)
3360 This gives us fast characteristic polynomial computations in
3361 product algebras, too::
3364 sage: J1 = JordanSpinEJA(2)
3365 sage: J2 = RealSymmetricEJA(3)
3366 sage: J = cartesian_product([J1,J2])
3367 sage: J.characteristic_polynomial_of().degree()
3373 def __init__(self
, algebras
, **kwargs
):
3374 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3376 self
._rational
_algebra
= None
3377 if self
.vector_space().base_field() is not QQ
:
3378 self
._rational
_algebra
= cartesian_product([
3379 r
._rational
_algebra
for r
in algebras
3383 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3385 random_eja
= ConcreteEJA
.random_instance