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 of the two non-matrix
605 simple algebras (whose matrix representations are columns)
606 back and forth faithfully::
608 sage: set_random_seed()
609 sage: J = HadamardEJA.random_instance()
610 sage: x = J.random_element()
611 sage: J(x.to_vector().column()) == x
613 sage: J = JordanSpinEJA.random_instance()
614 sage: x = J.random_element()
615 sage: J(x.to_vector().column()) == x
618 We cannot coerce elements between algebras just because their
619 matrix representations are compatible::
621 sage: J1 = HadamardEJA(3)
622 sage: J2 = JordanSpinEJA(3)
624 Traceback (most recent call last):
626 ValueError: not an element of this algebra
628 Traceback (most recent call last):
630 ValueError: not an element of this algebra
633 msg
= "not an element of this algebra"
634 if elt
in self
.base_ring():
635 # Ensure that no base ring -> algebra coercion is performed
636 # by this method. There's some stupidity in sage that would
637 # otherwise propagate to this method; for example, sage thinks
638 # that the integer 3 belongs to the space of 2-by-2 matrices.
639 raise ValueError(msg
)
642 # Try to convert a vector into a column-matrix...
644 except (AttributeError, TypeError):
645 # and ignore failure, because we weren't really expecting
646 # a vector as an argument anyway.
649 if elt
not in self
.matrix_space():
650 raise ValueError(msg
)
652 # Thanks for nothing! Matrix spaces aren't vector spaces in
653 # Sage, so we have to figure out its matrix-basis coordinates
654 # ourselves. We use the basis space's ring instead of the
655 # element's ring because the basis space might be an algebraic
656 # closure whereas the base ring of the 3-by-3 identity matrix
657 # could be QQ instead of QQbar.
659 # And, we also have to handle Cartesian product bases (when
660 # the matrix basis consists of tuples) here. The "good news"
661 # is that we're already converting everything to long vectors,
662 # and that strategy works for tuples as well.
664 # We pass check=False because the matrix basis is "guaranteed"
665 # to be linearly independent... right? Ha ha.
667 V
= VectorSpace(self
.base_ring(), len(elt
))
668 W
= V
.span_of_basis( (V(_all2list(s
)) for s
in self
.matrix_basis()),
672 coords
= W
.coordinate_vector(V(elt
))
673 except ArithmeticError: # vector is not in free module
674 raise ValueError(msg
)
676 return self
.from_vector(coords
)
680 Return a string representation of ``self``.
684 sage: from mjo.eja.eja_algebra import JordanSpinEJA
688 Ensure that it says what we think it says::
690 sage: JordanSpinEJA(2, field=AA)
691 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
692 sage: JordanSpinEJA(3, field=RDF)
693 Euclidean Jordan algebra of dimension 3 over Real Double Field
696 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
697 return fmt
.format(self
.dimension(), self
.base_ring())
701 def characteristic_polynomial_of(self
):
703 Return the algebra's "characteristic polynomial of" function,
704 which is itself a multivariate polynomial that, when evaluated
705 at the coordinates of some algebra element, returns that
706 element's characteristic polynomial.
708 The resulting polynomial has `n+1` variables, where `n` is the
709 dimension of this algebra. The first `n` variables correspond to
710 the coordinates of an algebra element: when evaluated at the
711 coordinates of an algebra element with respect to a certain
712 basis, the result is a univariate polynomial (in the one
713 remaining variable ``t``), namely the characteristic polynomial
718 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
722 The characteristic polynomial in the spin algebra is given in
723 Alizadeh, Example 11.11::
725 sage: J = JordanSpinEJA(3)
726 sage: p = J.characteristic_polynomial_of(); p
727 X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
728 sage: xvec = J.one().to_vector()
732 By definition, the characteristic polynomial is a monic
733 degree-zero polynomial in a rank-zero algebra. Note that
734 Cayley-Hamilton is indeed satisfied since the polynomial
735 ``1`` evaluates to the identity element of the algebra on
738 sage: J = TrivialEJA()
739 sage: J.characteristic_polynomial_of()
746 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
747 a
= self
._charpoly
_coefficients
()
749 # We go to a bit of trouble here to reorder the
750 # indeterminates, so that it's easier to evaluate the
751 # characteristic polynomial at x's coordinates and get back
752 # something in terms of t, which is what we want.
753 S
= PolynomialRing(self
.base_ring(),'t')
757 S
= PolynomialRing(S
, R
.variable_names())
760 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
762 def coordinate_polynomial_ring(self
):
764 The multivariate polynomial ring in which this algebra's
765 :meth:`characteristic_polynomial_of` lives.
769 sage: from mjo.eja.eja_algebra import (HadamardEJA,
770 ....: RealSymmetricEJA)
774 sage: J = HadamardEJA(2)
775 sage: J.coordinate_polynomial_ring()
776 Multivariate Polynomial Ring in X1, X2...
777 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
778 sage: J.coordinate_polynomial_ring()
779 Multivariate Polynomial Ring in X1, X2, X3, X4, X5, X6...
782 var_names
= tuple( "X%d" % z
for z
in range(1, self
.dimension()+1) )
783 return PolynomialRing(self
.base_ring(), var_names
)
785 def inner_product(self
, x
, y
):
787 The inner product associated with this Euclidean Jordan algebra.
789 Defaults to the trace inner product, but can be overridden by
790 subclasses if they are sure that the necessary properties are
795 sage: from mjo.eja.eja_algebra import (random_eja,
797 ....: BilinearFormEJA)
801 Our inner product is "associative," which means the following for
802 a symmetric bilinear form::
804 sage: set_random_seed()
805 sage: J = random_eja()
806 sage: x,y,z = J.random_elements(3)
807 sage: (x*y).inner_product(z) == y.inner_product(x*z)
812 Ensure that this is the usual inner product for the algebras
815 sage: set_random_seed()
816 sage: J = HadamardEJA.random_instance()
817 sage: x,y = J.random_elements(2)
818 sage: actual = x.inner_product(y)
819 sage: expected = x.to_vector().inner_product(y.to_vector())
820 sage: actual == expected
823 Ensure that this is one-half of the trace inner-product in a
824 BilinearFormEJA that isn't just the reals (when ``n`` isn't
825 one). This is in Faraut and Koranyi, and also my "On the
828 sage: set_random_seed()
829 sage: J = BilinearFormEJA.random_instance()
830 sage: n = J.dimension()
831 sage: x = J.random_element()
832 sage: y = J.random_element()
833 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
836 B
= self
._inner
_product
_matrix
837 return (B
*x
.to_vector()).inner_product(y
.to_vector())
840 def is_trivial(self
):
842 Return whether or not this algebra is trivial.
844 A trivial algebra contains only the zero element.
848 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
853 sage: J = ComplexHermitianEJA(3)
859 sage: J = TrivialEJA()
864 return self
.dimension() == 0
867 def multiplication_table(self
):
869 Return a visual representation of this algebra's multiplication
870 table (on basis elements).
874 sage: from mjo.eja.eja_algebra import JordanSpinEJA
878 sage: J = JordanSpinEJA(4)
879 sage: J.multiplication_table()
880 +----++----+----+----+----+
881 | * || e0 | e1 | e2 | e3 |
882 +====++====+====+====+====+
883 | e0 || e0 | e1 | e2 | e3 |
884 +----++----+----+----+----+
885 | e1 || e1 | e0 | 0 | 0 |
886 +----++----+----+----+----+
887 | e2 || e2 | 0 | e0 | 0 |
888 +----++----+----+----+----+
889 | e3 || e3 | 0 | 0 | e0 |
890 +----++----+----+----+----+
894 # Prepend the header row.
895 M
= [["*"] + list(self
.gens())]
897 # And to each subsequent row, prepend an entry that belongs to
898 # the left-side "header column."
899 M
+= [ [self
.gens()[i
]] + [ self
.product_on_basis(i
,j
)
903 return table(M
, header_row
=True, header_column
=True, frame
=True)
906 def matrix_basis(self
):
908 Return an (often more natural) representation of this algebras
909 basis as an ordered tuple of matrices.
911 Every finite-dimensional Euclidean Jordan Algebra is a, up to
912 Jordan isomorphism, a direct sum of five simple
913 algebras---four of which comprise Hermitian matrices. And the
914 last type of algebra can of course be thought of as `n`-by-`1`
915 column matrices (ambiguusly called column vectors) to avoid
916 special cases. As a result, matrices (and column vectors) are
917 a natural representation format for Euclidean Jordan algebra
920 But, when we construct an algebra from a basis of matrices,
921 those matrix representations are lost in favor of coordinate
922 vectors *with respect to* that basis. We could eventually
923 convert back if we tried hard enough, but having the original
924 representations handy is valuable enough that we simply store
925 them and return them from this method.
927 Why implement this for non-matrix algebras? Avoiding special
928 cases for the :class:`BilinearFormEJA` pays with simplicity in
929 its own right. But mainly, we would like to be able to assume
930 that elements of a :class:`CartesianProductEJA` can be displayed
931 nicely, without having to have special classes for direct sums
932 one of whose components was a matrix algebra.
936 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
937 ....: RealSymmetricEJA)
941 sage: J = RealSymmetricEJA(2)
943 Finite family {0: e0, 1: e1, 2: e2}
944 sage: J.matrix_basis()
946 [1 0] [ 0 0.7071067811865475?] [0 0]
947 [0 0], [0.7071067811865475? 0], [0 1]
952 sage: J = JordanSpinEJA(2)
954 Finite family {0: e0, 1: e1}
955 sage: J.matrix_basis()
961 return self
._matrix
_basis
964 def matrix_space(self
):
966 Return the matrix space in which this algebra's elements live, if
967 we think of them as matrices (including column vectors of the
970 "By default" this will be an `n`-by-`1` column-matrix space,
971 except when the algebra is trivial. There it's `n`-by-`n`
972 (where `n` is zero), to ensure that two elements of the matrix
973 space (empty matrices) can be multiplied. For algebras of
974 matrices, this returns the space in which their
975 real embeddings live.
979 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
981 ....: QuaternionHermitianEJA,
986 By default, the matrix representation is just a column-matrix
987 equivalent to the vector representation::
989 sage: J = JordanSpinEJA(3)
990 sage: J.matrix_space()
991 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
994 The matrix representation in the trivial algebra is
995 zero-by-zero instead of the usual `n`-by-one::
997 sage: J = TrivialEJA()
998 sage: J.matrix_space()
999 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
1002 The matrix space for complex/quaternion Hermitian matrix EJA
1003 is the space in which their real-embeddings live, not the
1004 original complex/quaternion matrix space::
1006 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1007 sage: J.matrix_space()
1008 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1009 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1010 sage: J.matrix_space()
1011 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1014 if self
.is_trivial():
1015 return MatrixSpace(self
.base_ring(), 0)
1017 return self
.matrix_basis()[0].parent()
1023 Return the unit element of this algebra.
1027 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1032 We can compute unit element in the Hadamard EJA::
1034 sage: J = HadamardEJA(5)
1036 e0 + e1 + e2 + e3 + e4
1038 The unit element in the Hadamard EJA is inherited in the
1039 subalgebras generated by its elements::
1041 sage: J = HadamardEJA(5)
1043 e0 + e1 + e2 + e3 + e4
1044 sage: x = sum(J.gens())
1045 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1048 sage: A.one().superalgebra_element()
1049 e0 + e1 + e2 + e3 + e4
1053 The identity element acts like the identity, regardless of
1054 whether or not we orthonormalize::
1056 sage: set_random_seed()
1057 sage: J = random_eja()
1058 sage: x = J.random_element()
1059 sage: J.one()*x == x and x*J.one() == x
1061 sage: A = x.subalgebra_generated_by()
1062 sage: y = A.random_element()
1063 sage: A.one()*y == y and y*A.one() == y
1068 sage: set_random_seed()
1069 sage: J = random_eja(field=QQ, orthonormalize=False)
1070 sage: x = J.random_element()
1071 sage: J.one()*x == x and x*J.one() == x
1073 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1074 sage: y = A.random_element()
1075 sage: A.one()*y == y and y*A.one() == y
1078 The matrix of the unit element's operator is the identity,
1079 regardless of the base field and whether or not we
1082 sage: set_random_seed()
1083 sage: J = random_eja()
1084 sage: actual = J.one().operator().matrix()
1085 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1086 sage: actual == expected
1088 sage: x = J.random_element()
1089 sage: A = x.subalgebra_generated_by()
1090 sage: actual = A.one().operator().matrix()
1091 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1092 sage: actual == expected
1097 sage: set_random_seed()
1098 sage: J = random_eja(field=QQ, orthonormalize=False)
1099 sage: actual = J.one().operator().matrix()
1100 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1101 sage: actual == expected
1103 sage: x = J.random_element()
1104 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1105 sage: actual = A.one().operator().matrix()
1106 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1107 sage: actual == expected
1110 Ensure that the cached unit element (often precomputed by
1111 hand) agrees with the computed one::
1113 sage: set_random_seed()
1114 sage: J = random_eja()
1115 sage: cached = J.one()
1116 sage: J.one.clear_cache()
1117 sage: J.one() == cached
1122 sage: set_random_seed()
1123 sage: J = random_eja(field=QQ, orthonormalize=False)
1124 sage: cached = J.one()
1125 sage: J.one.clear_cache()
1126 sage: J.one() == cached
1130 # We can brute-force compute the matrices of the operators
1131 # that correspond to the basis elements of this algebra.
1132 # If some linear combination of those basis elements is the
1133 # algebra identity, then the same linear combination of
1134 # their matrices has to be the identity matrix.
1136 # Of course, matrices aren't vectors in sage, so we have to
1137 # appeal to the "long vectors" isometry.
1138 oper_vecs
= [ _mat2vec(g
.operator().matrix()) for g
in self
.gens() ]
1140 # Now we use basic linear algebra to find the coefficients,
1141 # of the matrices-as-vectors-linear-combination, which should
1142 # work for the original algebra basis too.
1143 A
= matrix(self
.base_ring(), oper_vecs
)
1145 # We used the isometry on the left-hand side already, but we
1146 # still need to do it for the right-hand side. Recall that we
1147 # wanted something that summed to the identity matrix.
1148 b
= _mat2vec( matrix
.identity(self
.base_ring(), self
.dimension()) )
1150 # Now if there's an identity element in the algebra, this
1151 # should work. We solve on the left to avoid having to
1152 # transpose the matrix "A".
1153 return self
.from_vector(A
.solve_left(b
))
1156 def peirce_decomposition(self
, c
):
1158 The Peirce decomposition of this algebra relative to the
1161 In the future, this can be extended to a complete system of
1162 orthogonal idempotents.
1166 - ``c`` -- an idempotent of this algebra.
1170 A triple (J0, J5, J1) containing two subalgebras and one subspace
1173 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1174 corresponding to the eigenvalue zero.
1176 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1177 corresponding to the eigenvalue one-half.
1179 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1180 corresponding to the eigenvalue one.
1182 These are the only possible eigenspaces for that operator, and this
1183 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1184 orthogonal, and are subalgebras of this algebra with the appropriate
1189 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1193 The canonical example comes from the symmetric matrices, which
1194 decompose into diagonal and off-diagonal parts::
1196 sage: J = RealSymmetricEJA(3)
1197 sage: C = matrix(QQ, [ [1,0,0],
1201 sage: J0,J5,J1 = J.peirce_decomposition(c)
1203 Euclidean Jordan algebra of dimension 1...
1205 Vector space of degree 6 and dimension 2...
1207 Euclidean Jordan algebra of dimension 3...
1208 sage: J0.one().to_matrix()
1212 sage: orig_df = AA.options.display_format
1213 sage: AA.options.display_format = 'radical'
1214 sage: J.from_vector(J5.basis()[0]).to_matrix()
1218 sage: J.from_vector(J5.basis()[1]).to_matrix()
1222 sage: AA.options.display_format = orig_df
1223 sage: J1.one().to_matrix()
1230 Every algebra decomposes trivially with respect to its identity
1233 sage: set_random_seed()
1234 sage: J = random_eja()
1235 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1236 sage: J0.dimension() == 0 and J5.dimension() == 0
1238 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1241 The decomposition is into eigenspaces, and its components are
1242 therefore necessarily orthogonal. Moreover, the identity
1243 elements in the two subalgebras are the projections onto their
1244 respective subspaces of the superalgebra's identity element::
1246 sage: set_random_seed()
1247 sage: J = random_eja()
1248 sage: x = J.random_element()
1249 sage: if not J.is_trivial():
1250 ....: while x.is_nilpotent():
1251 ....: x = J.random_element()
1252 sage: c = x.subalgebra_idempotent()
1253 sage: J0,J5,J1 = J.peirce_decomposition(c)
1255 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1256 ....: w = w.superalgebra_element()
1257 ....: y = J.from_vector(y)
1258 ....: z = z.superalgebra_element()
1259 ....: ipsum += w.inner_product(y).abs()
1260 ....: ipsum += w.inner_product(z).abs()
1261 ....: ipsum += y.inner_product(z).abs()
1264 sage: J1(c) == J1.one()
1266 sage: J0(J.one() - c) == J0.one()
1270 if not c
.is_idempotent():
1271 raise ValueError("element is not idempotent: %s" % c
)
1273 # Default these to what they should be if they turn out to be
1274 # trivial, because eigenspaces_left() won't return eigenvalues
1275 # corresponding to trivial spaces (e.g. it returns only the
1276 # eigenspace corresponding to lambda=1 if you take the
1277 # decomposition relative to the identity element).
1278 trivial
= self
.subalgebra(())
1279 J0
= trivial
# eigenvalue zero
1280 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1281 J1
= trivial
# eigenvalue one
1283 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1284 if eigval
== ~
(self
.base_ring()(2)):
1287 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1288 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1294 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1299 def random_element(self
, thorough
=False):
1301 Return a random element of this algebra.
1303 Our algebra superclass method only returns a linear
1304 combination of at most two basis elements. We instead
1305 want the vector space "random element" method that
1306 returns a more diverse selection.
1310 - ``thorough`` -- (boolean; default False) whether or not we
1311 should generate irrational coefficients for the random
1312 element when our base ring is irrational; this slows the
1313 algebra operations to a crawl, but any truly random method
1317 # For a general base ring... maybe we can trust this to do the
1318 # right thing? Unlikely, but.
1319 V
= self
.vector_space()
1320 v
= V
.random_element()
1322 if self
.base_ring() is AA
:
1323 # The "random element" method of the algebraic reals is
1324 # stupid at the moment, and only returns integers between
1325 # -2 and 2, inclusive:
1327 # https://trac.sagemath.org/ticket/30875
1329 # Instead, we implement our own "random vector" method,
1330 # and then coerce that into the algebra. We use the vector
1331 # space degree here instead of the dimension because a
1332 # subalgebra could (for example) be spanned by only two
1333 # vectors, each with five coordinates. We need to
1334 # generate all five coordinates.
1336 v
*= QQbar
.random_element().real()
1338 v
*= QQ
.random_element()
1340 return self
.from_vector(V
.coordinate_vector(v
))
1342 def random_elements(self
, count
, thorough
=False):
1344 Return ``count`` random elements as a tuple.
1348 - ``thorough`` -- (boolean; default False) whether or not we
1349 should generate irrational coefficients for the random
1350 elements when our base ring is irrational; this slows the
1351 algebra operations to a crawl, but any truly random method
1356 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1360 sage: J = JordanSpinEJA(3)
1361 sage: x,y,z = J.random_elements(3)
1362 sage: all( [ x in J, y in J, z in J ])
1364 sage: len( J.random_elements(10) ) == 10
1368 return tuple( self
.random_element(thorough
)
1369 for idx
in range(count
) )
1373 def _charpoly_coefficients(self
):
1375 The `r` polynomial coefficients of the "characteristic polynomial
1380 sage: from mjo.eja.eja_algebra import random_eja
1384 The theory shows that these are all homogeneous polynomials of
1387 sage: set_random_seed()
1388 sage: J = random_eja()
1389 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1393 n
= self
.dimension()
1394 R
= self
.coordinate_polynomial_ring()
1396 F
= R
.fraction_field()
1399 # From a result in my book, these are the entries of the
1400 # basis representation of L_x.
1401 return sum( vars[k
]*self
.gens()[k
].operator().matrix()[i
,j
]
1404 L_x
= matrix(F
, n
, n
, L_x_i_j
)
1407 if self
.rank
.is_in_cache():
1409 # There's no need to pad the system with redundant
1410 # columns if we *know* they'll be redundant.
1413 # Compute an extra power in case the rank is equal to
1414 # the dimension (otherwise, we would stop at x^(r-1)).
1415 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1416 for k
in range(n
+1) ]
1417 A
= matrix
.column(F
, x_powers
[:n
])
1418 AE
= A
.extended_echelon_form()
1425 # The theory says that only the first "r" coefficients are
1426 # nonzero, and they actually live in the original polynomial
1427 # ring and not the fraction field. We negate them because in
1428 # the actual characteristic polynomial, they get moved to the
1429 # other side where x^r lives. We don't bother to trim A_rref
1430 # down to a square matrix and solve the resulting system,
1431 # because the upper-left r-by-r portion of A_rref is
1432 # guaranteed to be the identity matrix, so e.g.
1434 # A_rref.solve_right(Y)
1436 # would just be returning Y.
1437 return (-E
*b
)[:r
].change_ring(R
)
1442 Return the rank of this EJA.
1444 This is a cached method because we know the rank a priori for
1445 all of the algebras we can construct. Thus we can avoid the
1446 expensive ``_charpoly_coefficients()`` call unless we truly
1447 need to compute the whole characteristic polynomial.
1451 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1452 ....: JordanSpinEJA,
1453 ....: RealSymmetricEJA,
1454 ....: ComplexHermitianEJA,
1455 ....: QuaternionHermitianEJA,
1460 The rank of the Jordan spin algebra is always two::
1462 sage: JordanSpinEJA(2).rank()
1464 sage: JordanSpinEJA(3).rank()
1466 sage: JordanSpinEJA(4).rank()
1469 The rank of the `n`-by-`n` Hermitian real, complex, or
1470 quaternion matrices is `n`::
1472 sage: RealSymmetricEJA(4).rank()
1474 sage: ComplexHermitianEJA(3).rank()
1476 sage: QuaternionHermitianEJA(2).rank()
1481 Ensure that every EJA that we know how to construct has a
1482 positive integer rank, unless the algebra is trivial in
1483 which case its rank will be zero::
1485 sage: set_random_seed()
1486 sage: J = random_eja()
1490 sage: r > 0 or (r == 0 and J.is_trivial())
1493 Ensure that computing the rank actually works, since the ranks
1494 of all simple algebras are known and will be cached by default::
1496 sage: set_random_seed() # long time
1497 sage: J = random_eja() # long time
1498 sage: cached = J.rank() # long time
1499 sage: J.rank.clear_cache() # long time
1500 sage: J.rank() == cached # long time
1504 return len(self
._charpoly
_coefficients
())
1507 def subalgebra(self
, basis
, **kwargs
):
1509 Create a subalgebra of this algebra from the given basis.
1511 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1512 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1515 def vector_space(self
):
1517 Return the vector space that underlies this algebra.
1521 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1525 sage: J = RealSymmetricEJA(2)
1526 sage: J.vector_space()
1527 Vector space of dimension 3 over...
1530 return self
.zero().to_vector().parent().ambient_vector_space()
1534 class RationalBasisEJA(FiniteDimensionalEJA
):
1536 New class for algebras whose supplied basis elements have all rational entries.
1540 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1544 The supplied basis is orthonormalized by default::
1546 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1547 sage: J = BilinearFormEJA(B)
1548 sage: J.matrix_basis()
1565 # Abuse the check_field parameter to check that the entries of
1566 # out basis (in ambient coordinates) are in the field QQ.
1567 if not all( all(b_i
in QQ
for b_i
in b
.list()) for b
in basis
):
1568 raise TypeError("basis not rational")
1570 self
._rational
_algebra
= None
1572 # There's no point in constructing the extra algebra if this
1573 # one is already rational.
1575 # Note: the same Jordan and inner-products work here,
1576 # because they are necessarily defined with respect to
1577 # ambient coordinates and not any particular basis.
1578 self
._rational
_algebra
= FiniteDimensionalEJA(
1583 orthonormalize
=False,
1587 super().__init
__(basis
,
1591 check_field
=check_field
,
1595 def _charpoly_coefficients(self
):
1599 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1600 ....: JordanSpinEJA)
1604 The base ring of the resulting polynomial coefficients is what
1605 it should be, and not the rationals (unless the algebra was
1606 already over the rationals)::
1608 sage: J = JordanSpinEJA(3)
1609 sage: J._charpoly_coefficients()
1610 (X1^2 - X2^2 - X3^2, -2*X1)
1611 sage: a0 = J._charpoly_coefficients()[0]
1613 Algebraic Real Field
1614 sage: a0.base_ring()
1615 Algebraic Real Field
1618 if self
._rational
_algebra
is None:
1619 # There's no need to construct *another* algebra over the
1620 # rationals if this one is already over the
1621 # rationals. Likewise, if we never orthonormalized our
1622 # basis, we might as well just use the given one.
1623 return super()._charpoly
_coefficients
()
1625 # Do the computation over the rationals. The answer will be
1626 # the same, because all we've done is a change of basis.
1627 # Then, change back from QQ to our real base ring
1628 a
= ( a_i
.change_ring(self
.base_ring())
1629 for a_i
in self
._rational
_algebra
._charpoly
_coefficients
() )
1631 if self
._deortho
_matrix
is None:
1632 # This can happen if our base ring was, say, AA and we
1633 # chose not to (or didn't need to) orthonormalize. It's
1634 # still faster to do the computations over QQ even if
1635 # the numbers in the boxes stay the same.
1638 # Otherwise, convert the coordinate variables back to the
1639 # deorthonormalized ones.
1640 R
= self
.coordinate_polynomial_ring()
1641 from sage
.modules
.free_module_element
import vector
1642 X
= vector(R
, R
.gens())
1643 BX
= self
._deortho
_matrix
*X
1645 subs_dict
= { X[i]: BX[i] for i in range(len(X)) }
1646 return tuple( a_i
.subs(subs_dict
) for a_i
in a
)
1648 class ConcreteEJA(RationalBasisEJA
):
1650 A class for the Euclidean Jordan algebras that we know by name.
1652 These are the Jordan algebras whose basis, multiplication table,
1653 rank, and so on are known a priori. More to the point, they are
1654 the Euclidean Jordan algebras for which we are able to conjure up
1655 a "random instance."
1659 sage: from mjo.eja.eja_algebra import ConcreteEJA
1663 Our basis is normalized with respect to the algebra's inner
1664 product, unless we specify otherwise::
1666 sage: set_random_seed()
1667 sage: J = ConcreteEJA.random_instance()
1668 sage: all( b.norm() == 1 for b in J.gens() )
1671 Since our basis is orthonormal with respect to the algebra's inner
1672 product, and since we know that this algebra is an EJA, any
1673 left-multiplication operator's matrix will be symmetric because
1674 natural->EJA basis representation is an isometry and within the
1675 EJA the operator is self-adjoint by the Jordan axiom::
1677 sage: set_random_seed()
1678 sage: J = ConcreteEJA.random_instance()
1679 sage: x = J.random_element()
1680 sage: x.operator().is_self_adjoint()
1685 def _max_random_instance_size():
1687 Return an integer "size" that is an upper bound on the size of
1688 this algebra when it is used in a random test
1689 case. Unfortunately, the term "size" is ambiguous -- when
1690 dealing with `R^n` under either the Hadamard or Jordan spin
1691 product, the "size" refers to the dimension `n`. When dealing
1692 with a matrix algebra (real symmetric or complex/quaternion
1693 Hermitian), it refers to the size of the matrix, which is far
1694 less than the dimension of the underlying vector space.
1696 This method must be implemented in each subclass.
1698 raise NotImplementedError
1701 def random_instance(cls
, *args
, **kwargs
):
1703 Return a random instance of this type of algebra.
1705 This method should be implemented in each subclass.
1707 from sage
.misc
.prandom
import choice
1708 eja_class
= choice(cls
.__subclasses
__())
1710 # These all bubble up to the RationalBasisEJA superclass
1711 # constructor, so any (kw)args valid there are also valid
1713 return eja_class
.random_instance(*args
, **kwargs
)
1718 def dimension_over_reals():
1720 The dimension of this matrix's base ring over the reals.
1722 The reals are dimension one over themselves, obviously; that's
1723 just `\mathbb{R}^{1}`. Likewise, the complex numbers `a + bi`
1724 have dimension two. Finally, the quaternions have dimension
1725 four over the reals.
1727 This is used to determine the size of the matrix returned from
1728 :meth:`real_embed`, among other things.
1730 raise NotImplementedError
1733 def real_embed(cls
,M
):
1735 Embed the matrix ``M`` into a space of real matrices.
1737 The matrix ``M`` can have entries in any field at the moment:
1738 the real numbers, complex numbers, or quaternions. And although
1739 they are not a field, we can probably support octonions at some
1740 point, too. This function returns a real matrix that "acts like"
1741 the original with respect to matrix multiplication; i.e.
1743 real_embed(M*N) = real_embed(M)*real_embed(N)
1746 if M
.ncols() != M
.nrows():
1747 raise ValueError("the matrix 'M' must be square")
1752 def real_unembed(cls
,M
):
1754 The inverse of :meth:`real_embed`.
1756 if M
.ncols() != M
.nrows():
1757 raise ValueError("the matrix 'M' must be square")
1758 if not ZZ(M
.nrows()).mod(cls
.dimension_over_reals()).is_zero():
1759 raise ValueError("the matrix 'M' must be a real embedding")
1763 def jordan_product(X
,Y
):
1764 return (X
*Y
+ Y
*X
)/2
1767 def trace_inner_product(cls
,X
,Y
):
1769 Compute the trace inner-product of two real-embeddings.
1773 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
1774 ....: ComplexHermitianEJA,
1775 ....: QuaternionHermitianEJA)
1779 This gives the same answer as it would if we computed the trace
1780 from the unembedded (original) matrices::
1782 sage: set_random_seed()
1783 sage: J = RealSymmetricEJA.random_instance()
1784 sage: x,y = J.random_elements(2)
1785 sage: Xe = x.to_matrix()
1786 sage: Ye = y.to_matrix()
1787 sage: X = J.real_unembed(Xe)
1788 sage: Y = J.real_unembed(Ye)
1789 sage: expected = (X*Y).trace()
1790 sage: actual = J.trace_inner_product(Xe,Ye)
1791 sage: actual == expected
1796 sage: set_random_seed()
1797 sage: J = ComplexHermitianEJA.random_instance()
1798 sage: x,y = J.random_elements(2)
1799 sage: Xe = x.to_matrix()
1800 sage: Ye = y.to_matrix()
1801 sage: X = J.real_unembed(Xe)
1802 sage: Y = J.real_unembed(Ye)
1803 sage: expected = (X*Y).trace().real()
1804 sage: actual = J.trace_inner_product(Xe,Ye)
1805 sage: actual == expected
1810 sage: set_random_seed()
1811 sage: J = QuaternionHermitianEJA.random_instance()
1812 sage: x,y = J.random_elements(2)
1813 sage: Xe = x.to_matrix()
1814 sage: Ye = y.to_matrix()
1815 sage: X = J.real_unembed(Xe)
1816 sage: Y = J.real_unembed(Ye)
1817 sage: expected = (X*Y).trace().coefficient_tuple()[0]
1818 sage: actual = J.trace_inner_product(Xe,Ye)
1819 sage: actual == expected
1823 Xu
= cls
.real_unembed(X
)
1824 Yu
= cls
.real_unembed(Y
)
1825 tr
= (Xu
*Yu
).trace()
1828 # Works in QQ, AA, RDF, et cetera.
1830 except AttributeError:
1831 # A quaternion doesn't have a real() method, but does
1832 # have coefficient_tuple() method that returns the
1833 # coefficients of 1, i, j, and k -- in that order.
1834 return tr
.coefficient_tuple()[0]
1837 class RealMatrixEJA(MatrixEJA
):
1839 def dimension_over_reals():
1843 class RealSymmetricEJA(ConcreteEJA
, RealMatrixEJA
):
1845 The rank-n simple EJA consisting of real symmetric n-by-n
1846 matrices, the usual symmetric Jordan product, and the trace inner
1847 product. It has dimension `(n^2 + n)/2` over the reals.
1851 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1855 sage: J = RealSymmetricEJA(2)
1856 sage: e0, e1, e2 = J.gens()
1864 In theory, our "field" can be any subfield of the reals::
1866 sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
1867 Euclidean Jordan algebra of dimension 3 over Real Double Field
1868 sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
1869 Euclidean Jordan algebra of dimension 3 over Real Field with
1870 53 bits of precision
1874 The dimension of this algebra is `(n^2 + n) / 2`::
1876 sage: set_random_seed()
1877 sage: n_max = RealSymmetricEJA._max_random_instance_size()
1878 sage: n = ZZ.random_element(1, n_max)
1879 sage: J = RealSymmetricEJA(n)
1880 sage: J.dimension() == (n^2 + n)/2
1883 The Jordan multiplication is what we think it is::
1885 sage: set_random_seed()
1886 sage: J = RealSymmetricEJA.random_instance()
1887 sage: x,y = J.random_elements(2)
1888 sage: actual = (x*y).to_matrix()
1889 sage: X = x.to_matrix()
1890 sage: Y = y.to_matrix()
1891 sage: expected = (X*Y + Y*X)/2
1892 sage: actual == expected
1894 sage: J(expected) == x*y
1897 We can change the generator prefix::
1899 sage: RealSymmetricEJA(3, prefix='q').gens()
1900 (q0, q1, q2, q3, q4, q5)
1902 We can construct the (trivial) algebra of rank zero::
1904 sage: RealSymmetricEJA(0)
1905 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
1909 def _denormalized_basis(cls
, n
):
1911 Return a basis for the space of real symmetric n-by-n matrices.
1915 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1919 sage: set_random_seed()
1920 sage: n = ZZ.random_element(1,5)
1921 sage: B = RealSymmetricEJA._denormalized_basis(n)
1922 sage: all( M.is_symmetric() for M in B)
1926 # The basis of symmetric matrices, as matrices, in their R^(n-by-n)
1930 for j
in range(i
+1):
1931 Eij
= matrix(ZZ
, n
, lambda k
,l
: k
==i
and l
==j
)
1935 Sij
= Eij
+ Eij
.transpose()
1941 def _max_random_instance_size():
1942 return 4 # Dimension 10
1945 def random_instance(cls
, **kwargs
):
1947 Return a random instance of this type of algebra.
1949 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
1950 return cls(n
, **kwargs
)
1952 def __init__(self
, n
, **kwargs
):
1953 # We know this is a valid EJA, but will double-check
1954 # if the user passes check_axioms=True.
1955 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
1957 super(RealSymmetricEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
1958 self
.jordan_product
,
1959 self
.trace_inner_product
,
1962 # TODO: this could be factored out somehow, but is left here
1963 # because the MatrixEJA is not presently a subclass of the
1964 # FDEJA class that defines rank() and one().
1965 self
.rank
.set_cache(n
)
1966 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
1967 self
.one
.set_cache(self(idV
))
1971 class ComplexMatrixEJA(MatrixEJA
):
1972 # A manual dictionary-cache for the complex_extension() method,
1973 # since apparently @classmethods can't also be @cached_methods.
1974 _complex_extension
= {}
1977 def complex_extension(cls
,field
):
1979 The complex field that we embed/unembed, as an extension
1980 of the given ``field``.
1982 if field
in cls
._complex
_extension
:
1983 return cls
._complex
_extension
[field
]
1985 # Sage doesn't know how to adjoin the complex "i" (the root of
1986 # x^2 + 1) to a field in a general way. Here, we just enumerate
1987 # all of the cases that I have cared to support so far.
1989 # Sage doesn't know how to embed AA into QQbar, i.e. how
1990 # to adjoin sqrt(-1) to AA.
1992 elif not field
.is_exact():
1994 F
= field
.complex_field()
1996 # Works for QQ and... maybe some other fields.
1997 R
= PolynomialRing(field
, 'z')
1999 F
= field
.extension(z
**2 + 1, 'I', embedding
=CLF(-1).sqrt())
2001 cls
._complex
_extension
[field
] = F
2005 def dimension_over_reals():
2009 def real_embed(cls
,M
):
2011 Embed the n-by-n complex matrix ``M`` into the space of real
2012 matrices of size 2n-by-2n via the map the sends each entry `z = a +
2013 bi` to the block matrix ``[[a,b],[-b,a]]``.
2017 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2021 sage: F = QuadraticField(-1, 'I')
2022 sage: x1 = F(4 - 2*i)
2023 sage: x2 = F(1 + 2*i)
2026 sage: M = matrix(F,2,[[x1,x2],[x3,x4]])
2027 sage: ComplexMatrixEJA.real_embed(M)
2036 Embedding is a homomorphism (isomorphism, in fact)::
2038 sage: set_random_seed()
2039 sage: n = ZZ.random_element(3)
2040 sage: F = QuadraticField(-1, 'I')
2041 sage: X = random_matrix(F, n)
2042 sage: Y = random_matrix(F, n)
2043 sage: Xe = ComplexMatrixEJA.real_embed(X)
2044 sage: Ye = ComplexMatrixEJA.real_embed(Y)
2045 sage: XYe = ComplexMatrixEJA.real_embed(X*Y)
2050 super(ComplexMatrixEJA
,cls
).real_embed(M
)
2053 # We don't need any adjoined elements...
2054 field
= M
.base_ring().base_ring()
2060 blocks
.append(matrix(field
, 2, [ [ a
, b
],
2063 return matrix
.block(field
, n
, blocks
)
2067 def real_unembed(cls
,M
):
2069 The inverse of _embed_complex_matrix().
2073 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2077 sage: A = matrix(QQ,[ [ 1, 2, 3, 4],
2078 ....: [-2, 1, -4, 3],
2079 ....: [ 9, 10, 11, 12],
2080 ....: [-10, 9, -12, 11] ])
2081 sage: ComplexMatrixEJA.real_unembed(A)
2083 [ 10*I + 9 12*I + 11]
2087 Unembedding is the inverse of embedding::
2089 sage: set_random_seed()
2090 sage: F = QuadraticField(-1, 'I')
2091 sage: M = random_matrix(F, 3)
2092 sage: Me = ComplexMatrixEJA.real_embed(M)
2093 sage: ComplexMatrixEJA.real_unembed(Me) == M
2097 super(ComplexMatrixEJA
,cls
).real_unembed(M
)
2099 d
= cls
.dimension_over_reals()
2100 F
= cls
.complex_extension(M
.base_ring())
2103 # Go top-left to bottom-right (reading order), converting every
2104 # 2-by-2 block we see to a single complex element.
2106 for k
in range(n
/d
):
2107 for j
in range(n
/d
):
2108 submat
= M
[d
*k
:d
*k
+d
,d
*j
:d
*j
+d
]
2109 if submat
[0,0] != submat
[1,1]:
2110 raise ValueError('bad on-diagonal submatrix')
2111 if submat
[0,1] != -submat
[1,0]:
2112 raise ValueError('bad off-diagonal submatrix')
2113 z
= submat
[0,0] + submat
[0,1]*i
2116 return matrix(F
, n
/d
, elements
)
2119 class ComplexHermitianEJA(ConcreteEJA
, ComplexMatrixEJA
):
2121 The rank-n simple EJA consisting of complex Hermitian n-by-n
2122 matrices over the real numbers, the usual symmetric Jordan product,
2123 and the real-part-of-trace inner product. It has dimension `n^2` over
2128 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2132 In theory, our "field" can be any subfield of the reals::
2134 sage: ComplexHermitianEJA(2, field=RDF, check_axioms=True)
2135 Euclidean Jordan algebra of dimension 4 over Real Double Field
2136 sage: ComplexHermitianEJA(2, field=RR, check_axioms=True)
2137 Euclidean Jordan algebra of dimension 4 over Real Field with
2138 53 bits of precision
2142 The dimension of this algebra is `n^2`::
2144 sage: set_random_seed()
2145 sage: n_max = ComplexHermitianEJA._max_random_instance_size()
2146 sage: n = ZZ.random_element(1, n_max)
2147 sage: J = ComplexHermitianEJA(n)
2148 sage: J.dimension() == n^2
2151 The Jordan multiplication is what we think it is::
2153 sage: set_random_seed()
2154 sage: J = ComplexHermitianEJA.random_instance()
2155 sage: x,y = J.random_elements(2)
2156 sage: actual = (x*y).to_matrix()
2157 sage: X = x.to_matrix()
2158 sage: Y = y.to_matrix()
2159 sage: expected = (X*Y + Y*X)/2
2160 sage: actual == expected
2162 sage: J(expected) == x*y
2165 We can change the generator prefix::
2167 sage: ComplexHermitianEJA(2, prefix='z').gens()
2170 We can construct the (trivial) algebra of rank zero::
2172 sage: ComplexHermitianEJA(0)
2173 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2178 def _denormalized_basis(cls
, n
):
2180 Returns a basis for the space of complex Hermitian n-by-n matrices.
2182 Why do we embed these? Basically, because all of numerical linear
2183 algebra assumes that you're working with vectors consisting of `n`
2184 entries from a field and scalars from the same field. There's no way
2185 to tell SageMath that (for example) the vectors contain complex
2186 numbers, while the scalar field is real.
2190 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2194 sage: set_random_seed()
2195 sage: n = ZZ.random_element(1,5)
2196 sage: B = ComplexHermitianEJA._denormalized_basis(n)
2197 sage: all( M.is_symmetric() for M in B)
2202 R
= PolynomialRing(field
, 'z')
2204 F
= field
.extension(z
**2 + 1, 'I')
2207 # This is like the symmetric case, but we need to be careful:
2209 # * We want conjugate-symmetry, not just symmetry.
2210 # * The diagonal will (as a result) be real.
2213 Eij
= matrix
.zero(F
,n
)
2215 for j
in range(i
+1):
2219 Sij
= cls
.real_embed(Eij
)
2222 # The second one has a minus because it's conjugated.
2223 Eij
[j
,i
] = 1 # Eij = Eij + Eij.transpose()
2224 Sij_real
= cls
.real_embed(Eij
)
2226 # Eij = I*Eij - I*Eij.transpose()
2229 Sij_imag
= cls
.real_embed(Eij
)
2235 # Since we embedded these, we can drop back to the "field" that we
2236 # started with instead of the complex extension "F".
2237 return tuple( s
.change_ring(field
) for s
in S
)
2240 def __init__(self
, n
, **kwargs
):
2241 # We know this is a valid EJA, but will double-check
2242 # if the user passes check_axioms=True.
2243 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2245 super(ComplexHermitianEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
2246 self
.jordan_product
,
2247 self
.trace_inner_product
,
2249 # TODO: this could be factored out somehow, but is left here
2250 # because the MatrixEJA is not presently a subclass of the
2251 # FDEJA class that defines rank() and one().
2252 self
.rank
.set_cache(n
)
2253 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2254 self
.one
.set_cache(self(idV
))
2257 def _max_random_instance_size():
2258 return 3 # Dimension 9
2261 def random_instance(cls
, **kwargs
):
2263 Return a random instance of this type of algebra.
2265 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2266 return cls(n
, **kwargs
)
2268 class QuaternionMatrixEJA(MatrixEJA
):
2270 # A manual dictionary-cache for the quaternion_extension() method,
2271 # since apparently @classmethods can't also be @cached_methods.
2272 _quaternion_extension
= {}
2275 def quaternion_extension(cls
,field
):
2277 The quaternion field that we embed/unembed, as an extension
2278 of the given ``field``.
2280 if field
in cls
._quaternion
_extension
:
2281 return cls
._quaternion
_extension
[field
]
2283 Q
= QuaternionAlgebra(field
,-1,-1)
2285 cls
._quaternion
_extension
[field
] = Q
2289 def dimension_over_reals():
2293 def real_embed(cls
,M
):
2295 Embed the n-by-n quaternion matrix ``M`` into the space of real
2296 matrices of size 4n-by-4n by first sending each quaternion entry `z
2297 = a + bi + cj + dk` to the block-complex matrix ``[[a + bi,
2298 c+di],[-c + di, a-bi]]`, and then embedding those into a real
2303 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2307 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2308 sage: i,j,k = Q.gens()
2309 sage: x = 1 + 2*i + 3*j + 4*k
2310 sage: M = matrix(Q, 1, [[x]])
2311 sage: QuaternionMatrixEJA.real_embed(M)
2317 Embedding is a homomorphism (isomorphism, in fact)::
2319 sage: set_random_seed()
2320 sage: n = ZZ.random_element(2)
2321 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2322 sage: X = random_matrix(Q, n)
2323 sage: Y = random_matrix(Q, n)
2324 sage: Xe = QuaternionMatrixEJA.real_embed(X)
2325 sage: Ye = QuaternionMatrixEJA.real_embed(Y)
2326 sage: XYe = QuaternionMatrixEJA.real_embed(X*Y)
2331 super(QuaternionMatrixEJA
,cls
).real_embed(M
)
2332 quaternions
= M
.base_ring()
2335 F
= QuadraticField(-1, 'I')
2340 t
= z
.coefficient_tuple()
2345 cplxM
= matrix(F
, 2, [[ a
+ b
*i
, c
+ d
*i
],
2346 [-c
+ d
*i
, a
- b
*i
]])
2347 realM
= ComplexMatrixEJA
.real_embed(cplxM
)
2348 blocks
.append(realM
)
2350 # We should have real entries by now, so use the realest field
2351 # we've got for the return value.
2352 return matrix
.block(quaternions
.base_ring(), n
, blocks
)
2357 def real_unembed(cls
,M
):
2359 The inverse of _embed_quaternion_matrix().
2363 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2367 sage: M = matrix(QQ, [[ 1, 2, 3, 4],
2368 ....: [-2, 1, -4, 3],
2369 ....: [-3, 4, 1, -2],
2370 ....: [-4, -3, 2, 1]])
2371 sage: QuaternionMatrixEJA.real_unembed(M)
2372 [1 + 2*i + 3*j + 4*k]
2376 Unembedding is the inverse of embedding::
2378 sage: set_random_seed()
2379 sage: Q = QuaternionAlgebra(QQ, -1, -1)
2380 sage: M = random_matrix(Q, 3)
2381 sage: Me = QuaternionMatrixEJA.real_embed(M)
2382 sage: QuaternionMatrixEJA.real_unembed(Me) == M
2386 super(QuaternionMatrixEJA
,cls
).real_unembed(M
)
2388 d
= cls
.dimension_over_reals()
2390 # Use the base ring of the matrix to ensure that its entries can be
2391 # multiplied by elements of the quaternion algebra.
2392 Q
= cls
.quaternion_extension(M
.base_ring())
2395 # Go top-left to bottom-right (reading order), converting every
2396 # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
2399 for l
in range(n
/d
):
2400 for m
in range(n
/d
):
2401 submat
= ComplexMatrixEJA
.real_unembed(
2402 M
[d
*l
:d
*l
+d
,d
*m
:d
*m
+d
] )
2403 if submat
[0,0] != submat
[1,1].conjugate():
2404 raise ValueError('bad on-diagonal submatrix')
2405 if submat
[0,1] != -submat
[1,0].conjugate():
2406 raise ValueError('bad off-diagonal submatrix')
2407 z
= submat
[0,0].real()
2408 z
+= submat
[0,0].imag()*i
2409 z
+= submat
[0,1].real()*j
2410 z
+= submat
[0,1].imag()*k
2413 return matrix(Q
, n
/d
, elements
)
2416 class QuaternionHermitianEJA(ConcreteEJA
, QuaternionMatrixEJA
):
2418 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2419 matrices, the usual symmetric Jordan product, and the
2420 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2425 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2429 In theory, our "field" can be any subfield of the reals::
2431 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2432 Euclidean Jordan algebra of dimension 6 over Real Double Field
2433 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2434 Euclidean Jordan algebra of dimension 6 over Real Field with
2435 53 bits of precision
2439 The dimension of this algebra is `2*n^2 - n`::
2441 sage: set_random_seed()
2442 sage: n_max = QuaternionHermitianEJA._max_random_instance_size()
2443 sage: n = ZZ.random_element(1, n_max)
2444 sage: J = QuaternionHermitianEJA(n)
2445 sage: J.dimension() == 2*(n^2) - n
2448 The Jordan multiplication is what we think it is::
2450 sage: set_random_seed()
2451 sage: J = QuaternionHermitianEJA.random_instance()
2452 sage: x,y = J.random_elements(2)
2453 sage: actual = (x*y).to_matrix()
2454 sage: X = x.to_matrix()
2455 sage: Y = y.to_matrix()
2456 sage: expected = (X*Y + Y*X)/2
2457 sage: actual == expected
2459 sage: J(expected) == x*y
2462 We can change the generator prefix::
2464 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2465 (a0, a1, a2, a3, a4, a5)
2467 We can construct the (trivial) algebra of rank zero::
2469 sage: QuaternionHermitianEJA(0)
2470 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2474 def _denormalized_basis(cls
, n
):
2476 Returns a basis for the space of quaternion Hermitian n-by-n matrices.
2478 Why do we embed these? Basically, because all of numerical
2479 linear algebra assumes that you're working with vectors consisting
2480 of `n` entries from a field and scalars from the same field. There's
2481 no way to tell SageMath that (for example) the vectors contain
2482 complex numbers, while the scalar field is real.
2486 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2490 sage: set_random_seed()
2491 sage: n = ZZ.random_element(1,5)
2492 sage: B = QuaternionHermitianEJA._denormalized_basis(n)
2493 sage: all( M.is_symmetric() for M in B )
2498 Q
= QuaternionAlgebra(QQ
,-1,-1)
2501 # This is like the symmetric case, but we need to be careful:
2503 # * We want conjugate-symmetry, not just symmetry.
2504 # * The diagonal will (as a result) be real.
2507 Eij
= matrix
.zero(Q
,n
)
2509 for j
in range(i
+1):
2513 Sij
= cls
.real_embed(Eij
)
2516 # The second, third, and fourth ones have a minus
2517 # because they're conjugated.
2518 # Eij = Eij + Eij.transpose()
2520 Sij_real
= cls
.real_embed(Eij
)
2522 # Eij = I*(Eij - Eij.transpose())
2525 Sij_I
= cls
.real_embed(Eij
)
2527 # Eij = J*(Eij - Eij.transpose())
2530 Sij_J
= cls
.real_embed(Eij
)
2532 # Eij = K*(Eij - Eij.transpose())
2535 Sij_K
= cls
.real_embed(Eij
)
2541 # Since we embedded these, we can drop back to the "field" that we
2542 # started with instead of the quaternion algebra "Q".
2543 return tuple( s
.change_ring(field
) for s
in S
)
2546 def __init__(self
, n
, **kwargs
):
2547 # We know this is a valid EJA, but will double-check
2548 # if the user passes check_axioms=True.
2549 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2551 super(QuaternionHermitianEJA
, self
).__init
__(self
._denormalized
_basis
(n
),
2552 self
.jordan_product
,
2553 self
.trace_inner_product
,
2555 # TODO: this could be factored out somehow, but is left here
2556 # because the MatrixEJA is not presently a subclass of the
2557 # FDEJA class that defines rank() and one().
2558 self
.rank
.set_cache(n
)
2559 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2560 self
.one
.set_cache(self(idV
))
2564 def _max_random_instance_size():
2566 The maximum rank of a random QuaternionHermitianEJA.
2568 return 2 # Dimension 6
2571 def random_instance(cls
, **kwargs
):
2573 Return a random instance of this type of algebra.
2575 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2576 return cls(n
, **kwargs
)
2579 class HadamardEJA(ConcreteEJA
):
2581 Return the Euclidean Jordan Algebra corresponding to the set
2582 `R^n` under the Hadamard product.
2584 Note: this is nothing more than the Cartesian product of ``n``
2585 copies of the spin algebra. Once Cartesian product algebras
2586 are implemented, this can go.
2590 sage: from mjo.eja.eja_algebra import HadamardEJA
2594 This multiplication table can be verified by hand::
2596 sage: J = HadamardEJA(3)
2597 sage: e0,e1,e2 = J.gens()
2613 We can change the generator prefix::
2615 sage: HadamardEJA(3, prefix='r').gens()
2619 def __init__(self
, n
, **kwargs
):
2621 jordan_product
= lambda x
,y
: x
2622 inner_product
= lambda x
,y
: x
2624 def jordan_product(x
,y
):
2626 return P( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2628 def inner_product(x
,y
):
2631 # New defaults for keyword arguments. Don't orthonormalize
2632 # because our basis is already orthonormal with respect to our
2633 # inner-product. Don't check the axioms, because we know this
2634 # is a valid EJA... but do double-check if the user passes
2635 # check_axioms=True. Note: we DON'T override the "check_field"
2636 # default here, because the user can pass in a field!
2637 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2638 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2640 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2641 super().__init
__(column_basis
,
2646 self
.rank
.set_cache(n
)
2649 self
.one
.set_cache( self
.zero() )
2651 self
.one
.set_cache( sum(self
.gens()) )
2654 def _max_random_instance_size():
2656 The maximum dimension of a random HadamardEJA.
2661 def random_instance(cls
, **kwargs
):
2663 Return a random instance of this type of algebra.
2665 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2666 return cls(n
, **kwargs
)
2669 class BilinearFormEJA(ConcreteEJA
):
2671 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2672 with the half-trace inner product and jordan product ``x*y =
2673 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2674 a symmetric positive-definite "bilinear form" matrix. Its
2675 dimension is the size of `B`, and it has rank two in dimensions
2676 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2677 the identity matrix of order ``n``.
2679 We insist that the one-by-one upper-left identity block of `B` be
2680 passed in as well so that we can be passed a matrix of size zero
2681 to construct a trivial algebra.
2685 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2686 ....: JordanSpinEJA)
2690 When no bilinear form is specified, the identity matrix is used,
2691 and the resulting algebra is the Jordan spin algebra::
2693 sage: B = matrix.identity(AA,3)
2694 sage: J0 = BilinearFormEJA(B)
2695 sage: J1 = JordanSpinEJA(3)
2696 sage: J0.multiplication_table() == J0.multiplication_table()
2699 An error is raised if the matrix `B` does not correspond to a
2700 positive-definite bilinear form::
2702 sage: B = matrix.random(QQ,2,3)
2703 sage: J = BilinearFormEJA(B)
2704 Traceback (most recent call last):
2706 ValueError: bilinear form is not positive-definite
2707 sage: B = matrix.zero(QQ,3)
2708 sage: J = BilinearFormEJA(B)
2709 Traceback (most recent call last):
2711 ValueError: bilinear form is not positive-definite
2715 We can create a zero-dimensional algebra::
2717 sage: B = matrix.identity(AA,0)
2718 sage: J = BilinearFormEJA(B)
2722 We can check the multiplication condition given in the Jordan, von
2723 Neumann, and Wigner paper (and also discussed on my "On the
2724 symmetry..." paper). Note that this relies heavily on the standard
2725 choice of basis, as does anything utilizing the bilinear form
2726 matrix. We opt not to orthonormalize the basis, because if we
2727 did, we would have to normalize the `s_{i}` in a similar manner::
2729 sage: set_random_seed()
2730 sage: n = ZZ.random_element(5)
2731 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2732 sage: B11 = matrix.identity(QQ,1)
2733 sage: B22 = M.transpose()*M
2734 sage: B = block_matrix(2,2,[ [B11,0 ],
2736 sage: J = BilinearFormEJA(B, orthonormalize=False)
2737 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2738 sage: V = J.vector_space()
2739 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2740 ....: for ei in eis ]
2741 sage: actual = [ sis[i]*sis[j]
2742 ....: for i in range(n-1)
2743 ....: for j in range(n-1) ]
2744 sage: expected = [ J.one() if i == j else J.zero()
2745 ....: for i in range(n-1)
2746 ....: for j in range(n-1) ]
2747 sage: actual == expected
2751 def __init__(self
, B
, **kwargs
):
2752 # The matrix "B" is supplied by the user in most cases,
2753 # so it makes sense to check whether or not its positive-
2754 # definite unless we are specifically asked not to...
2755 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2756 if not B
.is_positive_definite():
2757 raise ValueError("bilinear form is not positive-definite")
2759 # However, all of the other data for this EJA is computed
2760 # by us in manner that guarantees the axioms are
2761 # satisfied. So, again, unless we are specifically asked to
2762 # verify things, we'll skip the rest of the checks.
2763 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2765 def inner_product(x
,y
):
2766 return (y
.T
*B
*x
)[0,0]
2768 def jordan_product(x
,y
):
2774 z0
= inner_product(y
,x
)
2775 zbar
= y0
*xbar
+ x0
*ybar
2776 return P([z0
] + zbar
.list())
2779 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2780 super(BilinearFormEJA
, self
).__init
__(column_basis
,
2785 # The rank of this algebra is two, unless we're in a
2786 # one-dimensional ambient space (because the rank is bounded
2787 # by the ambient dimension).
2788 self
.rank
.set_cache(min(n
,2))
2791 self
.one
.set_cache( self
.zero() )
2793 self
.one
.set_cache( self
.monomial(0) )
2796 def _max_random_instance_size():
2798 The maximum dimension of a random BilinearFormEJA.
2803 def random_instance(cls
, **kwargs
):
2805 Return a random instance of this algebra.
2807 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2809 B
= matrix
.identity(ZZ
, n
)
2810 return cls(B
, **kwargs
)
2812 B11
= matrix
.identity(ZZ
, 1)
2813 M
= matrix
.random(ZZ
, n
-1)
2814 I
= matrix
.identity(ZZ
, n
-1)
2816 while alpha
.is_zero():
2817 alpha
= ZZ
.random_element().abs()
2818 B22
= M
.transpose()*M
+ alpha
*I
2820 from sage
.matrix
.special
import block_matrix
2821 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2824 return cls(B
, **kwargs
)
2827 class JordanSpinEJA(BilinearFormEJA
):
2829 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2830 with the usual inner product and jordan product ``x*y =
2831 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2836 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2840 This multiplication table can be verified by hand::
2842 sage: J = JordanSpinEJA(4)
2843 sage: e0,e1,e2,e3 = J.gens()
2859 We can change the generator prefix::
2861 sage: JordanSpinEJA(2, prefix='B').gens()
2866 Ensure that we have the usual inner product on `R^n`::
2868 sage: set_random_seed()
2869 sage: J = JordanSpinEJA.random_instance()
2870 sage: x,y = J.random_elements(2)
2871 sage: actual = x.inner_product(y)
2872 sage: expected = x.to_vector().inner_product(y.to_vector())
2873 sage: actual == expected
2877 def __init__(self
, n
, **kwargs
):
2878 # This is a special case of the BilinearFormEJA with the
2879 # identity matrix as its bilinear form.
2880 B
= matrix
.identity(ZZ
, n
)
2882 # Don't orthonormalize because our basis is already
2883 # orthonormal with respect to our inner-product.
2884 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2886 # But also don't pass check_field=False here, because the user
2887 # can pass in a field!
2888 super(JordanSpinEJA
, self
).__init
__(B
, **kwargs
)
2891 def _max_random_instance_size():
2893 The maximum dimension of a random JordanSpinEJA.
2898 def random_instance(cls
, **kwargs
):
2900 Return a random instance of this type of algebra.
2902 Needed here to override the implementation for ``BilinearFormEJA``.
2904 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2905 return cls(n
, **kwargs
)
2908 class TrivialEJA(ConcreteEJA
):
2910 The trivial Euclidean Jordan algebra consisting of only a zero element.
2914 sage: from mjo.eja.eja_algebra import TrivialEJA
2918 sage: J = TrivialEJA()
2925 sage: 7*J.one()*12*J.one()
2927 sage: J.one().inner_product(J.one())
2929 sage: J.one().norm()
2931 sage: J.one().subalgebra_generated_by()
2932 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2937 def __init__(self
, **kwargs
):
2938 jordan_product
= lambda x
,y
: x
2939 inner_product
= lambda x
,y
: 0
2942 # New defaults for keyword arguments
2943 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2944 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2946 super(TrivialEJA
, self
).__init
__(basis
,
2950 # The rank is zero using my definition, namely the dimension of the
2951 # largest subalgebra generated by any element.
2952 self
.rank
.set_cache(0)
2953 self
.one
.set_cache( self
.zero() )
2956 def random_instance(cls
, **kwargs
):
2957 # We don't take a "size" argument so the superclass method is
2958 # inappropriate for us.
2959 return cls(**kwargs
)
2962 class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct
,
2963 FiniteDimensionalEJA
):
2965 The external (orthogonal) direct sum of two or more Euclidean
2966 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
2967 orthogonal direct sum of simple Euclidean Jordan algebras which is
2968 then isometric to a Cartesian product, so no generality is lost by
2969 providing only this construction.
2973 sage: from mjo.eja.eja_algebra import (random_eja,
2974 ....: CartesianProductEJA,
2976 ....: JordanSpinEJA,
2977 ....: RealSymmetricEJA)
2981 The Jordan product is inherited from our factors and implemented by
2982 our CombinatorialFreeModule Cartesian product superclass::
2984 sage: set_random_seed()
2985 sage: J1 = HadamardEJA(2)
2986 sage: J2 = RealSymmetricEJA(2)
2987 sage: J = cartesian_product([J1,J2])
2988 sage: x,y = J.random_elements(2)
2992 The ability to retrieve the original factors is implemented by our
2993 CombinatorialFreeModule Cartesian product superclass::
2995 sage: J1 = HadamardEJA(2, field=QQ)
2996 sage: J2 = JordanSpinEJA(3, field=QQ)
2997 sage: J = cartesian_product([J1,J2])
2998 sage: J.cartesian_factors()
2999 (Euclidean Jordan algebra of dimension 2 over Rational Field,
3000 Euclidean Jordan algebra of dimension 3 over Rational Field)
3002 You can provide more than two factors::
3004 sage: J1 = HadamardEJA(2)
3005 sage: J2 = JordanSpinEJA(3)
3006 sage: J3 = RealSymmetricEJA(3)
3007 sage: cartesian_product([J1,J2,J3])
3008 Euclidean Jordan algebra of dimension 2 over Algebraic Real
3009 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
3010 Real Field (+) Euclidean Jordan algebra of dimension 6 over
3011 Algebraic Real Field
3013 Rank is additive on a Cartesian product::
3015 sage: J1 = HadamardEJA(1)
3016 sage: J2 = RealSymmetricEJA(2)
3017 sage: J = cartesian_product([J1,J2])
3018 sage: J1.rank.clear_cache()
3019 sage: J2.rank.clear_cache()
3020 sage: J.rank.clear_cache()
3023 sage: J.rank() == J1.rank() + J2.rank()
3026 The same rank computation works over the rationals, with whatever
3029 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
3030 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
3031 sage: J = cartesian_product([J1,J2])
3032 sage: J1.rank.clear_cache()
3033 sage: J2.rank.clear_cache()
3034 sage: J.rank.clear_cache()
3037 sage: J.rank() == J1.rank() + J2.rank()
3040 The product algebra will be associative if and only if all of its
3041 components are associative::
3043 sage: J1 = HadamardEJA(2)
3044 sage: J1.is_associative()
3046 sage: J2 = HadamardEJA(3)
3047 sage: J2.is_associative()
3049 sage: J3 = RealSymmetricEJA(3)
3050 sage: J3.is_associative()
3052 sage: CP1 = cartesian_product([J1,J2])
3053 sage: CP1.is_associative()
3055 sage: CP2 = cartesian_product([J1,J3])
3056 sage: CP2.is_associative()
3061 All factors must share the same base field::
3063 sage: J1 = HadamardEJA(2, field=QQ)
3064 sage: J2 = RealSymmetricEJA(2)
3065 sage: CartesianProductEJA((J1,J2))
3066 Traceback (most recent call last):
3068 ValueError: all factors must share the same base field
3070 The cached unit element is the same one that would be computed::
3072 sage: set_random_seed() # long time
3073 sage: J1 = random_eja() # long time
3074 sage: J2 = random_eja() # long time
3075 sage: J = cartesian_product([J1,J2]) # long time
3076 sage: actual = J.one() # long time
3077 sage: J.one.clear_cache() # long time
3078 sage: expected = J.one() # long time
3079 sage: actual == expected # long time
3083 Element
= FiniteDimensionalEJAElement
3086 def __init__(self
, algebras
, **kwargs
):
3087 CombinatorialFreeModule_CartesianProduct
.__init
__(self
,
3090 field
= algebras
[0].base_ring()
3091 if not all( J
.base_ring() == field
for J
in algebras
):
3092 raise ValueError("all factors must share the same base field")
3094 associative
= all( m
.is_associative() for m
in algebras
)
3096 # The definition of matrix_space() and self.basis() relies
3097 # only on the stuff in the CFM_CartesianProduct class, which
3098 # we've already initialized.
3099 Js
= self
.cartesian_factors()
3101 MS
= self
.matrix_space()
3103 MS(tuple( self
.cartesian_projection(i
)(b
).to_matrix()
3104 for i
in range(m
) ))
3105 for b
in self
.basis()
3108 # Define jordan/inner products that operate on that matrix_basis.
3109 def jordan_product(x
,y
):
3111 (Js
[i
](x
[i
])*Js
[i
](y
[i
])).to_matrix() for i
in range(m
)
3114 def inner_product(x
, y
):
3116 Js
[i
](x
[i
]).inner_product(Js
[i
](y
[i
])) for i
in range(m
)
3119 # There's no need to check the field since it already came
3120 # from an EJA. Likewise the axioms are guaranteed to be
3121 # satisfied, unless the guy writing this class sucks.
3123 # If you want the basis to be orthonormalized, orthonormalize
3125 FiniteDimensionalEJA
.__init
__(self
,
3130 orthonormalize
=False,
3131 associative
=associative
,
3132 cartesian_product
=True,
3136 ones
= tuple(J
.one() for J
in algebras
)
3137 self
.one
.set_cache(self
._cartesian
_product
_of
_elements
(ones
))
3138 self
.rank
.set_cache(sum(J
.rank() for J
in algebras
))
3140 def matrix_space(self
):
3142 Return the space that our matrix basis lives in as a Cartesian
3147 sage: from mjo.eja.eja_algebra import (HadamardEJA,
3148 ....: RealSymmetricEJA)
3152 sage: J1 = HadamardEJA(1)
3153 sage: J2 = RealSymmetricEJA(2)
3154 sage: J = cartesian_product([J1,J2])
3155 sage: J.matrix_space()
3156 The Cartesian product of (Full MatrixSpace of 1 by 1 dense
3157 matrices over Algebraic Real Field, Full MatrixSpace of 2
3158 by 2 dense matrices over Algebraic Real Field)
3161 from sage
.categories
.cartesian_product
import cartesian_product
3162 return cartesian_product( [J
.matrix_space()
3163 for J
in self
.cartesian_factors()] )
3166 def cartesian_projection(self
, i
):
3170 sage: from mjo.eja.eja_algebra import (random_eja,
3171 ....: JordanSpinEJA,
3173 ....: RealSymmetricEJA,
3174 ....: ComplexHermitianEJA)
3178 The projection morphisms are Euclidean Jordan algebra
3181 sage: J1 = HadamardEJA(2)
3182 sage: J2 = RealSymmetricEJA(2)
3183 sage: J = cartesian_product([J1,J2])
3184 sage: J.cartesian_projection(0)
3185 Linear operator between finite-dimensional Euclidean Jordan
3186 algebras represented by the matrix:
3189 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3190 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3191 Algebraic Real Field
3192 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3194 sage: J.cartesian_projection(1)
3195 Linear operator between finite-dimensional Euclidean Jordan
3196 algebras represented by the matrix:
3200 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3201 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3202 Algebraic Real Field
3203 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3206 The projections work the way you'd expect on the vector
3207 representation of an element::
3209 sage: J1 = JordanSpinEJA(2)
3210 sage: J2 = ComplexHermitianEJA(2)
3211 sage: J = cartesian_product([J1,J2])
3212 sage: pi_left = J.cartesian_projection(0)
3213 sage: pi_right = J.cartesian_projection(1)
3214 sage: pi_left(J.one()).to_vector()
3216 sage: pi_right(J.one()).to_vector()
3218 sage: J.one().to_vector()
3223 The answer never changes::
3225 sage: set_random_seed()
3226 sage: J1 = random_eja()
3227 sage: J2 = random_eja()
3228 sage: J = cartesian_product([J1,J2])
3229 sage: P0 = J.cartesian_projection(0)
3230 sage: P1 = J.cartesian_projection(0)
3235 Ji
= self
.cartesian_factors()[i
]
3236 # Requires the fix on Trac 31421/31422 to work!
3237 Pi
= super().cartesian_projection(i
)
3238 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3241 def cartesian_embedding(self
, i
):
3245 sage: from mjo.eja.eja_algebra import (random_eja,
3246 ....: JordanSpinEJA,
3248 ....: RealSymmetricEJA)
3252 The embedding morphisms are Euclidean Jordan algebra
3255 sage: J1 = HadamardEJA(2)
3256 sage: J2 = RealSymmetricEJA(2)
3257 sage: J = cartesian_product([J1,J2])
3258 sage: J.cartesian_embedding(0)
3259 Linear operator between finite-dimensional Euclidean Jordan
3260 algebras represented by the matrix:
3266 Domain: Euclidean Jordan algebra of dimension 2 over
3267 Algebraic Real Field
3268 Codomain: Euclidean Jordan algebra of dimension 2 over
3269 Algebraic Real Field (+) Euclidean Jordan algebra of
3270 dimension 3 over Algebraic Real Field
3271 sage: J.cartesian_embedding(1)
3272 Linear operator between finite-dimensional Euclidean Jordan
3273 algebras represented by the matrix:
3279 Domain: Euclidean Jordan algebra of dimension 3 over
3280 Algebraic Real Field
3281 Codomain: Euclidean Jordan algebra of dimension 2 over
3282 Algebraic Real Field (+) Euclidean Jordan algebra of
3283 dimension 3 over Algebraic Real Field
3285 The embeddings work the way you'd expect on the vector
3286 representation of an element::
3288 sage: J1 = JordanSpinEJA(3)
3289 sage: J2 = RealSymmetricEJA(2)
3290 sage: J = cartesian_product([J1,J2])
3291 sage: iota_left = J.cartesian_embedding(0)
3292 sage: iota_right = J.cartesian_embedding(1)
3293 sage: iota_left(J1.zero()) == J.zero()
3295 sage: iota_right(J2.zero()) == J.zero()
3297 sage: J1.one().to_vector()
3299 sage: iota_left(J1.one()).to_vector()
3301 sage: J2.one().to_vector()
3303 sage: iota_right(J2.one()).to_vector()
3305 sage: J.one().to_vector()
3310 The answer never changes::
3312 sage: set_random_seed()
3313 sage: J1 = random_eja()
3314 sage: J2 = random_eja()
3315 sage: J = cartesian_product([J1,J2])
3316 sage: E0 = J.cartesian_embedding(0)
3317 sage: E1 = J.cartesian_embedding(0)
3321 Composing a projection with the corresponding inclusion should
3322 produce the identity map, and mismatching them should produce
3325 sage: set_random_seed()
3326 sage: J1 = random_eja()
3327 sage: J2 = random_eja()
3328 sage: J = cartesian_product([J1,J2])
3329 sage: iota_left = J.cartesian_embedding(0)
3330 sage: iota_right = J.cartesian_embedding(1)
3331 sage: pi_left = J.cartesian_projection(0)
3332 sage: pi_right = J.cartesian_projection(1)
3333 sage: pi_left*iota_left == J1.one().operator()
3335 sage: pi_right*iota_right == J2.one().operator()
3337 sage: (pi_left*iota_right).is_zero()
3339 sage: (pi_right*iota_left).is_zero()
3343 Ji
= self
.cartesian_factors()[i
]
3344 # Requires the fix on Trac 31421/31422 to work!
3345 Ei
= super().cartesian_embedding(i
)
3346 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3350 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3352 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3355 A separate class for products of algebras for which we know a
3360 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
3361 ....: RealSymmetricEJA)
3365 This gives us fast characteristic polynomial computations in
3366 product algebras, too::
3369 sage: J1 = JordanSpinEJA(2)
3370 sage: J2 = RealSymmetricEJA(3)
3371 sage: J = cartesian_product([J1,J2])
3372 sage: J.characteristic_polynomial_of().degree()
3378 def __init__(self
, algebras
, **kwargs
):
3379 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3381 self
._rational
_algebra
= None
3382 if self
.vector_space().base_field() is not QQ
:
3383 self
._rational
_algebra
= cartesian_product([
3384 r
._rational
_algebra
for r
in algebras
3388 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3390 random_eja
= ConcreteEJA
.random_instance