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 (JordanSpinEJA,
573 ....: RealSymmetricEJA)
577 The identity in `S^n` is converted to the identity in the EJA::
579 sage: J = RealSymmetricEJA(3)
580 sage: I = matrix.identity(QQ,3)
581 sage: J(I) == J.one()
584 This skew-symmetric matrix can't be represented in the EJA::
586 sage: J = RealSymmetricEJA(3)
587 sage: A = matrix(QQ,3, lambda i,j: i-j)
589 Traceback (most recent call last):
591 ValueError: not an element of this algebra
593 Tuples work as well, provided that the matrix basis for the
594 algebra consists of them::
596 sage: J1 = HadamardEJA(3)
597 sage: J2 = RealSymmetricEJA(2)
598 sage: J = cartesian_product([J1,J2])
599 sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
604 Ensure that we can convert any element back and forth
605 faithfully between its matrix and algebra representations::
607 sage: set_random_seed()
608 sage: J = random_eja()
609 sage: x = J.random_element()
610 sage: J(x.to_matrix()) == x
613 We cannot coerce elements between algebras just because their
614 matrix representations are compatible::
616 sage: J1 = HadamardEJA(3)
617 sage: J2 = JordanSpinEJA(3)
619 Traceback (most recent call last):
621 ValueError: not an element of this algebra
623 Traceback (most recent call last):
625 ValueError: not an element of this algebra
627 msg
= "not an element of this algebra"
628 if elt
in self
.base_ring():
629 # Ensure that no base ring -> algebra coercion is performed
630 # by this method. There's some stupidity in sage that would
631 # otherwise propagate to this method; for example, sage thinks
632 # that the integer 3 belongs to the space of 2-by-2 matrices.
633 raise ValueError(msg
)
636 # Try to convert a vector into a column-matrix...
638 except (AttributeError, TypeError):
639 # and ignore failure, because we weren't really expecting
640 # a vector as an argument anyway.
643 if elt
not in self
.matrix_space():
644 raise ValueError(msg
)
646 # Thanks for nothing! Matrix spaces aren't vector spaces in
647 # Sage, so we have to figure out its matrix-basis coordinates
648 # ourselves. We use the basis space's ring instead of the
649 # element's ring because the basis space might be an algebraic
650 # closure whereas the base ring of the 3-by-3 identity matrix
651 # could be QQ instead of QQbar.
653 # And, we also have to handle Cartesian product bases (when
654 # the matrix basis consists of tuples) here. The "good news"
655 # is that we're already converting everything to long vectors,
656 # and that strategy works for tuples as well.
658 # We pass check=False because the matrix basis is "guaranteed"
659 # to be linearly independent... right? Ha ha.
661 V
= VectorSpace(self
.base_ring(), len(elt
))
662 W
= V
.span_of_basis( (V(_all2list(s
)) for s
in self
.matrix_basis()),
666 coords
= W
.coordinate_vector(V(elt
))
667 except ArithmeticError: # vector is not in free module
668 raise ValueError(msg
)
670 return self
.from_vector(coords
)
674 Return a string representation of ``self``.
678 sage: from mjo.eja.eja_algebra import JordanSpinEJA
682 Ensure that it says what we think it says::
684 sage: JordanSpinEJA(2, field=AA)
685 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
686 sage: JordanSpinEJA(3, field=RDF)
687 Euclidean Jordan algebra of dimension 3 over Real Double Field
690 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
691 return fmt
.format(self
.dimension(), self
.base_ring())
695 def characteristic_polynomial_of(self
):
697 Return the algebra's "characteristic polynomial of" function,
698 which is itself a multivariate polynomial that, when evaluated
699 at the coordinates of some algebra element, returns that
700 element's characteristic polynomial.
702 The resulting polynomial has `n+1` variables, where `n` is the
703 dimension of this algebra. The first `n` variables correspond to
704 the coordinates of an algebra element: when evaluated at the
705 coordinates of an algebra element with respect to a certain
706 basis, the result is a univariate polynomial (in the one
707 remaining variable ``t``), namely the characteristic polynomial
712 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
716 The characteristic polynomial in the spin algebra is given in
717 Alizadeh, Example 11.11::
719 sage: J = JordanSpinEJA(3)
720 sage: p = J.characteristic_polynomial_of(); p
721 X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
722 sage: xvec = J.one().to_vector()
726 By definition, the characteristic polynomial is a monic
727 degree-zero polynomial in a rank-zero algebra. Note that
728 Cayley-Hamilton is indeed satisfied since the polynomial
729 ``1`` evaluates to the identity element of the algebra on
732 sage: J = TrivialEJA()
733 sage: J.characteristic_polynomial_of()
740 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
741 a
= self
._charpoly
_coefficients
()
743 # We go to a bit of trouble here to reorder the
744 # indeterminates, so that it's easier to evaluate the
745 # characteristic polynomial at x's coordinates and get back
746 # something in terms of t, which is what we want.
747 S
= PolynomialRing(self
.base_ring(),'t')
751 S
= PolynomialRing(S
, R
.variable_names())
754 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
756 def coordinate_polynomial_ring(self
):
758 The multivariate polynomial ring in which this algebra's
759 :meth:`characteristic_polynomial_of` lives.
763 sage: from mjo.eja.eja_algebra import (HadamardEJA,
764 ....: RealSymmetricEJA)
768 sage: J = HadamardEJA(2)
769 sage: J.coordinate_polynomial_ring()
770 Multivariate Polynomial Ring in X1, X2...
771 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
772 sage: J.coordinate_polynomial_ring()
773 Multivariate Polynomial Ring in X1, X2, X3, X4, X5, X6...
776 var_names
= tuple( "X%d" % z
for z
in range(1, self
.dimension()+1) )
777 return PolynomialRing(self
.base_ring(), var_names
)
779 def inner_product(self
, x
, y
):
781 The inner product associated with this Euclidean Jordan algebra.
783 Defaults to the trace inner product, but can be overridden by
784 subclasses if they are sure that the necessary properties are
789 sage: from mjo.eja.eja_algebra import (random_eja,
791 ....: BilinearFormEJA)
795 Our inner product is "associative," which means the following for
796 a symmetric bilinear form::
798 sage: set_random_seed()
799 sage: J = random_eja()
800 sage: x,y,z = J.random_elements(3)
801 sage: (x*y).inner_product(z) == y.inner_product(x*z)
806 Ensure that this is the usual inner product for the algebras
809 sage: set_random_seed()
810 sage: J = HadamardEJA.random_instance()
811 sage: x,y = J.random_elements(2)
812 sage: actual = x.inner_product(y)
813 sage: expected = x.to_vector().inner_product(y.to_vector())
814 sage: actual == expected
817 Ensure that this is one-half of the trace inner-product in a
818 BilinearFormEJA that isn't just the reals (when ``n`` isn't
819 one). This is in Faraut and Koranyi, and also my "On the
822 sage: set_random_seed()
823 sage: J = BilinearFormEJA.random_instance()
824 sage: n = J.dimension()
825 sage: x = J.random_element()
826 sage: y = J.random_element()
827 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
830 B
= self
._inner
_product
_matrix
831 return (B
*x
.to_vector()).inner_product(y
.to_vector())
834 def is_trivial(self
):
836 Return whether or not this algebra is trivial.
838 A trivial algebra contains only the zero element.
842 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
847 sage: J = ComplexHermitianEJA(3)
853 sage: J = TrivialEJA()
858 return self
.dimension() == 0
861 def multiplication_table(self
):
863 Return a visual representation of this algebra's multiplication
864 table (on basis elements).
868 sage: from mjo.eja.eja_algebra import JordanSpinEJA
872 sage: J = JordanSpinEJA(4)
873 sage: J.multiplication_table()
874 +----++----+----+----+----+
875 | * || e0 | e1 | e2 | e3 |
876 +====++====+====+====+====+
877 | e0 || e0 | e1 | e2 | e3 |
878 +----++----+----+----+----+
879 | e1 || e1 | e0 | 0 | 0 |
880 +----++----+----+----+----+
881 | e2 || e2 | 0 | e0 | 0 |
882 +----++----+----+----+----+
883 | e3 || e3 | 0 | 0 | e0 |
884 +----++----+----+----+----+
888 # Prepend the header row.
889 M
= [["*"] + list(self
.gens())]
891 # And to each subsequent row, prepend an entry that belongs to
892 # the left-side "header column."
893 M
+= [ [self
.gens()[i
]] + [ self
.product_on_basis(i
,j
)
897 return table(M
, header_row
=True, header_column
=True, frame
=True)
900 def matrix_basis(self
):
902 Return an (often more natural) representation of this algebras
903 basis as an ordered tuple of matrices.
905 Every finite-dimensional Euclidean Jordan Algebra is a, up to
906 Jordan isomorphism, a direct sum of five simple
907 algebras---four of which comprise Hermitian matrices. And the
908 last type of algebra can of course be thought of as `n`-by-`1`
909 column matrices (ambiguusly called column vectors) to avoid
910 special cases. As a result, matrices (and column vectors) are
911 a natural representation format for Euclidean Jordan algebra
914 But, when we construct an algebra from a basis of matrices,
915 those matrix representations are lost in favor of coordinate
916 vectors *with respect to* that basis. We could eventually
917 convert back if we tried hard enough, but having the original
918 representations handy is valuable enough that we simply store
919 them and return them from this method.
921 Why implement this for non-matrix algebras? Avoiding special
922 cases for the :class:`BilinearFormEJA` pays with simplicity in
923 its own right. But mainly, we would like to be able to assume
924 that elements of a :class:`CartesianProductEJA` can be displayed
925 nicely, without having to have special classes for direct sums
926 one of whose components was a matrix algebra.
930 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
931 ....: RealSymmetricEJA)
935 sage: J = RealSymmetricEJA(2)
937 Finite family {0: e0, 1: e1, 2: e2}
938 sage: J.matrix_basis()
940 [1 0] [ 0 0.7071067811865475?] [0 0]
941 [0 0], [0.7071067811865475? 0], [0 1]
946 sage: J = JordanSpinEJA(2)
948 Finite family {0: e0, 1: e1}
949 sage: J.matrix_basis()
955 return self
._matrix
_basis
958 def matrix_space(self
):
960 Return the matrix space in which this algebra's elements live, if
961 we think of them as matrices (including column vectors of the
964 "By default" this will be an `n`-by-`1` column-matrix space,
965 except when the algebra is trivial. There it's `n`-by-`n`
966 (where `n` is zero), to ensure that two elements of the matrix
967 space (empty matrices) can be multiplied. For algebras of
968 matrices, this returns the space in which their
969 real embeddings live.
973 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
975 ....: QuaternionHermitianEJA,
980 By default, the matrix representation is just a column-matrix
981 equivalent to the vector representation::
983 sage: J = JordanSpinEJA(3)
984 sage: J.matrix_space()
985 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
988 The matrix representation in the trivial algebra is
989 zero-by-zero instead of the usual `n`-by-one::
991 sage: J = TrivialEJA()
992 sage: J.matrix_space()
993 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
996 The matrix space for complex/quaternion Hermitian matrix EJA
997 is the space in which their real-embeddings live, not the
998 original complex/quaternion matrix space::
1000 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1001 sage: J.matrix_space()
1002 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1003 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1004 sage: J.matrix_space()
1005 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1008 if self
.is_trivial():
1009 return MatrixSpace(self
.base_ring(), 0)
1011 return self
.matrix_basis()[0].parent()
1017 Return the unit element of this algebra.
1021 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1026 We can compute unit element in the Hadamard EJA::
1028 sage: J = HadamardEJA(5)
1030 e0 + e1 + e2 + e3 + e4
1032 The unit element in the Hadamard EJA is inherited in the
1033 subalgebras generated by its elements::
1035 sage: J = HadamardEJA(5)
1037 e0 + e1 + e2 + e3 + e4
1038 sage: x = sum(J.gens())
1039 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1042 sage: A.one().superalgebra_element()
1043 e0 + e1 + e2 + e3 + e4
1047 The identity element acts like the identity, regardless of
1048 whether or not we orthonormalize::
1050 sage: set_random_seed()
1051 sage: J = random_eja()
1052 sage: x = J.random_element()
1053 sage: J.one()*x == x and x*J.one() == x
1055 sage: A = x.subalgebra_generated_by()
1056 sage: y = A.random_element()
1057 sage: A.one()*y == y and y*A.one() == y
1062 sage: set_random_seed()
1063 sage: J = random_eja(field=QQ, orthonormalize=False)
1064 sage: x = J.random_element()
1065 sage: J.one()*x == x and x*J.one() == x
1067 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1068 sage: y = A.random_element()
1069 sage: A.one()*y == y and y*A.one() == y
1072 The matrix of the unit element's operator is the identity,
1073 regardless of the base field and whether or not we
1076 sage: set_random_seed()
1077 sage: J = random_eja()
1078 sage: actual = J.one().operator().matrix()
1079 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1080 sage: actual == expected
1082 sage: x = J.random_element()
1083 sage: A = x.subalgebra_generated_by()
1084 sage: actual = A.one().operator().matrix()
1085 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1086 sage: actual == expected
1091 sage: set_random_seed()
1092 sage: J = random_eja(field=QQ, orthonormalize=False)
1093 sage: actual = J.one().operator().matrix()
1094 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1095 sage: actual == expected
1097 sage: x = J.random_element()
1098 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1099 sage: actual = A.one().operator().matrix()
1100 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1101 sage: actual == expected
1104 Ensure that the cached unit element (often precomputed by
1105 hand) agrees with the computed one::
1107 sage: set_random_seed()
1108 sage: J = random_eja()
1109 sage: cached = J.one()
1110 sage: J.one.clear_cache()
1111 sage: J.one() == cached
1116 sage: set_random_seed()
1117 sage: J = random_eja(field=QQ, orthonormalize=False)
1118 sage: cached = J.one()
1119 sage: J.one.clear_cache()
1120 sage: J.one() == cached
1124 # We can brute-force compute the matrices of the operators
1125 # that correspond to the basis elements of this algebra.
1126 # If some linear combination of those basis elements is the
1127 # algebra identity, then the same linear combination of
1128 # their matrices has to be the identity matrix.
1130 # Of course, matrices aren't vectors in sage, so we have to
1131 # appeal to the "long vectors" isometry.
1132 oper_vecs
= [ _mat2vec(g
.operator().matrix()) for g
in self
.gens() ]
1134 # Now we use basic linear algebra to find the coefficients,
1135 # of the matrices-as-vectors-linear-combination, which should
1136 # work for the original algebra basis too.
1137 A
= matrix(self
.base_ring(), oper_vecs
)
1139 # We used the isometry on the left-hand side already, but we
1140 # still need to do it for the right-hand side. Recall that we
1141 # wanted something that summed to the identity matrix.
1142 b
= _mat2vec( matrix
.identity(self
.base_ring(), self
.dimension()) )
1144 # Now if there's an identity element in the algebra, this
1145 # should work. We solve on the left to avoid having to
1146 # transpose the matrix "A".
1147 return self
.from_vector(A
.solve_left(b
))
1150 def peirce_decomposition(self
, c
):
1152 The Peirce decomposition of this algebra relative to the
1155 In the future, this can be extended to a complete system of
1156 orthogonal idempotents.
1160 - ``c`` -- an idempotent of this algebra.
1164 A triple (J0, J5, J1) containing two subalgebras and one subspace
1167 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1168 corresponding to the eigenvalue zero.
1170 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1171 corresponding to the eigenvalue one-half.
1173 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1174 corresponding to the eigenvalue one.
1176 These are the only possible eigenspaces for that operator, and this
1177 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1178 orthogonal, and are subalgebras of this algebra with the appropriate
1183 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1187 The canonical example comes from the symmetric matrices, which
1188 decompose into diagonal and off-diagonal parts::
1190 sage: J = RealSymmetricEJA(3)
1191 sage: C = matrix(QQ, [ [1,0,0],
1195 sage: J0,J5,J1 = J.peirce_decomposition(c)
1197 Euclidean Jordan algebra of dimension 1...
1199 Vector space of degree 6 and dimension 2...
1201 Euclidean Jordan algebra of dimension 3...
1202 sage: J0.one().to_matrix()
1206 sage: orig_df = AA.options.display_format
1207 sage: AA.options.display_format = 'radical'
1208 sage: J.from_vector(J5.basis()[0]).to_matrix()
1212 sage: J.from_vector(J5.basis()[1]).to_matrix()
1216 sage: AA.options.display_format = orig_df
1217 sage: J1.one().to_matrix()
1224 Every algebra decomposes trivially with respect to its identity
1227 sage: set_random_seed()
1228 sage: J = random_eja()
1229 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1230 sage: J0.dimension() == 0 and J5.dimension() == 0
1232 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1235 The decomposition is into eigenspaces, and its components are
1236 therefore necessarily orthogonal. Moreover, the identity
1237 elements in the two subalgebras are the projections onto their
1238 respective subspaces of the superalgebra's identity element::
1240 sage: set_random_seed()
1241 sage: J = random_eja()
1242 sage: x = J.random_element()
1243 sage: if not J.is_trivial():
1244 ....: while x.is_nilpotent():
1245 ....: x = J.random_element()
1246 sage: c = x.subalgebra_idempotent()
1247 sage: J0,J5,J1 = J.peirce_decomposition(c)
1249 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1250 ....: w = w.superalgebra_element()
1251 ....: y = J.from_vector(y)
1252 ....: z = z.superalgebra_element()
1253 ....: ipsum += w.inner_product(y).abs()
1254 ....: ipsum += w.inner_product(z).abs()
1255 ....: ipsum += y.inner_product(z).abs()
1258 sage: J1(c) == J1.one()
1260 sage: J0(J.one() - c) == J0.one()
1264 if not c
.is_idempotent():
1265 raise ValueError("element is not idempotent: %s" % c
)
1267 # Default these to what they should be if they turn out to be
1268 # trivial, because eigenspaces_left() won't return eigenvalues
1269 # corresponding to trivial spaces (e.g. it returns only the
1270 # eigenspace corresponding to lambda=1 if you take the
1271 # decomposition relative to the identity element).
1272 trivial
= self
.subalgebra(())
1273 J0
= trivial
# eigenvalue zero
1274 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1275 J1
= trivial
# eigenvalue one
1277 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1278 if eigval
== ~
(self
.base_ring()(2)):
1281 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1282 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1288 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1293 def random_element(self
, thorough
=False):
1295 Return a random element of this algebra.
1297 Our algebra superclass method only returns a linear
1298 combination of at most two basis elements. We instead
1299 want the vector space "random element" method that
1300 returns a more diverse selection.
1304 - ``thorough`` -- (boolean; default False) whether or not we
1305 should generate irrational coefficients for the random
1306 element when our base ring is irrational; this slows the
1307 algebra operations to a crawl, but any truly random method
1311 # For a general base ring... maybe we can trust this to do the
1312 # right thing? Unlikely, but.
1313 V
= self
.vector_space()
1314 v
= V
.random_element()
1316 if self
.base_ring() is AA
:
1317 # The "random element" method of the algebraic reals is
1318 # stupid at the moment, and only returns integers between
1319 # -2 and 2, inclusive:
1321 # https://trac.sagemath.org/ticket/30875
1323 # Instead, we implement our own "random vector" method,
1324 # and then coerce that into the algebra. We use the vector
1325 # space degree here instead of the dimension because a
1326 # subalgebra could (for example) be spanned by only two
1327 # vectors, each with five coordinates. We need to
1328 # generate all five coordinates.
1330 v
*= QQbar
.random_element().real()
1332 v
*= QQ
.random_element()
1334 return self
.from_vector(V
.coordinate_vector(v
))
1336 def random_elements(self
, count
, thorough
=False):
1338 Return ``count`` random elements as a tuple.
1342 - ``thorough`` -- (boolean; default False) whether or not we
1343 should generate irrational coefficients for the random
1344 elements when our base ring is irrational; this slows the
1345 algebra operations to a crawl, but any truly random method
1350 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1354 sage: J = JordanSpinEJA(3)
1355 sage: x,y,z = J.random_elements(3)
1356 sage: all( [ x in J, y in J, z in J ])
1358 sage: len( J.random_elements(10) ) == 10
1362 return tuple( self
.random_element(thorough
)
1363 for idx
in range(count
) )
1367 def _charpoly_coefficients(self
):
1369 The `r` polynomial coefficients of the "characteristic polynomial
1374 sage: from mjo.eja.eja_algebra import random_eja
1378 The theory shows that these are all homogeneous polynomials of
1381 sage: set_random_seed()
1382 sage: J = random_eja()
1383 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1387 n
= self
.dimension()
1388 R
= self
.coordinate_polynomial_ring()
1390 F
= R
.fraction_field()
1393 # From a result in my book, these are the entries of the
1394 # basis representation of L_x.
1395 return sum( vars[k
]*self
.gens()[k
].operator().matrix()[i
,j
]
1398 L_x
= matrix(F
, n
, n
, L_x_i_j
)
1401 if self
.rank
.is_in_cache():
1403 # There's no need to pad the system with redundant
1404 # columns if we *know* they'll be redundant.
1407 # Compute an extra power in case the rank is equal to
1408 # the dimension (otherwise, we would stop at x^(r-1)).
1409 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1410 for k
in range(n
+1) ]
1411 A
= matrix
.column(F
, x_powers
[:n
])
1412 AE
= A
.extended_echelon_form()
1419 # The theory says that only the first "r" coefficients are
1420 # nonzero, and they actually live in the original polynomial
1421 # ring and not the fraction field. We negate them because in
1422 # the actual characteristic polynomial, they get moved to the
1423 # other side where x^r lives. We don't bother to trim A_rref
1424 # down to a square matrix and solve the resulting system,
1425 # because the upper-left r-by-r portion of A_rref is
1426 # guaranteed to be the identity matrix, so e.g.
1428 # A_rref.solve_right(Y)
1430 # would just be returning Y.
1431 return (-E
*b
)[:r
].change_ring(R
)
1436 Return the rank of this EJA.
1438 This is a cached method because we know the rank a priori for
1439 all of the algebras we can construct. Thus we can avoid the
1440 expensive ``_charpoly_coefficients()`` call unless we truly
1441 need to compute the whole characteristic polynomial.
1445 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1446 ....: JordanSpinEJA,
1447 ....: RealSymmetricEJA,
1448 ....: ComplexHermitianEJA,
1449 ....: QuaternionHermitianEJA,
1454 The rank of the Jordan spin algebra is always two::
1456 sage: JordanSpinEJA(2).rank()
1458 sage: JordanSpinEJA(3).rank()
1460 sage: JordanSpinEJA(4).rank()
1463 The rank of the `n`-by-`n` Hermitian real, complex, or
1464 quaternion matrices is `n`::
1466 sage: RealSymmetricEJA(4).rank()
1468 sage: ComplexHermitianEJA(3).rank()
1470 sage: QuaternionHermitianEJA(2).rank()
1475 Ensure that every EJA that we know how to construct has a
1476 positive integer rank, unless the algebra is trivial in
1477 which case its rank will be zero::
1479 sage: set_random_seed()
1480 sage: J = random_eja()
1484 sage: r > 0 or (r == 0 and J.is_trivial())
1487 Ensure that computing the rank actually works, since the ranks
1488 of all simple algebras are known and will be cached by default::
1490 sage: set_random_seed() # long time
1491 sage: J = random_eja() # long time
1492 sage: cached = J.rank() # long time
1493 sage: J.rank.clear_cache() # long time
1494 sage: J.rank() == cached # long time
1498 return len(self
._charpoly
_coefficients
())
1501 def subalgebra(self
, basis
, **kwargs
):
1503 Create a subalgebra of this algebra from the given basis.
1505 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1506 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1509 def vector_space(self
):
1511 Return the vector space that underlies this algebra.
1515 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1519 sage: J = RealSymmetricEJA(2)
1520 sage: J.vector_space()
1521 Vector space of dimension 3 over...
1524 return self
.zero().to_vector().parent().ambient_vector_space()
1528 class RationalBasisEJA(FiniteDimensionalEJA
):
1530 New class for algebras whose supplied basis elements have all rational entries.
1534 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1538 The supplied basis is orthonormalized by default::
1540 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1541 sage: J = BilinearFormEJA(B)
1542 sage: J.matrix_basis()
1559 # Abuse the check_field parameter to check that the entries of
1560 # out basis (in ambient coordinates) are in the field QQ.
1561 if not all( all(b_i
in QQ
for b_i
in b
.list()) for b
in basis
):
1562 raise TypeError("basis not rational")
1564 self
._rational
_algebra
= None
1566 # There's no point in constructing the extra algebra if this
1567 # one is already rational.
1569 # Note: the same Jordan and inner-products work here,
1570 # because they are necessarily defined with respect to
1571 # ambient coordinates and not any particular basis.
1572 self
._rational
_algebra
= FiniteDimensionalEJA(
1577 orthonormalize
=False,
1581 super().__init
__(basis
,
1585 check_field
=check_field
,
1589 def _charpoly_coefficients(self
):
1593 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1594 ....: JordanSpinEJA)
1598 The base ring of the resulting polynomial coefficients is what
1599 it should be, and not the rationals (unless the algebra was
1600 already over the rationals)::
1602 sage: J = JordanSpinEJA(3)
1603 sage: J._charpoly_coefficients()
1604 (X1^2 - X2^2 - X3^2, -2*X1)
1605 sage: a0 = J._charpoly_coefficients()[0]
1607 Algebraic Real Field
1608 sage: a0.base_ring()
1609 Algebraic Real Field
1612 if self
._rational
_algebra
is None:
1613 # There's no need to construct *another* algebra over the
1614 # rationals if this one is already over the
1615 # rationals. Likewise, if we never orthonormalized our
1616 # basis, we might as well just use the given one.
1617 return super()._charpoly
_coefficients
()
1619 # Do the computation over the rationals. The answer will be
1620 # the same, because all we've done is a change of basis.
1621 # Then, change back from QQ to our real base ring
1622 a
= ( a_i
.change_ring(self
.base_ring())
1623 for a_i
in self
._rational
_algebra
._charpoly
_coefficients
() )
1625 if self
._deortho
_matrix
is None:
1626 # This can happen if our base ring was, say, AA and we
1627 # chose not to (or didn't need to) orthonormalize. It's
1628 # still faster to do the computations over QQ even if
1629 # the numbers in the boxes stay the same.
1632 # Otherwise, convert the coordinate variables back to the
1633 # deorthonormalized ones.
1634 R
= self
.coordinate_polynomial_ring()
1635 from sage
.modules
.free_module_element
import vector
1636 X
= vector(R
, R
.gens())
1637 BX
= self
._deortho
_matrix
*X
1639 subs_dict
= { X[i]: BX[i] for i in range(len(X)) }
1640 return tuple( a_i
.subs(subs_dict
) for a_i
in a
)
1642 class ConcreteEJA(RationalBasisEJA
):
1644 A class for the Euclidean Jordan algebras that we know by name.
1646 These are the Jordan algebras whose basis, multiplication table,
1647 rank, and so on are known a priori. More to the point, they are
1648 the Euclidean Jordan algebras for which we are able to conjure up
1649 a "random instance."
1653 sage: from mjo.eja.eja_algebra import ConcreteEJA
1657 Our basis is normalized with respect to the algebra's inner
1658 product, unless we specify otherwise::
1660 sage: set_random_seed()
1661 sage: J = ConcreteEJA.random_instance()
1662 sage: all( b.norm() == 1 for b in J.gens() )
1665 Since our basis is orthonormal with respect to the algebra's inner
1666 product, and since we know that this algebra is an EJA, any
1667 left-multiplication operator's matrix will be symmetric because
1668 natural->EJA basis representation is an isometry and within the
1669 EJA the operator is self-adjoint by the Jordan axiom::
1671 sage: set_random_seed()
1672 sage: J = ConcreteEJA.random_instance()
1673 sage: x = J.random_element()
1674 sage: x.operator().is_self_adjoint()
1679 def _max_random_instance_size():
1681 Return an integer "size" that is an upper bound on the size of
1682 this algebra when it is used in a random test
1683 case. Unfortunately, the term "size" is ambiguous -- when
1684 dealing with `R^n` under either the Hadamard or Jordan spin
1685 product, the "size" refers to the dimension `n`. When dealing
1686 with a matrix algebra (real symmetric or complex/quaternion
1687 Hermitian), it refers to the size of the matrix, which is far
1688 less than the dimension of the underlying vector space.
1690 This method must be implemented in each subclass.
1692 raise NotImplementedError
1695 def random_instance(cls
, *args
, **kwargs
):
1697 Return a random instance of this type of algebra.
1699 This method should be implemented in each subclass.
1701 from sage
.misc
.prandom
import choice
1702 eja_class
= choice(cls
.__subclasses
__())
1704 # These all bubble up to the RationalBasisEJA superclass
1705 # constructor, so any (kw)args valid there are also valid
1707 return eja_class
.random_instance(*args
, **kwargs
)
1712 def dimension_over_reals():
1714 The dimension of this matrix's base ring over the reals.
1716 The reals are dimension one over themselves, obviously; that's
1717 just `\mathbb{R}^{1}`. Likewise, the complex numbers `a + bi`
1718 have dimension two. Finally, the quaternions have dimension
1719 four over the reals.
1721 This is used to determine the size of the matrix returned from
1722 :meth:`real_embed`, among other things.
1724 raise NotImplementedError
1727 def real_embed(cls
,M
):
1729 Embed the matrix ``M`` into a space of real matrices.
1731 The matrix ``M`` can have entries in any field at the moment:
1732 the real numbers, complex numbers, or quaternions. And although
1733 they are not a field, we can probably support octonions at some
1734 point, too. This function returns a real matrix that "acts like"
1735 the original with respect to matrix multiplication; i.e.
1737 real_embed(M*N) = real_embed(M)*real_embed(N)
1740 if M
.ncols() != M
.nrows():
1741 raise ValueError("the matrix 'M' must be square")
1746 def real_unembed(cls
,M
):
1748 The inverse of :meth:`real_embed`.
1750 if M
.ncols() != M
.nrows():
1751 raise ValueError("the matrix 'M' must be square")
1752 if not ZZ(M
.nrows()).mod(cls
.dimension_over_reals()).is_zero():
1753 raise ValueError("the matrix 'M' must be a real embedding")
1757 def jordan_product(X
,Y
):
1758 return (X
*Y
+ Y
*X
)/2
1761 def trace_inner_product(cls
,X
,Y
):
1763 Compute the trace inner-product of two real-embeddings.
1767 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
1768 ....: ComplexHermitianEJA,
1769 ....: QuaternionHermitianEJA)
1773 This gives the same answer as it would if we computed the trace
1774 from the unembedded (original) matrices::
1776 sage: set_random_seed()
1777 sage: J = RealSymmetricEJA.random_instance()
1778 sage: x,y = J.random_elements(2)
1779 sage: Xe = x.to_matrix()
1780 sage: Ye = y.to_matrix()
1781 sage: X = J.real_unembed(Xe)
1782 sage: Y = J.real_unembed(Ye)
1783 sage: expected = (X*Y).trace()
1784 sage: actual = J.trace_inner_product(Xe,Ye)
1785 sage: actual == expected
1790 sage: set_random_seed()
1791 sage: J = ComplexHermitianEJA.random_instance()
1792 sage: x,y = J.random_elements(2)
1793 sage: Xe = x.to_matrix()
1794 sage: Ye = y.to_matrix()
1795 sage: X = J.real_unembed(Xe)
1796 sage: Y = J.real_unembed(Ye)
1797 sage: expected = (X*Y).trace().real()
1798 sage: actual = J.trace_inner_product(Xe,Ye)
1799 sage: actual == expected
1804 sage: set_random_seed()
1805 sage: J = QuaternionHermitianEJA.random_instance()
1806 sage: x,y = J.random_elements(2)
1807 sage: Xe = x.to_matrix()
1808 sage: Ye = y.to_matrix()
1809 sage: X = J.real_unembed(Xe)
1810 sage: Y = J.real_unembed(Ye)
1811 sage: expected = (X*Y).trace().coefficient_tuple()[0]
1812 sage: actual = J.trace_inner_product(Xe,Ye)
1813 sage: actual == expected
1817 Xu
= cls
.real_unembed(X
)
1818 Yu
= cls
.real_unembed(Y
)
1819 tr
= (Xu
*Yu
).trace()
1822 # Works in QQ, AA, RDF, et cetera.
1824 except AttributeError:
1825 # A quaternion doesn't have a real() method, but does
1826 # have coefficient_tuple() method that returns the
1827 # coefficients of 1, i, j, and k -- in that order.
1828 return tr
.coefficient_tuple()[0]
1831 class RealMatrixEJA(MatrixEJA
):
1833 def dimension_over_reals():
1837 class RealSymmetricEJA(ConcreteEJA
, RealMatrixEJA
):
1839 The rank-n simple EJA consisting of real symmetric n-by-n
1840 matrices, the usual symmetric Jordan product, and the trace inner
1841 product. It has dimension `(n^2 + n)/2` over the reals.
1845 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1849 sage: J = RealSymmetricEJA(2)
1850 sage: e0, e1, e2 = J.gens()
1858 In theory, our "field" can be any subfield of the reals::
1860 sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
1861 Euclidean Jordan algebra of dimension 3 over Real Double Field
1862 sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
1863 Euclidean Jordan algebra of dimension 3 over Real Field with
1864 53 bits of precision
1868 The dimension of this algebra is `(n^2 + n) / 2`::
1870 sage: set_random_seed()
1871 sage: n_max = RealSymmetricEJA._max_random_instance_size()
1872 sage: n = ZZ.random_element(1, n_max)
1873 sage: J = RealSymmetricEJA(n)
1874 sage: J.dimension() == (n^2 + n)/2
1877 The Jordan multiplication is what we think it is::
1879 sage: set_random_seed()
1880 sage: J = RealSymmetricEJA.random_instance()
1881 sage: x,y = J.random_elements(2)
1882 sage: actual = (x*y).to_matrix()
1883 sage: X = x.to_matrix()
1884 sage: Y = y.to_matrix()
1885 sage: expected = (X*Y + Y*X)/2
1886 sage: actual == expected
1888 sage: J(expected) == x*y
1891 We can change the generator prefix::
1893 sage: RealSymmetricEJA(3, prefix='q').gens()
1894 (q0, q1, q2, q3, q4, q5)
1896 We can construct the (trivial) algebra of rank zero::
1898 sage: RealSymmetricEJA(0)
1899 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
1903 def _denormalized_basis(cls
, n
):
1905 Return a basis for the space of real symmetric n-by-n matrices.
1909 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1913 sage: set_random_seed()
1914 sage: n = ZZ.random_element(1,5)
1915 sage: B = RealSymmetricEJA._denormalized_basis(n)
1916 sage: all( M.is_symmetric() for M in B)
1920 # The basis of symmetric matrices, as matrices, in their R^(n-by-n)
1924 for j
in range(i
+1):
1925 Eij
= matrix(ZZ
, n
, lambda k
,l
: k
==i
and l
==j
)
1929 Sij
= Eij
+ Eij
.transpose()
1935 def _max_random_instance_size():
1936 return 4 # Dimension 10
1939 def random_instance(cls
, **kwargs
):
1941 Return a random instance of this type of algebra.
1943 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
1944 return cls(n
, **kwargs
)
1946 def __init__(self
, n
, **kwargs
):
1947 # We know this is a valid EJA, but will double-check
1948 # if the user passes check_axioms=True.
1949 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
1951 super(RealSymmetricEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
1952 self
.jordan_product
,
1953 self
.trace_inner_product
,
1956 # TODO: this could be factored out somehow, but is left here
1957 # because the MatrixEJA is not presently a subclass of the
1958 # FDEJA class that defines rank() and one().
1959 self
.rank
.set_cache(n
)
1960 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
1961 self
.one
.set_cache(self(idV
))
1965 class ComplexMatrixEJA(MatrixEJA
):
1966 # A manual dictionary-cache for the complex_extension() method,
1967 # since apparently @classmethods can't also be @cached_methods.
1968 _complex_extension
= {}
1971 def complex_extension(cls
,field
):
1973 The complex field that we embed/unembed, as an extension
1974 of the given ``field``.
1976 if field
in cls
._complex
_extension
:
1977 return cls
._complex
_extension
[field
]
1979 # Sage doesn't know how to adjoin the complex "i" (the root of
1980 # x^2 + 1) to a field in a general way. Here, we just enumerate
1981 # all of the cases that I have cared to support so far.
1983 # Sage doesn't know how to embed AA into QQbar, i.e. how
1984 # to adjoin sqrt(-1) to AA.
1986 elif not field
.is_exact():
1988 F
= field
.complex_field()
1990 # Works for QQ and... maybe some other fields.
1991 R
= PolynomialRing(field
, 'z')
1993 F
= field
.extension(z
**2 + 1, 'I', embedding
=CLF(-1).sqrt())
1995 cls
._complex
_extension
[field
] = F
1999 def dimension_over_reals():
2003 def real_embed(cls
,M
):
2005 Embed the n-by-n complex matrix ``M`` into the space of real
2006 matrices of size 2n-by-2n via the map the sends each entry `z = a +
2007 bi` to the block matrix ``[[a,b],[-b,a]]``.
2011 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2015 sage: F = QuadraticField(-1, 'I')
2016 sage: x1 = F(4 - 2*i)
2017 sage: x2 = F(1 + 2*i)
2020 sage: M = matrix(F,2,[[x1,x2],[x3,x4]])
2021 sage: ComplexMatrixEJA.real_embed(M)
2030 Embedding is a homomorphism (isomorphism, in fact)::
2032 sage: set_random_seed()
2033 sage: n = ZZ.random_element(3)
2034 sage: F = QuadraticField(-1, 'I')
2035 sage: X = random_matrix(F, n)
2036 sage: Y = random_matrix(F, n)
2037 sage: Xe = ComplexMatrixEJA.real_embed(X)
2038 sage: Ye = ComplexMatrixEJA.real_embed(Y)
2039 sage: XYe = ComplexMatrixEJA.real_embed(X*Y)
2044 super(ComplexMatrixEJA
,cls
).real_embed(M
)
2047 # We don't need any adjoined elements...
2048 field
= M
.base_ring().base_ring()
2054 blocks
.append(matrix(field
, 2, [ [ a
, b
],
2057 return matrix
.block(field
, n
, blocks
)
2061 def real_unembed(cls
,M
):
2063 The inverse of _embed_complex_matrix().
2067 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2071 sage: A = matrix(QQ,[ [ 1, 2, 3, 4],
2072 ....: [-2, 1, -4, 3],
2073 ....: [ 9, 10, 11, 12],
2074 ....: [-10, 9, -12, 11] ])
2075 sage: ComplexMatrixEJA.real_unembed(A)
2077 [ 10*I + 9 12*I + 11]
2081 Unembedding is the inverse of embedding::
2083 sage: set_random_seed()
2084 sage: F = QuadraticField(-1, 'I')
2085 sage: M = random_matrix(F, 3)
2086 sage: Me = ComplexMatrixEJA.real_embed(M)
2087 sage: ComplexMatrixEJA.real_unembed(Me) == M
2091 super(ComplexMatrixEJA
,cls
).real_unembed(M
)
2093 d
= cls
.dimension_over_reals()
2094 F
= cls
.complex_extension(M
.base_ring())
2097 # Go top-left to bottom-right (reading order), converting every
2098 # 2-by-2 block we see to a single complex element.
2100 for k
in range(n
/d
):
2101 for j
in range(n
/d
):
2102 submat
= M
[d
*k
:d
*k
+d
,d
*j
:d
*j
+d
]
2103 if submat
[0,0] != submat
[1,1]:
2104 raise ValueError('bad on-diagonal submatrix')
2105 if submat
[0,1] != -submat
[1,0]:
2106 raise ValueError('bad off-diagonal submatrix')
2107 z
= submat
[0,0] + submat
[0,1]*i
2110 return matrix(F
, n
/d
, elements
)
2113 class ComplexHermitianEJA(ConcreteEJA
, ComplexMatrixEJA
):
2115 The rank-n simple EJA consisting of complex Hermitian n-by-n
2116 matrices over the real numbers, the usual symmetric Jordan product,
2117 and the real-part-of-trace inner product. It has dimension `n^2` over
2122 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2126 In theory, our "field" can be any subfield of the reals::
2128 sage: ComplexHermitianEJA(2, field=RDF, check_axioms=True)
2129 Euclidean Jordan algebra of dimension 4 over Real Double Field
2130 sage: ComplexHermitianEJA(2, field=RR, check_axioms=True)
2131 Euclidean Jordan algebra of dimension 4 over Real Field with
2132 53 bits of precision
2136 The dimension of this algebra is `n^2`::
2138 sage: set_random_seed()
2139 sage: n_max = ComplexHermitianEJA._max_random_instance_size()
2140 sage: n = ZZ.random_element(1, n_max)
2141 sage: J = ComplexHermitianEJA(n)
2142 sage: J.dimension() == n^2
2145 The Jordan multiplication is what we think it is::
2147 sage: set_random_seed()
2148 sage: J = ComplexHermitianEJA.random_instance()
2149 sage: x,y = J.random_elements(2)
2150 sage: actual = (x*y).to_matrix()
2151 sage: X = x.to_matrix()
2152 sage: Y = y.to_matrix()
2153 sage: expected = (X*Y + Y*X)/2
2154 sage: actual == expected
2156 sage: J(expected) == x*y
2159 We can change the generator prefix::
2161 sage: ComplexHermitianEJA(2, prefix='z').gens()
2164 We can construct the (trivial) algebra of rank zero::
2166 sage: ComplexHermitianEJA(0)
2167 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2172 def _denormalized_basis(cls
, n
):
2174 Returns a basis for the space of complex Hermitian n-by-n matrices.
2176 Why do we embed these? Basically, because all of numerical linear
2177 algebra assumes that you're working with vectors consisting of `n`
2178 entries from a field and scalars from the same field. There's no way
2179 to tell SageMath that (for example) the vectors contain complex
2180 numbers, while the scalar field is real.
2184 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2188 sage: set_random_seed()
2189 sage: n = ZZ.random_element(1,5)
2190 sage: B = ComplexHermitianEJA._denormalized_basis(n)
2191 sage: all( M.is_symmetric() for M in B)
2196 R
= PolynomialRing(field
, 'z')
2198 F
= field
.extension(z
**2 + 1, 'I')
2201 # This is like the symmetric case, but we need to be careful:
2203 # * We want conjugate-symmetry, not just symmetry.
2204 # * The diagonal will (as a result) be real.
2207 Eij
= matrix
.zero(F
,n
)
2209 for j
in range(i
+1):
2213 Sij
= cls
.real_embed(Eij
)
2216 # The second one has a minus because it's conjugated.
2217 Eij
[j
,i
] = 1 # Eij = Eij + Eij.transpose()
2218 Sij_real
= cls
.real_embed(Eij
)
2220 # Eij = I*Eij - I*Eij.transpose()
2223 Sij_imag
= cls
.real_embed(Eij
)
2229 # Since we embedded these, we can drop back to the "field" that we
2230 # started with instead of the complex extension "F".
2231 return tuple( s
.change_ring(field
) for s
in S
)
2234 def __init__(self
, n
, **kwargs
):
2235 # We know this is a valid EJA, but will double-check
2236 # if the user passes check_axioms=True.
2237 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2239 super(ComplexHermitianEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
2240 self
.jordan_product
,
2241 self
.trace_inner_product
,
2243 # TODO: this could be factored out somehow, but is left here
2244 # because the MatrixEJA is not presently a subclass of the
2245 # FDEJA class that defines rank() and one().
2246 self
.rank
.set_cache(n
)
2247 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2248 self
.one
.set_cache(self(idV
))
2251 def _max_random_instance_size():
2252 return 3 # Dimension 9
2255 def random_instance(cls
, **kwargs
):
2257 Return a random instance of this type of algebra.
2259 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2260 return cls(n
, **kwargs
)
2262 class QuaternionMatrixEJA(MatrixEJA
):
2264 # A manual dictionary-cache for the quaternion_extension() method,
2265 # since apparently @classmethods can't also be @cached_methods.
2266 _quaternion_extension
= {}
2269 def quaternion_extension(cls
,field
):
2271 The quaternion field that we embed/unembed, as an extension
2272 of the given ``field``.
2274 if field
in cls
._quaternion
_extension
:
2275 return cls
._quaternion
_extension
[field
]
2277 Q
= QuaternionAlgebra(field
,-1,-1)
2279 cls
._quaternion
_extension
[field
] = Q
2283 def dimension_over_reals():
2287 def real_embed(cls
,M
):
2289 Embed the n-by-n quaternion matrix ``M`` into the space of real
2290 matrices of size 4n-by-4n by first sending each quaternion entry `z
2291 = a + bi + cj + dk` to the block-complex matrix ``[[a + bi,
2292 c+di],[-c + di, a-bi]]`, and then embedding those into a real
2297 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2301 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2302 sage: i,j,k = Q.gens()
2303 sage: x = 1 + 2*i + 3*j + 4*k
2304 sage: M = matrix(Q, 1, [[x]])
2305 sage: QuaternionMatrixEJA.real_embed(M)
2311 Embedding is a homomorphism (isomorphism, in fact)::
2313 sage: set_random_seed()
2314 sage: n = ZZ.random_element(2)
2315 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2316 sage: X = random_matrix(Q, n)
2317 sage: Y = random_matrix(Q, n)
2318 sage: Xe = QuaternionMatrixEJA.real_embed(X)
2319 sage: Ye = QuaternionMatrixEJA.real_embed(Y)
2320 sage: XYe = QuaternionMatrixEJA.real_embed(X*Y)
2325 super(QuaternionMatrixEJA
,cls
).real_embed(M
)
2326 quaternions
= M
.base_ring()
2329 F
= QuadraticField(-1, 'I')
2334 t
= z
.coefficient_tuple()
2339 cplxM
= matrix(F
, 2, [[ a
+ b
*i
, c
+ d
*i
],
2340 [-c
+ d
*i
, a
- b
*i
]])
2341 realM
= ComplexMatrixEJA
.real_embed(cplxM
)
2342 blocks
.append(realM
)
2344 # We should have real entries by now, so use the realest field
2345 # we've got for the return value.
2346 return matrix
.block(quaternions
.base_ring(), n
, blocks
)
2351 def real_unembed(cls
,M
):
2353 The inverse of _embed_quaternion_matrix().
2357 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2361 sage: M = matrix(QQ, [[ 1, 2, 3, 4],
2362 ....: [-2, 1, -4, 3],
2363 ....: [-3, 4, 1, -2],
2364 ....: [-4, -3, 2, 1]])
2365 sage: QuaternionMatrixEJA.real_unembed(M)
2366 [1 + 2*i + 3*j + 4*k]
2370 Unembedding is the inverse of embedding::
2372 sage: set_random_seed()
2373 sage: Q = QuaternionAlgebra(QQ, -1, -1)
2374 sage: M = random_matrix(Q, 3)
2375 sage: Me = QuaternionMatrixEJA.real_embed(M)
2376 sage: QuaternionMatrixEJA.real_unembed(Me) == M
2380 super(QuaternionMatrixEJA
,cls
).real_unembed(M
)
2382 d
= cls
.dimension_over_reals()
2384 # Use the base ring of the matrix to ensure that its entries can be
2385 # multiplied by elements of the quaternion algebra.
2386 Q
= cls
.quaternion_extension(M
.base_ring())
2389 # Go top-left to bottom-right (reading order), converting every
2390 # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
2393 for l
in range(n
/d
):
2394 for m
in range(n
/d
):
2395 submat
= ComplexMatrixEJA
.real_unembed(
2396 M
[d
*l
:d
*l
+d
,d
*m
:d
*m
+d
] )
2397 if submat
[0,0] != submat
[1,1].conjugate():
2398 raise ValueError('bad on-diagonal submatrix')
2399 if submat
[0,1] != -submat
[1,0].conjugate():
2400 raise ValueError('bad off-diagonal submatrix')
2401 z
= submat
[0,0].real()
2402 z
+= submat
[0,0].imag()*i
2403 z
+= submat
[0,1].real()*j
2404 z
+= submat
[0,1].imag()*k
2407 return matrix(Q
, n
/d
, elements
)
2410 class QuaternionHermitianEJA(ConcreteEJA
, QuaternionMatrixEJA
):
2412 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2413 matrices, the usual symmetric Jordan product, and the
2414 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2419 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2423 In theory, our "field" can be any subfield of the reals::
2425 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2426 Euclidean Jordan algebra of dimension 6 over Real Double Field
2427 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2428 Euclidean Jordan algebra of dimension 6 over Real Field with
2429 53 bits of precision
2433 The dimension of this algebra is `2*n^2 - n`::
2435 sage: set_random_seed()
2436 sage: n_max = QuaternionHermitianEJA._max_random_instance_size()
2437 sage: n = ZZ.random_element(1, n_max)
2438 sage: J = QuaternionHermitianEJA(n)
2439 sage: J.dimension() == 2*(n^2) - n
2442 The Jordan multiplication is what we think it is::
2444 sage: set_random_seed()
2445 sage: J = QuaternionHermitianEJA.random_instance()
2446 sage: x,y = J.random_elements(2)
2447 sage: actual = (x*y).to_matrix()
2448 sage: X = x.to_matrix()
2449 sage: Y = y.to_matrix()
2450 sage: expected = (X*Y + Y*X)/2
2451 sage: actual == expected
2453 sage: J(expected) == x*y
2456 We can change the generator prefix::
2458 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2459 (a0, a1, a2, a3, a4, a5)
2461 We can construct the (trivial) algebra of rank zero::
2463 sage: QuaternionHermitianEJA(0)
2464 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2468 def _denormalized_basis(cls
, n
):
2470 Returns a basis for the space of quaternion Hermitian n-by-n matrices.
2472 Why do we embed these? Basically, because all of numerical
2473 linear algebra assumes that you're working with vectors consisting
2474 of `n` entries from a field and scalars from the same field. There's
2475 no way to tell SageMath that (for example) the vectors contain
2476 complex numbers, while the scalar field is real.
2480 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2484 sage: set_random_seed()
2485 sage: n = ZZ.random_element(1,5)
2486 sage: B = QuaternionHermitianEJA._denormalized_basis(n)
2487 sage: all( M.is_symmetric() for M in B )
2492 Q
= QuaternionAlgebra(QQ
,-1,-1)
2495 # This is like the symmetric case, but we need to be careful:
2497 # * We want conjugate-symmetry, not just symmetry.
2498 # * The diagonal will (as a result) be real.
2501 Eij
= matrix
.zero(Q
,n
)
2503 for j
in range(i
+1):
2507 Sij
= cls
.real_embed(Eij
)
2510 # The second, third, and fourth ones have a minus
2511 # because they're conjugated.
2512 # Eij = Eij + Eij.transpose()
2514 Sij_real
= cls
.real_embed(Eij
)
2516 # Eij = I*(Eij - Eij.transpose())
2519 Sij_I
= cls
.real_embed(Eij
)
2521 # Eij = J*(Eij - Eij.transpose())
2524 Sij_J
= cls
.real_embed(Eij
)
2526 # Eij = K*(Eij - Eij.transpose())
2529 Sij_K
= cls
.real_embed(Eij
)
2535 # Since we embedded these, we can drop back to the "field" that we
2536 # started with instead of the quaternion algebra "Q".
2537 return tuple( s
.change_ring(field
) for s
in S
)
2540 def __init__(self
, n
, **kwargs
):
2541 # We know this is a valid EJA, but will double-check
2542 # if the user passes check_axioms=True.
2543 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2545 super(QuaternionHermitianEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
2546 self
.jordan_product
,
2547 self
.trace_inner_product
,
2549 # TODO: this could be factored out somehow, but is left here
2550 # because the MatrixEJA is not presently a subclass of the
2551 # FDEJA class that defines rank() and one().
2552 self
.rank
.set_cache(n
)
2553 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2554 self
.one
.set_cache(self(idV
))
2558 def _max_random_instance_size():
2560 The maximum rank of a random QuaternionHermitianEJA.
2562 return 2 # Dimension 6
2565 def random_instance(cls
, **kwargs
):
2567 Return a random instance of this type of algebra.
2569 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2570 return cls(n
, **kwargs
)
2573 class HadamardEJA(ConcreteEJA
):
2575 Return the Euclidean Jordan Algebra corresponding to the set
2576 `R^n` under the Hadamard product.
2578 Note: this is nothing more than the Cartesian product of ``n``
2579 copies of the spin algebra. Once Cartesian product algebras
2580 are implemented, this can go.
2584 sage: from mjo.eja.eja_algebra import HadamardEJA
2588 This multiplication table can be verified by hand::
2590 sage: J = HadamardEJA(3)
2591 sage: e0,e1,e2 = J.gens()
2607 We can change the generator prefix::
2609 sage: HadamardEJA(3, prefix='r').gens()
2613 def __init__(self
, n
, **kwargs
):
2615 jordan_product
= lambda x
,y
: x
2616 inner_product
= lambda x
,y
: x
2618 def jordan_product(x
,y
):
2620 return P( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2622 def inner_product(x
,y
):
2625 # New defaults for keyword arguments. Don't orthonormalize
2626 # because our basis is already orthonormal with respect to our
2627 # inner-product. Don't check the axioms, because we know this
2628 # is a valid EJA... but do double-check if the user passes
2629 # check_axioms=True. Note: we DON'T override the "check_field"
2630 # default here, because the user can pass in a field!
2631 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2632 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2634 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2635 super().__init
__(column_basis
,
2640 self
.rank
.set_cache(n
)
2643 self
.one
.set_cache( self
.zero() )
2645 self
.one
.set_cache( sum(self
.gens()) )
2648 def _max_random_instance_size():
2650 The maximum dimension of a random HadamardEJA.
2655 def random_instance(cls
, **kwargs
):
2657 Return a random instance of this type of algebra.
2659 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2660 return cls(n
, **kwargs
)
2663 class BilinearFormEJA(ConcreteEJA
):
2665 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2666 with the half-trace inner product and jordan product ``x*y =
2667 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2668 a symmetric positive-definite "bilinear form" matrix. Its
2669 dimension is the size of `B`, and it has rank two in dimensions
2670 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2671 the identity matrix of order ``n``.
2673 We insist that the one-by-one upper-left identity block of `B` be
2674 passed in as well so that we can be passed a matrix of size zero
2675 to construct a trivial algebra.
2679 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2680 ....: JordanSpinEJA)
2684 When no bilinear form is specified, the identity matrix is used,
2685 and the resulting algebra is the Jordan spin algebra::
2687 sage: B = matrix.identity(AA,3)
2688 sage: J0 = BilinearFormEJA(B)
2689 sage: J1 = JordanSpinEJA(3)
2690 sage: J0.multiplication_table() == J0.multiplication_table()
2693 An error is raised if the matrix `B` does not correspond to a
2694 positive-definite bilinear form::
2696 sage: B = matrix.random(QQ,2,3)
2697 sage: J = BilinearFormEJA(B)
2698 Traceback (most recent call last):
2700 ValueError: bilinear form is not positive-definite
2701 sage: B = matrix.zero(QQ,3)
2702 sage: J = BilinearFormEJA(B)
2703 Traceback (most recent call last):
2705 ValueError: bilinear form is not positive-definite
2709 We can create a zero-dimensional algebra::
2711 sage: B = matrix.identity(AA,0)
2712 sage: J = BilinearFormEJA(B)
2716 We can check the multiplication condition given in the Jordan, von
2717 Neumann, and Wigner paper (and also discussed on my "On the
2718 symmetry..." paper). Note that this relies heavily on the standard
2719 choice of basis, as does anything utilizing the bilinear form
2720 matrix. We opt not to orthonormalize the basis, because if we
2721 did, we would have to normalize the `s_{i}` in a similar manner::
2723 sage: set_random_seed()
2724 sage: n = ZZ.random_element(5)
2725 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2726 sage: B11 = matrix.identity(QQ,1)
2727 sage: B22 = M.transpose()*M
2728 sage: B = block_matrix(2,2,[ [B11,0 ],
2730 sage: J = BilinearFormEJA(B, orthonormalize=False)
2731 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2732 sage: V = J.vector_space()
2733 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2734 ....: for ei in eis ]
2735 sage: actual = [ sis[i]*sis[j]
2736 ....: for i in range(n-1)
2737 ....: for j in range(n-1) ]
2738 sage: expected = [ J.one() if i == j else J.zero()
2739 ....: for i in range(n-1)
2740 ....: for j in range(n-1) ]
2741 sage: actual == expected
2745 def __init__(self
, B
, **kwargs
):
2746 # The matrix "B" is supplied by the user in most cases,
2747 # so it makes sense to check whether or not its positive-
2748 # definite unless we are specifically asked not to...
2749 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2750 if not B
.is_positive_definite():
2751 raise ValueError("bilinear form is not positive-definite")
2753 # However, all of the other data for this EJA is computed
2754 # by us in manner that guarantees the axioms are
2755 # satisfied. So, again, unless we are specifically asked to
2756 # verify things, we'll skip the rest of the checks.
2757 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2759 def inner_product(x
,y
):
2760 return (y
.T
*B
*x
)[0,0]
2762 def jordan_product(x
,y
):
2768 z0
= inner_product(y
,x
)
2769 zbar
= y0
*xbar
+ x0
*ybar
2770 return P([z0
] + zbar
.list())
2773 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2774 super(BilinearFormEJA
, self
).__init
__(column_basis
,
2779 # The rank of this algebra is two, unless we're in a
2780 # one-dimensional ambient space (because the rank is bounded
2781 # by the ambient dimension).
2782 self
.rank
.set_cache(min(n
,2))
2785 self
.one
.set_cache( self
.zero() )
2787 self
.one
.set_cache( self
.monomial(0) )
2790 def _max_random_instance_size():
2792 The maximum dimension of a random BilinearFormEJA.
2797 def random_instance(cls
, **kwargs
):
2799 Return a random instance of this algebra.
2801 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2803 B
= matrix
.identity(ZZ
, n
)
2804 return cls(B
, **kwargs
)
2806 B11
= matrix
.identity(ZZ
, 1)
2807 M
= matrix
.random(ZZ
, n
-1)
2808 I
= matrix
.identity(ZZ
, n
-1)
2810 while alpha
.is_zero():
2811 alpha
= ZZ
.random_element().abs()
2812 B22
= M
.transpose()*M
+ alpha
*I
2814 from sage
.matrix
.special
import block_matrix
2815 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2818 return cls(B
, **kwargs
)
2821 class JordanSpinEJA(BilinearFormEJA
):
2823 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2824 with the usual inner product and jordan product ``x*y =
2825 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2830 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2834 This multiplication table can be verified by hand::
2836 sage: J = JordanSpinEJA(4)
2837 sage: e0,e1,e2,e3 = J.gens()
2853 We can change the generator prefix::
2855 sage: JordanSpinEJA(2, prefix='B').gens()
2860 Ensure that we have the usual inner product on `R^n`::
2862 sage: set_random_seed()
2863 sage: J = JordanSpinEJA.random_instance()
2864 sage: x,y = J.random_elements(2)
2865 sage: actual = x.inner_product(y)
2866 sage: expected = x.to_vector().inner_product(y.to_vector())
2867 sage: actual == expected
2871 def __init__(self
, n
, **kwargs
):
2872 # This is a special case of the BilinearFormEJA with the
2873 # identity matrix as its bilinear form.
2874 B
= matrix
.identity(ZZ
, n
)
2876 # Don't orthonormalize because our basis is already
2877 # orthonormal with respect to our inner-product.
2878 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2880 # But also don't pass check_field=False here, because the user
2881 # can pass in a field!
2882 super(JordanSpinEJA
, self
).__init
__(B
, **kwargs
)
2885 def _max_random_instance_size():
2887 The maximum dimension of a random JordanSpinEJA.
2892 def random_instance(cls
, **kwargs
):
2894 Return a random instance of this type of algebra.
2896 Needed here to override the implementation for ``BilinearFormEJA``.
2898 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2899 return cls(n
, **kwargs
)
2902 class TrivialEJA(ConcreteEJA
):
2904 The trivial Euclidean Jordan algebra consisting of only a zero element.
2908 sage: from mjo.eja.eja_algebra import TrivialEJA
2912 sage: J = TrivialEJA()
2919 sage: 7*J.one()*12*J.one()
2921 sage: J.one().inner_product(J.one())
2923 sage: J.one().norm()
2925 sage: J.one().subalgebra_generated_by()
2926 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2931 def __init__(self
, **kwargs
):
2932 jordan_product
= lambda x
,y
: x
2933 inner_product
= lambda x
,y
: 0
2936 # New defaults for keyword arguments
2937 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2938 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2940 super(TrivialEJA
, self
).__init
__(basis
,
2944 # The rank is zero using my definition, namely the dimension of the
2945 # largest subalgebra generated by any element.
2946 self
.rank
.set_cache(0)
2947 self
.one
.set_cache( self
.zero() )
2950 def random_instance(cls
, **kwargs
):
2951 # We don't take a "size" argument so the superclass method is
2952 # inappropriate for us.
2953 return cls(**kwargs
)
2956 class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct
,
2957 FiniteDimensionalEJA
):
2959 The external (orthogonal) direct sum of two or more Euclidean
2960 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
2961 orthogonal direct sum of simple Euclidean Jordan algebras which is
2962 then isometric to a Cartesian product, so no generality is lost by
2963 providing only this construction.
2967 sage: from mjo.eja.eja_algebra import (random_eja,
2968 ....: CartesianProductEJA,
2970 ....: JordanSpinEJA,
2971 ....: RealSymmetricEJA)
2975 The Jordan product is inherited from our factors and implemented by
2976 our CombinatorialFreeModule Cartesian product superclass::
2978 sage: set_random_seed()
2979 sage: J1 = HadamardEJA(2)
2980 sage: J2 = RealSymmetricEJA(2)
2981 sage: J = cartesian_product([J1,J2])
2982 sage: x,y = J.random_elements(2)
2986 The ability to retrieve the original factors is implemented by our
2987 CombinatorialFreeModule Cartesian product superclass::
2989 sage: J1 = HadamardEJA(2, field=QQ)
2990 sage: J2 = JordanSpinEJA(3, field=QQ)
2991 sage: J = cartesian_product([J1,J2])
2992 sage: J.cartesian_factors()
2993 (Euclidean Jordan algebra of dimension 2 over Rational Field,
2994 Euclidean Jordan algebra of dimension 3 over Rational Field)
2996 You can provide more than two factors::
2998 sage: J1 = HadamardEJA(2)
2999 sage: J2 = JordanSpinEJA(3)
3000 sage: J3 = RealSymmetricEJA(3)
3001 sage: cartesian_product([J1,J2,J3])
3002 Euclidean Jordan algebra of dimension 2 over Algebraic Real
3003 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
3004 Real Field (+) Euclidean Jordan algebra of dimension 6 over
3005 Algebraic Real Field
3007 Rank is additive on a Cartesian product::
3009 sage: J1 = HadamardEJA(1)
3010 sage: J2 = RealSymmetricEJA(2)
3011 sage: J = cartesian_product([J1,J2])
3012 sage: J1.rank.clear_cache()
3013 sage: J2.rank.clear_cache()
3014 sage: J.rank.clear_cache()
3017 sage: J.rank() == J1.rank() + J2.rank()
3020 The same rank computation works over the rationals, with whatever
3023 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
3024 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
3025 sage: J = cartesian_product([J1,J2])
3026 sage: J1.rank.clear_cache()
3027 sage: J2.rank.clear_cache()
3028 sage: J.rank.clear_cache()
3031 sage: J.rank() == J1.rank() + J2.rank()
3034 The product algebra will be associative if and only if all of its
3035 components are associative::
3037 sage: J1 = HadamardEJA(2)
3038 sage: J1.is_associative()
3040 sage: J2 = HadamardEJA(3)
3041 sage: J2.is_associative()
3043 sage: J3 = RealSymmetricEJA(3)
3044 sage: J3.is_associative()
3046 sage: CP1 = cartesian_product([J1,J2])
3047 sage: CP1.is_associative()
3049 sage: CP2 = cartesian_product([J1,J3])
3050 sage: CP2.is_associative()
3055 All factors must share the same base field::
3057 sage: J1 = HadamardEJA(2, field=QQ)
3058 sage: J2 = RealSymmetricEJA(2)
3059 sage: CartesianProductEJA((J1,J2))
3060 Traceback (most recent call last):
3062 ValueError: all factors must share the same base field
3064 The cached unit element is the same one that would be computed::
3066 sage: set_random_seed() # long time
3067 sage: J1 = random_eja() # long time
3068 sage: J2 = random_eja() # long time
3069 sage: J = cartesian_product([J1,J2]) # long time
3070 sage: actual = J.one() # long time
3071 sage: J.one.clear_cache() # long time
3072 sage: expected = J.one() # long time
3073 sage: actual == expected # long time
3077 Element
= FiniteDimensionalEJAElement
3080 def __init__(self
, algebras
, **kwargs
):
3081 CombinatorialFreeModule_CartesianProduct
.__init
__(self
,
3084 field
= algebras
[0].base_ring()
3085 if not all( J
.base_ring() == field
for J
in algebras
):
3086 raise ValueError("all factors must share the same base field")
3088 associative
= all( m
.is_associative() for m
in algebras
)
3090 # The definition of matrix_space() and self.basis() relies
3091 # only on the stuff in the CFM_CartesianProduct class, which
3092 # we've already initialized.
3093 Js
= self
.cartesian_factors()
3095 MS
= self
.matrix_space()
3097 MS(tuple( self
.cartesian_projection(i
)(b
).to_matrix()
3098 for i
in range(m
) ))
3099 for b
in self
.basis()
3102 # Define jordan/inner products that operate on that matrix_basis.
3103 def jordan_product(x
,y
):
3105 (Js
[i
](x
[i
])*Js
[i
](y
[i
])).to_matrix() for i
in range(m
)
3108 def inner_product(x
, y
):
3110 Js
[i
](x
[i
]).inner_product(Js
[i
](y
[i
])) for i
in range(m
)
3113 # There's no need to check the field since it already came
3114 # from an EJA. Likewise the axioms are guaranteed to be
3115 # satisfied, unless the guy writing this class sucks.
3117 # If you want the basis to be orthonormalized, orthonormalize
3119 FiniteDimensionalEJA
.__init
__(self
,
3124 orthonormalize
=False,
3125 associative
=associative
,
3126 cartesian_product
=True,
3130 ones
= tuple(J
.one() for J
in algebras
)
3131 self
.one
.set_cache(self
._cartesian
_product
_of
_elements
(ones
))
3132 self
.rank
.set_cache(sum(J
.rank() for J
in algebras
))
3134 def matrix_space(self
):
3136 Return the space that our matrix basis lives in as a Cartesian
3141 sage: from mjo.eja.eja_algebra import (HadamardEJA,
3142 ....: RealSymmetricEJA)
3146 sage: J1 = HadamardEJA(1)
3147 sage: J2 = RealSymmetricEJA(2)
3148 sage: J = cartesian_product([J1,J2])
3149 sage: J.matrix_space()
3150 The Cartesian product of (Full MatrixSpace of 1 by 1 dense
3151 matrices over Algebraic Real Field, Full MatrixSpace of 2
3152 by 2 dense matrices over Algebraic Real Field)
3155 from sage
.categories
.cartesian_product
import cartesian_product
3156 return cartesian_product( [J
.matrix_space()
3157 for J
in self
.cartesian_factors()] )
3160 def cartesian_projection(self
, i
):
3164 sage: from mjo.eja.eja_algebra import (random_eja,
3165 ....: JordanSpinEJA,
3167 ....: RealSymmetricEJA,
3168 ....: ComplexHermitianEJA)
3172 The projection morphisms are Euclidean Jordan algebra
3175 sage: J1 = HadamardEJA(2)
3176 sage: J2 = RealSymmetricEJA(2)
3177 sage: J = cartesian_product([J1,J2])
3178 sage: J.cartesian_projection(0)
3179 Linear operator between finite-dimensional Euclidean Jordan
3180 algebras represented by the matrix:
3183 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3184 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3185 Algebraic Real Field
3186 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3188 sage: J.cartesian_projection(1)
3189 Linear operator between finite-dimensional Euclidean Jordan
3190 algebras represented by the matrix:
3194 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3195 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3196 Algebraic Real Field
3197 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3200 The projections work the way you'd expect on the vector
3201 representation of an element::
3203 sage: J1 = JordanSpinEJA(2)
3204 sage: J2 = ComplexHermitianEJA(2)
3205 sage: J = cartesian_product([J1,J2])
3206 sage: pi_left = J.cartesian_projection(0)
3207 sage: pi_right = J.cartesian_projection(1)
3208 sage: pi_left(J.one()).to_vector()
3210 sage: pi_right(J.one()).to_vector()
3212 sage: J.one().to_vector()
3217 The answer never changes::
3219 sage: set_random_seed()
3220 sage: J1 = random_eja()
3221 sage: J2 = random_eja()
3222 sage: J = cartesian_product([J1,J2])
3223 sage: P0 = J.cartesian_projection(0)
3224 sage: P1 = J.cartesian_projection(0)
3229 Ji
= self
.cartesian_factors()[i
]
3230 # Requires the fix on Trac 31421/31422 to work!
3231 Pi
= super().cartesian_projection(i
)
3232 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3235 def cartesian_embedding(self
, i
):
3239 sage: from mjo.eja.eja_algebra import (random_eja,
3240 ....: JordanSpinEJA,
3242 ....: RealSymmetricEJA)
3246 The embedding morphisms are Euclidean Jordan algebra
3249 sage: J1 = HadamardEJA(2)
3250 sage: J2 = RealSymmetricEJA(2)
3251 sage: J = cartesian_product([J1,J2])
3252 sage: J.cartesian_embedding(0)
3253 Linear operator between finite-dimensional Euclidean Jordan
3254 algebras represented by the matrix:
3260 Domain: Euclidean Jordan algebra of dimension 2 over
3261 Algebraic Real Field
3262 Codomain: Euclidean Jordan algebra of dimension 2 over
3263 Algebraic Real Field (+) Euclidean Jordan algebra of
3264 dimension 3 over Algebraic Real Field
3265 sage: J.cartesian_embedding(1)
3266 Linear operator between finite-dimensional Euclidean Jordan
3267 algebras represented by the matrix:
3273 Domain: Euclidean Jordan algebra of dimension 3 over
3274 Algebraic Real Field
3275 Codomain: Euclidean Jordan algebra of dimension 2 over
3276 Algebraic Real Field (+) Euclidean Jordan algebra of
3277 dimension 3 over Algebraic Real Field
3279 The embeddings work the way you'd expect on the vector
3280 representation of an element::
3282 sage: J1 = JordanSpinEJA(3)
3283 sage: J2 = RealSymmetricEJA(2)
3284 sage: J = cartesian_product([J1,J2])
3285 sage: iota_left = J.cartesian_embedding(0)
3286 sage: iota_right = J.cartesian_embedding(1)
3287 sage: iota_left(J1.zero()) == J.zero()
3289 sage: iota_right(J2.zero()) == J.zero()
3291 sage: J1.one().to_vector()
3293 sage: iota_left(J1.one()).to_vector()
3295 sage: J2.one().to_vector()
3297 sage: iota_right(J2.one()).to_vector()
3299 sage: J.one().to_vector()
3304 The answer never changes::
3306 sage: set_random_seed()
3307 sage: J1 = random_eja()
3308 sage: J2 = random_eja()
3309 sage: J = cartesian_product([J1,J2])
3310 sage: E0 = J.cartesian_embedding(0)
3311 sage: E1 = J.cartesian_embedding(0)
3315 Composing a projection with the corresponding inclusion should
3316 produce the identity map, and mismatching them should produce
3319 sage: set_random_seed()
3320 sage: J1 = random_eja()
3321 sage: J2 = random_eja()
3322 sage: J = cartesian_product([J1,J2])
3323 sage: iota_left = J.cartesian_embedding(0)
3324 sage: iota_right = J.cartesian_embedding(1)
3325 sage: pi_left = J.cartesian_projection(0)
3326 sage: pi_right = J.cartesian_projection(1)
3327 sage: pi_left*iota_left == J1.one().operator()
3329 sage: pi_right*iota_right == J2.one().operator()
3331 sage: (pi_left*iota_right).is_zero()
3333 sage: (pi_right*iota_left).is_zero()
3337 Ji
= self
.cartesian_factors()[i
]
3338 # Requires the fix on Trac 31421/31422 to work!
3339 Ei
= super().cartesian_embedding(i
)
3340 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3344 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3346 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3349 A separate class for products of algebras for which we know a
3354 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
3355 ....: RealSymmetricEJA)
3359 This gives us fast characteristic polynomial computations in
3360 product algebras, too::
3363 sage: J1 = JordanSpinEJA(2)
3364 sage: J2 = RealSymmetricEJA(3)
3365 sage: J = cartesian_product([J1,J2])
3366 sage: J.characteristic_polynomial_of().degree()
3372 def __init__(self
, algebras
, **kwargs
):
3373 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3375 self
._rational
_algebra
= None
3376 if self
.vector_space().base_field() is not QQ
:
3377 self
._rational
_algebra
= cartesian_product([
3378 r
._rational
_algebra
for r
in algebras
3382 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3384 random_eja
= ConcreteEJA
.random_instance