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.
117 sage: from mjo.eja.eja_algebra import random_eja
121 We should compute that an element subalgebra is associative even
122 if we circumvent the element method::
124 sage: set_random_seed()
125 sage: J = random_eja(field=QQ,orthonormalize=False)
126 sage: x = J.random_element()
127 sage: A = x.subalgebra_generated_by(orthonormalize=False)
128 sage: basis = tuple(b.superalgebra_element() for b in A.basis())
129 sage: J.subalgebra(basis, orthonormalize=False).is_associative()
133 Element
= FiniteDimensionalEJAElement
142 cartesian_product
=False,
147 # Keep track of whether or not the matrix basis consists of
148 # tuples, since we need special cases for them damned near
149 # everywhere. This is INDEPENDENT of whether or not the
150 # algebra is a cartesian product, since a subalgebra of a
151 # cartesian product will have a basis of tuples, but will not
152 # in general itself be a cartesian product algebra.
153 self
._matrix
_basis
_is
_cartesian
= False
156 if hasattr(basis
[0], 'cartesian_factors'):
157 self
._matrix
_basis
_is
_cartesian
= True
160 if not field
.is_subring(RR
):
161 # Note: this does return true for the real algebraic
162 # field, the rationals, and any quadratic field where
163 # we've specified a real embedding.
164 raise ValueError("scalar field is not real")
166 # If the basis given to us wasn't over the field that it's
167 # supposed to be over, fix that. Or, you know, crash.
168 if not cartesian_product
:
169 # The field for a cartesian product algebra comes from one
170 # of its factors and is the same for all factors, so
171 # there's no need to "reapply" it on product algebras.
172 if self
._matrix
_basis
_is
_cartesian
:
173 # OK since if n == 0, the basis does not consist of tuples.
174 P
= basis
[0].parent()
175 basis
= tuple( P(tuple(b_i
.change_ring(field
) for b_i
in b
))
178 basis
= tuple( b
.change_ring(field
) for b
in basis
)
182 # Check commutativity of the Jordan and inner-products.
183 # This has to be done before we build the multiplication
184 # and inner-product tables/matrices, because we take
185 # advantage of symmetry in the process.
186 if not all( jordan_product(bi
,bj
) == jordan_product(bj
,bi
)
189 raise ValueError("Jordan product is not commutative")
191 if not all( inner_product(bi
,bj
) == inner_product(bj
,bi
)
194 raise ValueError("inner-product is not commutative")
197 category
= MagmaticAlgebras(field
).FiniteDimensional()
198 category
= category
.WithBasis().Unital().Commutative()
200 if associative
is None:
201 # We should figure it out. As with check_axioms, we have to do
202 # this without the help of the _jordan_product_is_associative()
203 # method because we need to know the category before we
204 # initialize the algebra.
205 associative
= all( jordan_product(jordan_product(bi
,bj
),bk
)
207 jordan_product(bi
,jordan_product(bj
,bk
))
213 # Element subalgebras can take advantage of this.
214 category
= category
.Associative()
215 if cartesian_product
:
216 category
= category
.CartesianProducts()
218 # Call the superclass constructor so that we can use its from_vector()
219 # method to build our multiplication table.
220 CombinatorialFreeModule
.__init
__(self
,
227 # Now comes all of the hard work. We'll be constructing an
228 # ambient vector space V that our (vectorized) basis lives in,
229 # as well as a subspace W of V spanned by those (vectorized)
230 # basis elements. The W-coordinates are the coefficients that
231 # we see in things like x = 1*e1 + 2*e2.
236 degree
= len(_all2list(basis
[0]))
238 # Build an ambient space that fits our matrix basis when
239 # written out as "long vectors."
240 V
= VectorSpace(field
, degree
)
242 # The matrix that will hole the orthonormal -> unorthonormal
243 # coordinate transformation.
244 self
._deortho
_matrix
= None
247 # Save a copy of the un-orthonormalized basis for later.
248 # Convert it to ambient V (vector) coordinates while we're
249 # at it, because we'd have to do it later anyway.
250 deortho_vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
252 from mjo
.eja
.eja_utils
import gram_schmidt
253 basis
= tuple(gram_schmidt(basis
, inner_product
))
255 # Save the (possibly orthonormalized) matrix basis for
257 self
._matrix
_basis
= basis
259 # Now create the vector space for the algebra, which will have
260 # its own set of non-ambient coordinates (in terms of the
262 vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
263 W
= V
.span_of_basis( vector_basis
, check
=check_axioms
)
266 # Now "W" is the vector space of our algebra coordinates. The
267 # variables "X1", "X2",... refer to the entries of vectors in
268 # W. Thus to convert back and forth between the orthonormal
269 # coordinates and the given ones, we need to stick the original
271 U
= V
.span_of_basis( deortho_vector_basis
, check
=check_axioms
)
272 self
._deortho
_matrix
= matrix( U
.coordinate_vector(q
)
273 for q
in vector_basis
)
276 # Now we actually compute the multiplication and inner-product
277 # tables/matrices using the possibly-orthonormalized basis.
278 self
._inner
_product
_matrix
= matrix
.identity(field
, n
)
279 self
._multiplication
_table
= [ [0 for j
in range(i
+1)]
282 # Note: the Jordan and inner-products are defined in terms
283 # of the ambient basis. It's important that their arguments
284 # are in ambient coordinates as well.
287 # ortho basis w.r.t. ambient coords
291 # The jordan product returns a matrixy answer, so we
292 # have to convert it to the algebra coordinates.
293 elt
= jordan_product(q_i
, q_j
)
294 elt
= W
.coordinate_vector(V(_all2list(elt
)))
295 self
._multiplication
_table
[i
][j
] = self
.from_vector(elt
)
297 if not orthonormalize
:
298 # If we're orthonormalizing the basis with respect
299 # to an inner-product, then the inner-product
300 # matrix with respect to the resulting basis is
301 # just going to be the identity.
302 ip
= inner_product(q_i
, q_j
)
303 self
._inner
_product
_matrix
[i
,j
] = ip
304 self
._inner
_product
_matrix
[j
,i
] = ip
306 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
307 self
._inner
_product
_matrix
.set_immutable()
310 if not self
._is
_jordanian
():
311 raise ValueError("Jordan identity does not hold")
312 if not self
._inner
_product
_is
_associative
():
313 raise ValueError("inner product is not associative")
316 def _coerce_map_from_base_ring(self
):
318 Disable the map from the base ring into the algebra.
320 Performing a nonsense conversion like this automatically
321 is counterpedagogical. The fallback is to try the usual
322 element constructor, which should also fail.
326 sage: from mjo.eja.eja_algebra import random_eja
330 sage: set_random_seed()
331 sage: J = random_eja()
333 Traceback (most recent call last):
335 ValueError: not an element of this algebra
341 def product_on_basis(self
, i
, j
):
343 Returns the Jordan product of the `i` and `j`th basis elements.
345 This completely defines the Jordan product on the algebra, and
346 is used direclty by our superclass machinery to implement
351 sage: from mjo.eja.eja_algebra import random_eja
355 sage: set_random_seed()
356 sage: J = random_eja()
357 sage: n = J.dimension()
360 sage: ei_ej = J.zero()*J.zero()
362 ....: i = ZZ.random_element(n)
363 ....: j = ZZ.random_element(n)
364 ....: ei = J.gens()[i]
365 ....: ej = J.gens()[j]
366 ....: ei_ej = J.product_on_basis(i,j)
371 # We only stored the lower-triangular portion of the
372 # multiplication table.
374 return self
._multiplication
_table
[i
][j
]
376 return self
._multiplication
_table
[j
][i
]
378 def inner_product(self
, x
, y
):
380 The inner product associated with this Euclidean Jordan algebra.
382 Defaults to the trace inner product, but can be overridden by
383 subclasses if they are sure that the necessary properties are
388 sage: from mjo.eja.eja_algebra import (random_eja,
390 ....: BilinearFormEJA)
394 Our inner product is "associative," which means the following for
395 a symmetric bilinear form::
397 sage: set_random_seed()
398 sage: J = random_eja()
399 sage: x,y,z = J.random_elements(3)
400 sage: (x*y).inner_product(z) == y.inner_product(x*z)
405 Ensure that this is the usual inner product for the algebras
408 sage: set_random_seed()
409 sage: J = HadamardEJA.random_instance()
410 sage: x,y = J.random_elements(2)
411 sage: actual = x.inner_product(y)
412 sage: expected = x.to_vector().inner_product(y.to_vector())
413 sage: actual == expected
416 Ensure that this is one-half of the trace inner-product in a
417 BilinearFormEJA that isn't just the reals (when ``n`` isn't
418 one). This is in Faraut and Koranyi, and also my "On the
421 sage: set_random_seed()
422 sage: J = BilinearFormEJA.random_instance()
423 sage: n = J.dimension()
424 sage: x = J.random_element()
425 sage: y = J.random_element()
426 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
430 B
= self
._inner
_product
_matrix
431 return (B
*x
.to_vector()).inner_product(y
.to_vector())
434 def is_associative(self
):
436 Return whether or not this algebra's Jordan product is associative.
440 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
444 sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False)
445 sage: J.is_associative()
447 sage: x = sum(J.gens())
448 sage: A = x.subalgebra_generated_by(orthonormalize=False)
449 sage: A.is_associative()
453 return "Associative" in self
.category().axioms()
455 def _is_commutative(self
):
457 Whether or not this algebra's multiplication table is commutative.
459 This method should of course always return ``True``, unless
460 this algebra was constructed with ``check_axioms=False`` and
461 passed an invalid multiplication table.
463 return all( self
.product_on_basis(i
,j
) == self
.product_on_basis(i
,j
)
464 for i
in range(self
.dimension())
465 for j
in range(self
.dimension()) )
467 def _is_jordanian(self
):
469 Whether or not this algebra's multiplication table respects the
470 Jordan identity `(x^{2})(xy) = x(x^{2}y)`.
472 We only check one arrangement of `x` and `y`, so for a
473 ``True`` result to be truly true, you should also check
474 :meth:`_is_commutative`. This method should of course always
475 return ``True``, unless this algebra was constructed with
476 ``check_axioms=False`` and passed an invalid multiplication table.
478 return all( (self
.gens()[i
]**2)*(self
.gens()[i
]*self
.gens()[j
])
480 (self
.gens()[i
])*((self
.gens()[i
]**2)*self
.gens()[j
])
481 for i
in range(self
.dimension())
482 for j
in range(self
.dimension()) )
484 def _jordan_product_is_associative(self
):
486 Return whether or not this algebra's Jordan product is
487 associative; that is, whether or not `x*(y*z) = (x*y)*z`
490 This method should agree with :meth:`is_associative` unless
491 you lied about the value of the ``associative`` parameter
492 when you constructed the algebra.
496 sage: from mjo.eja.eja_algebra import (random_eja,
497 ....: RealSymmetricEJA,
498 ....: ComplexHermitianEJA,
499 ....: QuaternionHermitianEJA)
503 sage: J = RealSymmetricEJA(4, orthonormalize=False)
504 sage: J._jordan_product_is_associative()
506 sage: x = sum(J.gens())
507 sage: A = x.subalgebra_generated_by()
508 sage: A._jordan_product_is_associative()
513 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
514 sage: J._jordan_product_is_associative()
516 sage: x = sum(J.gens())
517 sage: A = x.subalgebra_generated_by(orthonormalize=False)
518 sage: A._jordan_product_is_associative()
523 sage: J = QuaternionHermitianEJA(2)
524 sage: J._jordan_product_is_associative()
526 sage: x = sum(J.gens())
527 sage: A = x.subalgebra_generated_by()
528 sage: A._jordan_product_is_associative()
533 The values we've presupplied to the constructors agree with
536 sage: set_random_seed()
537 sage: J = random_eja()
538 sage: J.is_associative() == J._jordan_product_is_associative()
544 # Used to check whether or not something is zero.
547 # I don't know of any examples that make this magnitude
548 # necessary because I don't know how to make an
549 # associative algebra when the element subalgebra
550 # construction is unreliable (as it is over RDF; we can't
551 # find the degree of an element because we can't compute
552 # the rank of a matrix). But even multiplication of floats
553 # is non-associative, so *some* epsilon is needed... let's
554 # just take the one from _inner_product_is_associative?
557 for i
in range(self
.dimension()):
558 for j
in range(self
.dimension()):
559 for k
in range(self
.dimension()):
563 diff
= (x
*y
)*z
- x
*(y
*z
)
565 if diff
.norm() > epsilon
:
570 def _inner_product_is_associative(self
):
572 Return whether or not this algebra's inner product `B` is
573 associative; that is, whether or not `B(xy,z) = B(x,yz)`.
575 This method should of course always return ``True``, unless
576 this algebra was constructed with ``check_axioms=False`` and
577 passed an invalid Jordan or inner-product.
581 # Used to check whether or not something is zero.
584 # This choice is sufficient to allow the construction of
585 # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
588 for i
in range(self
.dimension()):
589 for j
in range(self
.dimension()):
590 for k
in range(self
.dimension()):
594 diff
= (x
*y
).inner_product(z
) - x
.inner_product(y
*z
)
596 if diff
.abs() > epsilon
:
601 def _element_constructor_(self
, elt
):
603 Construct an element of this algebra from its vector or matrix
606 This gets called only after the parent element _call_ method
607 fails to find a coercion for the argument.
611 sage: from mjo.eja.eja_algebra import (random_eja,
614 ....: RealSymmetricEJA)
618 The identity in `S^n` is converted to the identity in the EJA::
620 sage: J = RealSymmetricEJA(3)
621 sage: I = matrix.identity(QQ,3)
622 sage: J(I) == J.one()
625 This skew-symmetric matrix can't be represented in the EJA::
627 sage: J = RealSymmetricEJA(3)
628 sage: A = matrix(QQ,3, lambda i,j: i-j)
630 Traceback (most recent call last):
632 ValueError: not an element of this algebra
634 Tuples work as well, provided that the matrix basis for the
635 algebra consists of them::
637 sage: J1 = HadamardEJA(3)
638 sage: J2 = RealSymmetricEJA(2)
639 sage: J = cartesian_product([J1,J2])
640 sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
645 Ensure that we can convert any element back and forth
646 faithfully between its matrix and algebra representations::
648 sage: set_random_seed()
649 sage: J = random_eja()
650 sage: x = J.random_element()
651 sage: J(x.to_matrix()) == x
654 We cannot coerce elements between algebras just because their
655 matrix representations are compatible::
657 sage: J1 = HadamardEJA(3)
658 sage: J2 = JordanSpinEJA(3)
660 Traceback (most recent call last):
662 ValueError: not an element of this algebra
664 Traceback (most recent call last):
666 ValueError: not an element of this algebra
668 msg
= "not an element of this algebra"
669 if elt
in self
.base_ring():
670 # Ensure that no base ring -> algebra coercion is performed
671 # by this method. There's some stupidity in sage that would
672 # otherwise propagate to this method; for example, sage thinks
673 # that the integer 3 belongs to the space of 2-by-2 matrices.
674 raise ValueError(msg
)
677 # Try to convert a vector into a column-matrix...
679 except (AttributeError, TypeError):
680 # and ignore failure, because we weren't really expecting
681 # a vector as an argument anyway.
684 if elt
not in self
.matrix_space():
685 raise ValueError(msg
)
687 # Thanks for nothing! Matrix spaces aren't vector spaces in
688 # Sage, so we have to figure out its matrix-basis coordinates
689 # ourselves. We use the basis space's ring instead of the
690 # element's ring because the basis space might be an algebraic
691 # closure whereas the base ring of the 3-by-3 identity matrix
692 # could be QQ instead of QQbar.
694 # And, we also have to handle Cartesian product bases (when
695 # the matrix basis consists of tuples) here. The "good news"
696 # is that we're already converting everything to long vectors,
697 # and that strategy works for tuples as well.
699 # We pass check=False because the matrix basis is "guaranteed"
700 # to be linearly independent... right? Ha ha.
702 V
= VectorSpace(self
.base_ring(), len(elt
))
703 W
= V
.span_of_basis( (V(_all2list(s
)) for s
in self
.matrix_basis()),
707 coords
= W
.coordinate_vector(V(elt
))
708 except ArithmeticError: # vector is not in free module
709 raise ValueError(msg
)
711 return self
.from_vector(coords
)
715 Return a string representation of ``self``.
719 sage: from mjo.eja.eja_algebra import JordanSpinEJA
723 Ensure that it says what we think it says::
725 sage: JordanSpinEJA(2, field=AA)
726 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
727 sage: JordanSpinEJA(3, field=RDF)
728 Euclidean Jordan algebra of dimension 3 over Real Double Field
731 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
732 return fmt
.format(self
.dimension(), self
.base_ring())
736 def characteristic_polynomial_of(self
):
738 Return the algebra's "characteristic polynomial of" function,
739 which is itself a multivariate polynomial that, when evaluated
740 at the coordinates of some algebra element, returns that
741 element's characteristic polynomial.
743 The resulting polynomial has `n+1` variables, where `n` is the
744 dimension of this algebra. The first `n` variables correspond to
745 the coordinates of an algebra element: when evaluated at the
746 coordinates of an algebra element with respect to a certain
747 basis, the result is a univariate polynomial (in the one
748 remaining variable ``t``), namely the characteristic polynomial
753 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
757 The characteristic polynomial in the spin algebra is given in
758 Alizadeh, Example 11.11::
760 sage: J = JordanSpinEJA(3)
761 sage: p = J.characteristic_polynomial_of(); p
762 X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
763 sage: xvec = J.one().to_vector()
767 By definition, the characteristic polynomial is a monic
768 degree-zero polynomial in a rank-zero algebra. Note that
769 Cayley-Hamilton is indeed satisfied since the polynomial
770 ``1`` evaluates to the identity element of the algebra on
773 sage: J = TrivialEJA()
774 sage: J.characteristic_polynomial_of()
781 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
782 a
= self
._charpoly
_coefficients
()
784 # We go to a bit of trouble here to reorder the
785 # indeterminates, so that it's easier to evaluate the
786 # characteristic polynomial at x's coordinates and get back
787 # something in terms of t, which is what we want.
788 S
= PolynomialRing(self
.base_ring(),'t')
792 S
= PolynomialRing(S
, R
.variable_names())
795 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
797 def coordinate_polynomial_ring(self
):
799 The multivariate polynomial ring in which this algebra's
800 :meth:`characteristic_polynomial_of` lives.
804 sage: from mjo.eja.eja_algebra import (HadamardEJA,
805 ....: RealSymmetricEJA)
809 sage: J = HadamardEJA(2)
810 sage: J.coordinate_polynomial_ring()
811 Multivariate Polynomial Ring in X1, X2...
812 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
813 sage: J.coordinate_polynomial_ring()
814 Multivariate Polynomial Ring in X1, X2, X3, X4, X5, X6...
817 var_names
= tuple( "X%d" % z
for z
in range(1, self
.dimension()+1) )
818 return PolynomialRing(self
.base_ring(), var_names
)
820 def inner_product(self
, x
, y
):
822 The inner product associated with this Euclidean Jordan algebra.
824 Defaults to the trace inner product, but can be overridden by
825 subclasses if they are sure that the necessary properties are
830 sage: from mjo.eja.eja_algebra import (random_eja,
832 ....: BilinearFormEJA)
836 Our inner product is "associative," which means the following for
837 a symmetric bilinear form::
839 sage: set_random_seed()
840 sage: J = random_eja()
841 sage: x,y,z = J.random_elements(3)
842 sage: (x*y).inner_product(z) == y.inner_product(x*z)
847 Ensure that this is the usual inner product for the algebras
850 sage: set_random_seed()
851 sage: J = HadamardEJA.random_instance()
852 sage: x,y = J.random_elements(2)
853 sage: actual = x.inner_product(y)
854 sage: expected = x.to_vector().inner_product(y.to_vector())
855 sage: actual == expected
858 Ensure that this is one-half of the trace inner-product in a
859 BilinearFormEJA that isn't just the reals (when ``n`` isn't
860 one). This is in Faraut and Koranyi, and also my "On the
863 sage: set_random_seed()
864 sage: J = BilinearFormEJA.random_instance()
865 sage: n = J.dimension()
866 sage: x = J.random_element()
867 sage: y = J.random_element()
868 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
871 B
= self
._inner
_product
_matrix
872 return (B
*x
.to_vector()).inner_product(y
.to_vector())
875 def is_trivial(self
):
877 Return whether or not this algebra is trivial.
879 A trivial algebra contains only the zero element.
883 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
888 sage: J = ComplexHermitianEJA(3)
894 sage: J = TrivialEJA()
899 return self
.dimension() == 0
902 def multiplication_table(self
):
904 Return a visual representation of this algebra's multiplication
905 table (on basis elements).
909 sage: from mjo.eja.eja_algebra import JordanSpinEJA
913 sage: J = JordanSpinEJA(4)
914 sage: J.multiplication_table()
915 +----++----+----+----+----+
916 | * || e0 | e1 | e2 | e3 |
917 +====++====+====+====+====+
918 | e0 || e0 | e1 | e2 | e3 |
919 +----++----+----+----+----+
920 | e1 || e1 | e0 | 0 | 0 |
921 +----++----+----+----+----+
922 | e2 || e2 | 0 | e0 | 0 |
923 +----++----+----+----+----+
924 | e3 || e3 | 0 | 0 | e0 |
925 +----++----+----+----+----+
929 # Prepend the header row.
930 M
= [["*"] + list(self
.gens())]
932 # And to each subsequent row, prepend an entry that belongs to
933 # the left-side "header column."
934 M
+= [ [self
.gens()[i
]] + [ self
.product_on_basis(i
,j
)
938 return table(M
, header_row
=True, header_column
=True, frame
=True)
941 def matrix_basis(self
):
943 Return an (often more natural) representation of this algebras
944 basis as an ordered tuple of matrices.
946 Every finite-dimensional Euclidean Jordan Algebra is a, up to
947 Jordan isomorphism, a direct sum of five simple
948 algebras---four of which comprise Hermitian matrices. And the
949 last type of algebra can of course be thought of as `n`-by-`1`
950 column matrices (ambiguusly called column vectors) to avoid
951 special cases. As a result, matrices (and column vectors) are
952 a natural representation format for Euclidean Jordan algebra
955 But, when we construct an algebra from a basis of matrices,
956 those matrix representations are lost in favor of coordinate
957 vectors *with respect to* that basis. We could eventually
958 convert back if we tried hard enough, but having the original
959 representations handy is valuable enough that we simply store
960 them and return them from this method.
962 Why implement this for non-matrix algebras? Avoiding special
963 cases for the :class:`BilinearFormEJA` pays with simplicity in
964 its own right. But mainly, we would like to be able to assume
965 that elements of a :class:`CartesianProductEJA` can be displayed
966 nicely, without having to have special classes for direct sums
967 one of whose components was a matrix algebra.
971 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
972 ....: RealSymmetricEJA)
976 sage: J = RealSymmetricEJA(2)
978 Finite family {0: e0, 1: e1, 2: e2}
979 sage: J.matrix_basis()
981 [1 0] [ 0 0.7071067811865475?] [0 0]
982 [0 0], [0.7071067811865475? 0], [0 1]
987 sage: J = JordanSpinEJA(2)
989 Finite family {0: e0, 1: e1}
990 sage: J.matrix_basis()
996 return self
._matrix
_basis
999 def matrix_space(self
):
1001 Return the matrix space in which this algebra's elements live, if
1002 we think of them as matrices (including column vectors of the
1005 "By default" this will be an `n`-by-`1` column-matrix space,
1006 except when the algebra is trivial. There it's `n`-by-`n`
1007 (where `n` is zero), to ensure that two elements of the matrix
1008 space (empty matrices) can be multiplied. For algebras of
1009 matrices, this returns the space in which their
1010 real embeddings live.
1014 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
1015 ....: JordanSpinEJA,
1016 ....: QuaternionHermitianEJA,
1021 By default, the matrix representation is just a column-matrix
1022 equivalent to the vector representation::
1024 sage: J = JordanSpinEJA(3)
1025 sage: J.matrix_space()
1026 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
1029 The matrix representation in the trivial algebra is
1030 zero-by-zero instead of the usual `n`-by-one::
1032 sage: J = TrivialEJA()
1033 sage: J.matrix_space()
1034 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
1037 The matrix space for complex/quaternion Hermitian matrix EJA
1038 is the space in which their real-embeddings live, not the
1039 original complex/quaternion matrix space::
1041 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1042 sage: J.matrix_space()
1043 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1044 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1045 sage: J.matrix_space()
1046 Full MatrixSpace of 4 by 4 dense matrices over Rational Field
1049 if self
.is_trivial():
1050 return MatrixSpace(self
.base_ring(), 0)
1052 return self
.matrix_basis()[0].parent()
1058 Return the unit element of this algebra.
1062 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1067 We can compute unit element in the Hadamard EJA::
1069 sage: J = HadamardEJA(5)
1071 e0 + e1 + e2 + e3 + e4
1073 The unit element in the Hadamard EJA is inherited in the
1074 subalgebras generated by its elements::
1076 sage: J = HadamardEJA(5)
1078 e0 + e1 + e2 + e3 + e4
1079 sage: x = sum(J.gens())
1080 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1083 sage: A.one().superalgebra_element()
1084 e0 + e1 + e2 + e3 + e4
1088 The identity element acts like the identity, regardless of
1089 whether or not we orthonormalize::
1091 sage: set_random_seed()
1092 sage: J = random_eja()
1093 sage: x = J.random_element()
1094 sage: J.one()*x == x and x*J.one() == x
1096 sage: A = x.subalgebra_generated_by()
1097 sage: y = A.random_element()
1098 sage: A.one()*y == y and y*A.one() == y
1103 sage: set_random_seed()
1104 sage: J = random_eja(field=QQ, orthonormalize=False)
1105 sage: x = J.random_element()
1106 sage: J.one()*x == x and x*J.one() == x
1108 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1109 sage: y = A.random_element()
1110 sage: A.one()*y == y and y*A.one() == y
1113 The matrix of the unit element's operator is the identity,
1114 regardless of the base field and whether or not we
1117 sage: set_random_seed()
1118 sage: J = random_eja()
1119 sage: actual = J.one().operator().matrix()
1120 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1121 sage: actual == expected
1123 sage: x = J.random_element()
1124 sage: A = x.subalgebra_generated_by()
1125 sage: actual = A.one().operator().matrix()
1126 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1127 sage: actual == expected
1132 sage: set_random_seed()
1133 sage: J = random_eja(field=QQ, orthonormalize=False)
1134 sage: actual = J.one().operator().matrix()
1135 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1136 sage: actual == expected
1138 sage: x = J.random_element()
1139 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1140 sage: actual = A.one().operator().matrix()
1141 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1142 sage: actual == expected
1145 Ensure that the cached unit element (often precomputed by
1146 hand) agrees with the computed one::
1148 sage: set_random_seed()
1149 sage: J = random_eja()
1150 sage: cached = J.one()
1151 sage: J.one.clear_cache()
1152 sage: J.one() == cached
1157 sage: set_random_seed()
1158 sage: J = random_eja(field=QQ, orthonormalize=False)
1159 sage: cached = J.one()
1160 sage: J.one.clear_cache()
1161 sage: J.one() == cached
1165 # We can brute-force compute the matrices of the operators
1166 # that correspond to the basis elements of this algebra.
1167 # If some linear combination of those basis elements is the
1168 # algebra identity, then the same linear combination of
1169 # their matrices has to be the identity matrix.
1171 # Of course, matrices aren't vectors in sage, so we have to
1172 # appeal to the "long vectors" isometry.
1173 oper_vecs
= [ _mat2vec(g
.operator().matrix()) for g
in self
.gens() ]
1175 # Now we use basic linear algebra to find the coefficients,
1176 # of the matrices-as-vectors-linear-combination, which should
1177 # work for the original algebra basis too.
1178 A
= matrix(self
.base_ring(), oper_vecs
)
1180 # We used the isometry on the left-hand side already, but we
1181 # still need to do it for the right-hand side. Recall that we
1182 # wanted something that summed to the identity matrix.
1183 b
= _mat2vec( matrix
.identity(self
.base_ring(), self
.dimension()) )
1185 # Now if there's an identity element in the algebra, this
1186 # should work. We solve on the left to avoid having to
1187 # transpose the matrix "A".
1188 return self
.from_vector(A
.solve_left(b
))
1191 def peirce_decomposition(self
, c
):
1193 The Peirce decomposition of this algebra relative to the
1196 In the future, this can be extended to a complete system of
1197 orthogonal idempotents.
1201 - ``c`` -- an idempotent of this algebra.
1205 A triple (J0, J5, J1) containing two subalgebras and one subspace
1208 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1209 corresponding to the eigenvalue zero.
1211 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1212 corresponding to the eigenvalue one-half.
1214 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1215 corresponding to the eigenvalue one.
1217 These are the only possible eigenspaces for that operator, and this
1218 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1219 orthogonal, and are subalgebras of this algebra with the appropriate
1224 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1228 The canonical example comes from the symmetric matrices, which
1229 decompose into diagonal and off-diagonal parts::
1231 sage: J = RealSymmetricEJA(3)
1232 sage: C = matrix(QQ, [ [1,0,0],
1236 sage: J0,J5,J1 = J.peirce_decomposition(c)
1238 Euclidean Jordan algebra of dimension 1...
1240 Vector space of degree 6 and dimension 2...
1242 Euclidean Jordan algebra of dimension 3...
1243 sage: J0.one().to_matrix()
1247 sage: orig_df = AA.options.display_format
1248 sage: AA.options.display_format = 'radical'
1249 sage: J.from_vector(J5.basis()[0]).to_matrix()
1253 sage: J.from_vector(J5.basis()[1]).to_matrix()
1257 sage: AA.options.display_format = orig_df
1258 sage: J1.one().to_matrix()
1265 Every algebra decomposes trivially with respect to its identity
1268 sage: set_random_seed()
1269 sage: J = random_eja()
1270 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1271 sage: J0.dimension() == 0 and J5.dimension() == 0
1273 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1276 The decomposition is into eigenspaces, and its components are
1277 therefore necessarily orthogonal. Moreover, the identity
1278 elements in the two subalgebras are the projections onto their
1279 respective subspaces of the superalgebra's identity element::
1281 sage: set_random_seed()
1282 sage: J = random_eja()
1283 sage: x = J.random_element()
1284 sage: if not J.is_trivial():
1285 ....: while x.is_nilpotent():
1286 ....: x = J.random_element()
1287 sage: c = x.subalgebra_idempotent()
1288 sage: J0,J5,J1 = J.peirce_decomposition(c)
1290 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1291 ....: w = w.superalgebra_element()
1292 ....: y = J.from_vector(y)
1293 ....: z = z.superalgebra_element()
1294 ....: ipsum += w.inner_product(y).abs()
1295 ....: ipsum += w.inner_product(z).abs()
1296 ....: ipsum += y.inner_product(z).abs()
1299 sage: J1(c) == J1.one()
1301 sage: J0(J.one() - c) == J0.one()
1305 if not c
.is_idempotent():
1306 raise ValueError("element is not idempotent: %s" % c
)
1308 # Default these to what they should be if they turn out to be
1309 # trivial, because eigenspaces_left() won't return eigenvalues
1310 # corresponding to trivial spaces (e.g. it returns only the
1311 # eigenspace corresponding to lambda=1 if you take the
1312 # decomposition relative to the identity element).
1313 trivial
= self
.subalgebra(())
1314 J0
= trivial
# eigenvalue zero
1315 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1316 J1
= trivial
# eigenvalue one
1318 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1319 if eigval
== ~
(self
.base_ring()(2)):
1322 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1323 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1329 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1334 def random_element(self
, thorough
=False):
1336 Return a random element of this algebra.
1338 Our algebra superclass method only returns a linear
1339 combination of at most two basis elements. We instead
1340 want the vector space "random element" method that
1341 returns a more diverse selection.
1345 - ``thorough`` -- (boolean; default False) whether or not we
1346 should generate irrational coefficients for the random
1347 element when our base ring is irrational; this slows the
1348 algebra operations to a crawl, but any truly random method
1352 # For a general base ring... maybe we can trust this to do the
1353 # right thing? Unlikely, but.
1354 V
= self
.vector_space()
1355 v
= V
.random_element()
1357 if self
.base_ring() is AA
:
1358 # The "random element" method of the algebraic reals is
1359 # stupid at the moment, and only returns integers between
1360 # -2 and 2, inclusive:
1362 # https://trac.sagemath.org/ticket/30875
1364 # Instead, we implement our own "random vector" method,
1365 # and then coerce that into the algebra. We use the vector
1366 # space degree here instead of the dimension because a
1367 # subalgebra could (for example) be spanned by only two
1368 # vectors, each with five coordinates. We need to
1369 # generate all five coordinates.
1371 v
*= QQbar
.random_element().real()
1373 v
*= QQ
.random_element()
1375 return self
.from_vector(V
.coordinate_vector(v
))
1377 def random_elements(self
, count
, thorough
=False):
1379 Return ``count`` random elements as a tuple.
1383 - ``thorough`` -- (boolean; default False) whether or not we
1384 should generate irrational coefficients for the random
1385 elements when our base ring is irrational; this slows the
1386 algebra operations to a crawl, but any truly random method
1391 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1395 sage: J = JordanSpinEJA(3)
1396 sage: x,y,z = J.random_elements(3)
1397 sage: all( [ x in J, y in J, z in J ])
1399 sage: len( J.random_elements(10) ) == 10
1403 return tuple( self
.random_element(thorough
)
1404 for idx
in range(count
) )
1408 def _charpoly_coefficients(self
):
1410 The `r` polynomial coefficients of the "characteristic polynomial
1415 sage: from mjo.eja.eja_algebra import random_eja
1419 The theory shows that these are all homogeneous polynomials of
1422 sage: set_random_seed()
1423 sage: J = random_eja()
1424 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1428 n
= self
.dimension()
1429 R
= self
.coordinate_polynomial_ring()
1431 F
= R
.fraction_field()
1434 # From a result in my book, these are the entries of the
1435 # basis representation of L_x.
1436 return sum( vars[k
]*self
.gens()[k
].operator().matrix()[i
,j
]
1439 L_x
= matrix(F
, n
, n
, L_x_i_j
)
1442 if self
.rank
.is_in_cache():
1444 # There's no need to pad the system with redundant
1445 # columns if we *know* they'll be redundant.
1448 # Compute an extra power in case the rank is equal to
1449 # the dimension (otherwise, we would stop at x^(r-1)).
1450 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1451 for k
in range(n
+1) ]
1452 A
= matrix
.column(F
, x_powers
[:n
])
1453 AE
= A
.extended_echelon_form()
1460 # The theory says that only the first "r" coefficients are
1461 # nonzero, and they actually live in the original polynomial
1462 # ring and not the fraction field. We negate them because in
1463 # the actual characteristic polynomial, they get moved to the
1464 # other side where x^r lives. We don't bother to trim A_rref
1465 # down to a square matrix and solve the resulting system,
1466 # because the upper-left r-by-r portion of A_rref is
1467 # guaranteed to be the identity matrix, so e.g.
1469 # A_rref.solve_right(Y)
1471 # would just be returning Y.
1472 return (-E
*b
)[:r
].change_ring(R
)
1477 Return the rank of this EJA.
1479 This is a cached method because we know the rank a priori for
1480 all of the algebras we can construct. Thus we can avoid the
1481 expensive ``_charpoly_coefficients()`` call unless we truly
1482 need to compute the whole characteristic polynomial.
1486 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1487 ....: JordanSpinEJA,
1488 ....: RealSymmetricEJA,
1489 ....: ComplexHermitianEJA,
1490 ....: QuaternionHermitianEJA,
1495 The rank of the Jordan spin algebra is always two::
1497 sage: JordanSpinEJA(2).rank()
1499 sage: JordanSpinEJA(3).rank()
1501 sage: JordanSpinEJA(4).rank()
1504 The rank of the `n`-by-`n` Hermitian real, complex, or
1505 quaternion matrices is `n`::
1507 sage: RealSymmetricEJA(4).rank()
1509 sage: ComplexHermitianEJA(3).rank()
1511 sage: QuaternionHermitianEJA(2).rank()
1516 Ensure that every EJA that we know how to construct has a
1517 positive integer rank, unless the algebra is trivial in
1518 which case its rank will be zero::
1520 sage: set_random_seed()
1521 sage: J = random_eja()
1525 sage: r > 0 or (r == 0 and J.is_trivial())
1528 Ensure that computing the rank actually works, since the ranks
1529 of all simple algebras are known and will be cached by default::
1531 sage: set_random_seed() # long time
1532 sage: J = random_eja() # long time
1533 sage: cached = J.rank() # long time
1534 sage: J.rank.clear_cache() # long time
1535 sage: J.rank() == cached # long time
1539 return len(self
._charpoly
_coefficients
())
1542 def subalgebra(self
, basis
, **kwargs
):
1544 Create a subalgebra of this algebra from the given basis.
1546 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1547 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1550 def vector_space(self
):
1552 Return the vector space that underlies this algebra.
1556 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1560 sage: J = RealSymmetricEJA(2)
1561 sage: J.vector_space()
1562 Vector space of dimension 3 over...
1565 return self
.zero().to_vector().parent().ambient_vector_space()
1569 class RationalBasisEJA(FiniteDimensionalEJA
):
1571 New class for algebras whose supplied basis elements have all rational entries.
1575 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1579 The supplied basis is orthonormalized by default::
1581 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1582 sage: J = BilinearFormEJA(B)
1583 sage: J.matrix_basis()
1600 # Abuse the check_field parameter to check that the entries of
1601 # out basis (in ambient coordinates) are in the field QQ.
1602 if not all( all(b_i
in QQ
for b_i
in b
.list()) for b
in basis
):
1603 raise TypeError("basis not rational")
1605 super().__init
__(basis
,
1609 check_field
=check_field
,
1612 self
._rational
_algebra
= None
1614 # There's no point in constructing the extra algebra if this
1615 # one is already rational.
1617 # Note: the same Jordan and inner-products work here,
1618 # because they are necessarily defined with respect to
1619 # ambient coordinates and not any particular basis.
1620 self
._rational
_algebra
= FiniteDimensionalEJA(
1625 associative
=self
.is_associative(),
1626 orthonormalize
=False,
1631 def _charpoly_coefficients(self
):
1635 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1636 ....: JordanSpinEJA)
1640 The base ring of the resulting polynomial coefficients is what
1641 it should be, and not the rationals (unless the algebra was
1642 already over the rationals)::
1644 sage: J = JordanSpinEJA(3)
1645 sage: J._charpoly_coefficients()
1646 (X1^2 - X2^2 - X3^2, -2*X1)
1647 sage: a0 = J._charpoly_coefficients()[0]
1649 Algebraic Real Field
1650 sage: a0.base_ring()
1651 Algebraic Real Field
1654 if self
._rational
_algebra
is None:
1655 # There's no need to construct *another* algebra over the
1656 # rationals if this one is already over the
1657 # rationals. Likewise, if we never orthonormalized our
1658 # basis, we might as well just use the given one.
1659 return super()._charpoly
_coefficients
()
1661 # Do the computation over the rationals. The answer will be
1662 # the same, because all we've done is a change of basis.
1663 # Then, change back from QQ to our real base ring
1664 a
= ( a_i
.change_ring(self
.base_ring())
1665 for a_i
in self
._rational
_algebra
._charpoly
_coefficients
() )
1667 if self
._deortho
_matrix
is None:
1668 # This can happen if our base ring was, say, AA and we
1669 # chose not to (or didn't need to) orthonormalize. It's
1670 # still faster to do the computations over QQ even if
1671 # the numbers in the boxes stay the same.
1674 # Otherwise, convert the coordinate variables back to the
1675 # deorthonormalized ones.
1676 R
= self
.coordinate_polynomial_ring()
1677 from sage
.modules
.free_module_element
import vector
1678 X
= vector(R
, R
.gens())
1679 BX
= self
._deortho
_matrix
*X
1681 subs_dict
= { X[i]: BX[i] for i in range(len(X)) }
1682 return tuple( a_i
.subs(subs_dict
) for a_i
in a
)
1684 class ConcreteEJA(RationalBasisEJA
):
1686 A class for the Euclidean Jordan algebras that we know by name.
1688 These are the Jordan algebras whose basis, multiplication table,
1689 rank, and so on are known a priori. More to the point, they are
1690 the Euclidean Jordan algebras for which we are able to conjure up
1691 a "random instance."
1695 sage: from mjo.eja.eja_algebra import ConcreteEJA
1699 Our basis is normalized with respect to the algebra's inner
1700 product, unless we specify otherwise::
1702 sage: set_random_seed()
1703 sage: J = ConcreteEJA.random_instance()
1704 sage: all( b.norm() == 1 for b in J.gens() )
1707 Since our basis is orthonormal with respect to the algebra's inner
1708 product, and since we know that this algebra is an EJA, any
1709 left-multiplication operator's matrix will be symmetric because
1710 natural->EJA basis representation is an isometry and within the
1711 EJA the operator is self-adjoint by the Jordan axiom::
1713 sage: set_random_seed()
1714 sage: J = ConcreteEJA.random_instance()
1715 sage: x = J.random_element()
1716 sage: x.operator().is_self_adjoint()
1721 def _max_random_instance_size():
1723 Return an integer "size" that is an upper bound on the size of
1724 this algebra when it is used in a random test
1725 case. Unfortunately, the term "size" is ambiguous -- when
1726 dealing with `R^n` under either the Hadamard or Jordan spin
1727 product, the "size" refers to the dimension `n`. When dealing
1728 with a matrix algebra (real symmetric or complex/quaternion
1729 Hermitian), it refers to the size of the matrix, which is far
1730 less than the dimension of the underlying vector space.
1732 This method must be implemented in each subclass.
1734 raise NotImplementedError
1737 def random_instance(cls
, *args
, **kwargs
):
1739 Return a random instance of this type of algebra.
1741 This method should be implemented in each subclass.
1743 from sage
.misc
.prandom
import choice
1744 eja_class
= choice(cls
.__subclasses
__())
1746 # These all bubble up to the RationalBasisEJA superclass
1747 # constructor, so any (kw)args valid there are also valid
1749 return eja_class
.random_instance(*args
, **kwargs
)
1754 def dimension_over_reals():
1756 The dimension of this matrix's base ring over the reals.
1758 The reals are dimension one over themselves, obviously; that's
1759 just `\mathbb{R}^{1}`. Likewise, the complex numbers `a + bi`
1760 have dimension two. Finally, the quaternions have dimension
1761 four over the reals.
1763 This is used to determine the size of the matrix returned from
1764 :meth:`real_embed`, among other things.
1766 raise NotImplementedError
1769 def real_embed(cls
,M
):
1771 Embed the matrix ``M`` into a space of real matrices.
1773 The matrix ``M`` can have entries in any field at the moment:
1774 the real numbers, complex numbers, or quaternions. And although
1775 they are not a field, we can probably support octonions at some
1776 point, too. This function returns a real matrix that "acts like"
1777 the original with respect to matrix multiplication; i.e.
1779 real_embed(M*N) = real_embed(M)*real_embed(N)
1782 if M
.ncols() != M
.nrows():
1783 raise ValueError("the matrix 'M' must be square")
1788 def real_unembed(cls
,M
):
1790 The inverse of :meth:`real_embed`.
1792 if M
.ncols() != M
.nrows():
1793 raise ValueError("the matrix 'M' must be square")
1794 if not ZZ(M
.nrows()).mod(cls
.dimension_over_reals()).is_zero():
1795 raise ValueError("the matrix 'M' must be a real embedding")
1799 def jordan_product(X
,Y
):
1800 return (X
*Y
+ Y
*X
)/2
1803 def trace_inner_product(cls
,X
,Y
):
1805 Compute the trace inner-product of two real-embeddings.
1809 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
1810 ....: ComplexHermitianEJA,
1811 ....: QuaternionHermitianEJA)
1815 This gives the same answer as it would if we computed the trace
1816 from the unembedded (original) matrices::
1818 sage: set_random_seed()
1819 sage: J = RealSymmetricEJA.random_instance()
1820 sage: x,y = J.random_elements(2)
1821 sage: Xe = x.to_matrix()
1822 sage: Ye = y.to_matrix()
1823 sage: X = J.real_unembed(Xe)
1824 sage: Y = J.real_unembed(Ye)
1825 sage: expected = (X*Y).trace()
1826 sage: actual = J.trace_inner_product(Xe,Ye)
1827 sage: actual == expected
1832 sage: set_random_seed()
1833 sage: J = ComplexHermitianEJA.random_instance()
1834 sage: x,y = J.random_elements(2)
1835 sage: Xe = x.to_matrix()
1836 sage: Ye = y.to_matrix()
1837 sage: X = J.real_unembed(Xe)
1838 sage: Y = J.real_unembed(Ye)
1839 sage: expected = (X*Y).trace().real()
1840 sage: actual = J.trace_inner_product(Xe,Ye)
1841 sage: actual == expected
1846 sage: set_random_seed()
1847 sage: J = QuaternionHermitianEJA.random_instance()
1848 sage: x,y = J.random_elements(2)
1849 sage: Xe = x.to_matrix()
1850 sage: Ye = y.to_matrix()
1851 sage: X = J.real_unembed(Xe)
1852 sage: Y = J.real_unembed(Ye)
1853 sage: expected = (X*Y).trace().coefficient_tuple()[0]
1854 sage: actual = J.trace_inner_product(Xe,Ye)
1855 sage: actual == expected
1859 Xu
= cls
.real_unembed(X
)
1860 Yu
= cls
.real_unembed(Y
)
1861 tr
= (Xu
*Yu
).trace()
1864 # Works in QQ, AA, RDF, et cetera.
1866 except AttributeError:
1867 # A quaternion doesn't have a real() method, but does
1868 # have coefficient_tuple() method that returns the
1869 # coefficients of 1, i, j, and k -- in that order.
1870 return tr
.coefficient_tuple()[0]
1873 class RealMatrixEJA(MatrixEJA
):
1875 def dimension_over_reals():
1879 class RealSymmetricEJA(ConcreteEJA
, RealMatrixEJA
):
1881 The rank-n simple EJA consisting of real symmetric n-by-n
1882 matrices, the usual symmetric Jordan product, and the trace inner
1883 product. It has dimension `(n^2 + n)/2` over the reals.
1887 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1891 sage: J = RealSymmetricEJA(2)
1892 sage: e0, e1, e2 = J.gens()
1900 In theory, our "field" can be any subfield of the reals::
1902 sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
1903 Euclidean Jordan algebra of dimension 3 over Real Double Field
1904 sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
1905 Euclidean Jordan algebra of dimension 3 over Real Field with
1906 53 bits of precision
1910 The dimension of this algebra is `(n^2 + n) / 2`::
1912 sage: set_random_seed()
1913 sage: n_max = RealSymmetricEJA._max_random_instance_size()
1914 sage: n = ZZ.random_element(1, n_max)
1915 sage: J = RealSymmetricEJA(n)
1916 sage: J.dimension() == (n^2 + n)/2
1919 The Jordan multiplication is what we think it is::
1921 sage: set_random_seed()
1922 sage: J = RealSymmetricEJA.random_instance()
1923 sage: x,y = J.random_elements(2)
1924 sage: actual = (x*y).to_matrix()
1925 sage: X = x.to_matrix()
1926 sage: Y = y.to_matrix()
1927 sage: expected = (X*Y + Y*X)/2
1928 sage: actual == expected
1930 sage: J(expected) == x*y
1933 We can change the generator prefix::
1935 sage: RealSymmetricEJA(3, prefix='q').gens()
1936 (q0, q1, q2, q3, q4, q5)
1938 We can construct the (trivial) algebra of rank zero::
1940 sage: RealSymmetricEJA(0)
1941 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
1945 def _denormalized_basis(cls
, n
):
1947 Return a basis for the space of real symmetric n-by-n matrices.
1951 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1955 sage: set_random_seed()
1956 sage: n = ZZ.random_element(1,5)
1957 sage: B = RealSymmetricEJA._denormalized_basis(n)
1958 sage: all( M.is_symmetric() for M in B)
1962 # The basis of symmetric matrices, as matrices, in their R^(n-by-n)
1966 for j
in range(i
+1):
1967 Eij
= matrix(ZZ
, n
, lambda k
,l
: k
==i
and l
==j
)
1971 Sij
= Eij
+ Eij
.transpose()
1977 def _max_random_instance_size():
1978 return 4 # Dimension 10
1981 def random_instance(cls
, **kwargs
):
1983 Return a random instance of this type of algebra.
1985 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
1986 return cls(n
, **kwargs
)
1988 def __init__(self
, n
, **kwargs
):
1989 # We know this is a valid EJA, but will double-check
1990 # if the user passes check_axioms=True.
1991 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
1997 super().__init
__(self
._denormalized
_basis
(n
),
1998 self
.jordan_product
,
1999 self
.trace_inner_product
,
2000 associative
=associative
,
2003 # TODO: this could be factored out somehow, but is left here
2004 # because the MatrixEJA is not presently a subclass of the
2005 # FDEJA class that defines rank() and one().
2006 self
.rank
.set_cache(n
)
2007 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2008 self
.one
.set_cache(self(idV
))
2012 class ComplexMatrixEJA(MatrixEJA
):
2013 # A manual dictionary-cache for the complex_extension() method,
2014 # since apparently @classmethods can't also be @cached_methods.
2015 _complex_extension
= {}
2018 def complex_extension(cls
,field
):
2020 The complex field that we embed/unembed, as an extension
2021 of the given ``field``.
2023 if field
in cls
._complex
_extension
:
2024 return cls
._complex
_extension
[field
]
2026 # Sage doesn't know how to adjoin the complex "i" (the root of
2027 # x^2 + 1) to a field in a general way. Here, we just enumerate
2028 # all of the cases that I have cared to support so far.
2030 # Sage doesn't know how to embed AA into QQbar, i.e. how
2031 # to adjoin sqrt(-1) to AA.
2033 elif not field
.is_exact():
2035 F
= field
.complex_field()
2037 # Works for QQ and... maybe some other fields.
2038 R
= PolynomialRing(field
, 'z')
2040 F
= field
.extension(z
**2 + 1, 'I', embedding
=CLF(-1).sqrt())
2042 cls
._complex
_extension
[field
] = F
2046 def dimension_over_reals():
2050 def real_embed(cls
,M
):
2052 Embed the n-by-n complex matrix ``M`` into the space of real
2053 matrices of size 2n-by-2n via the map the sends each entry `z = a +
2054 bi` to the block matrix ``[[a,b],[-b,a]]``.
2058 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2062 sage: F = QuadraticField(-1, 'I')
2063 sage: x1 = F(4 - 2*i)
2064 sage: x2 = F(1 + 2*i)
2067 sage: M = matrix(F,2,[[x1,x2],[x3,x4]])
2068 sage: ComplexMatrixEJA.real_embed(M)
2077 Embedding is a homomorphism (isomorphism, in fact)::
2079 sage: set_random_seed()
2080 sage: n = ZZ.random_element(3)
2081 sage: F = QuadraticField(-1, 'I')
2082 sage: X = random_matrix(F, n)
2083 sage: Y = random_matrix(F, n)
2084 sage: Xe = ComplexMatrixEJA.real_embed(X)
2085 sage: Ye = ComplexMatrixEJA.real_embed(Y)
2086 sage: XYe = ComplexMatrixEJA.real_embed(X*Y)
2091 super().real_embed(M
)
2094 # We don't need any adjoined elements...
2095 field
= M
.base_ring().base_ring()
2101 blocks
.append(matrix(field
, 2, [ [ a
, b
],
2104 return matrix
.block(field
, n
, blocks
)
2108 def real_unembed(cls
,M
):
2110 The inverse of _embed_complex_matrix().
2114 sage: from mjo.eja.eja_algebra import ComplexMatrixEJA
2118 sage: A = matrix(QQ,[ [ 1, 2, 3, 4],
2119 ....: [-2, 1, -4, 3],
2120 ....: [ 9, 10, 11, 12],
2121 ....: [-10, 9, -12, 11] ])
2122 sage: ComplexMatrixEJA.real_unembed(A)
2124 [ 10*I + 9 12*I + 11]
2128 Unembedding is the inverse of embedding::
2130 sage: set_random_seed()
2131 sage: F = QuadraticField(-1, 'I')
2132 sage: M = random_matrix(F, 3)
2133 sage: Me = ComplexMatrixEJA.real_embed(M)
2134 sage: ComplexMatrixEJA.real_unembed(Me) == M
2138 super().real_unembed(M
)
2140 d
= cls
.dimension_over_reals()
2141 F
= cls
.complex_extension(M
.base_ring())
2144 # Go top-left to bottom-right (reading order), converting every
2145 # 2-by-2 block we see to a single complex element.
2147 for k
in range(n
/d
):
2148 for j
in range(n
/d
):
2149 submat
= M
[d
*k
:d
*k
+d
,d
*j
:d
*j
+d
]
2150 if submat
[0,0] != submat
[1,1]:
2151 raise ValueError('bad on-diagonal submatrix')
2152 if submat
[0,1] != -submat
[1,0]:
2153 raise ValueError('bad off-diagonal submatrix')
2154 z
= submat
[0,0] + submat
[0,1]*i
2157 return matrix(F
, n
/d
, elements
)
2160 class ComplexHermitianEJA(ConcreteEJA
, ComplexMatrixEJA
):
2162 The rank-n simple EJA consisting of complex Hermitian n-by-n
2163 matrices over the real numbers, the usual symmetric Jordan product,
2164 and the real-part-of-trace inner product. It has dimension `n^2` over
2169 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2173 In theory, our "field" can be any subfield of the reals::
2175 sage: ComplexHermitianEJA(2, field=RDF, check_axioms=True)
2176 Euclidean Jordan algebra of dimension 4 over Real Double Field
2177 sage: ComplexHermitianEJA(2, field=RR, check_axioms=True)
2178 Euclidean Jordan algebra of dimension 4 over Real Field with
2179 53 bits of precision
2183 The dimension of this algebra is `n^2`::
2185 sage: set_random_seed()
2186 sage: n_max = ComplexHermitianEJA._max_random_instance_size()
2187 sage: n = ZZ.random_element(1, n_max)
2188 sage: J = ComplexHermitianEJA(n)
2189 sage: J.dimension() == n^2
2192 The Jordan multiplication is what we think it is::
2194 sage: set_random_seed()
2195 sage: J = ComplexHermitianEJA.random_instance()
2196 sage: x,y = J.random_elements(2)
2197 sage: actual = (x*y).to_matrix()
2198 sage: X = x.to_matrix()
2199 sage: Y = y.to_matrix()
2200 sage: expected = (X*Y + Y*X)/2
2201 sage: actual == expected
2203 sage: J(expected) == x*y
2206 We can change the generator prefix::
2208 sage: ComplexHermitianEJA(2, prefix='z').gens()
2211 We can construct the (trivial) algebra of rank zero::
2213 sage: ComplexHermitianEJA(0)
2214 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2219 def _denormalized_basis(cls
, n
):
2221 Returns a basis for the space of complex Hermitian n-by-n matrices.
2223 Why do we embed these? Basically, because all of numerical linear
2224 algebra assumes that you're working with vectors consisting of `n`
2225 entries from a field and scalars from the same field. There's no way
2226 to tell SageMath that (for example) the vectors contain complex
2227 numbers, while the scalar field is real.
2231 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2235 sage: set_random_seed()
2236 sage: n = ZZ.random_element(1,5)
2237 sage: B = ComplexHermitianEJA._denormalized_basis(n)
2238 sage: all( M.is_symmetric() for M in B)
2243 R
= PolynomialRing(field
, 'z')
2245 F
= field
.extension(z
**2 + 1, 'I')
2248 # This is like the symmetric case, but we need to be careful:
2250 # * We want conjugate-symmetry, not just symmetry.
2251 # * The diagonal will (as a result) be real.
2254 Eij
= matrix
.zero(F
,n
)
2256 for j
in range(i
+1):
2260 Sij
= cls
.real_embed(Eij
)
2263 # The second one has a minus because it's conjugated.
2264 Eij
[j
,i
] = 1 # Eij = Eij + Eij.transpose()
2265 Sij_real
= cls
.real_embed(Eij
)
2267 # Eij = I*Eij - I*Eij.transpose()
2270 Sij_imag
= cls
.real_embed(Eij
)
2276 # Since we embedded these, we can drop back to the "field" that we
2277 # started with instead of the complex extension "F".
2278 return tuple( s
.change_ring(field
) for s
in S
)
2281 def __init__(self
, n
, **kwargs
):
2282 # We know this is a valid EJA, but will double-check
2283 # if the user passes check_axioms=True.
2284 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2290 super().__init
__(self
._denormalized
_basis
(n
),
2291 self
.jordan_product
,
2292 self
.trace_inner_product
,
2293 associative
=associative
,
2295 # TODO: this could be factored out somehow, but is left here
2296 # because the MatrixEJA is not presently a subclass of the
2297 # FDEJA class that defines rank() and one().
2298 self
.rank
.set_cache(n
)
2299 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2300 self
.one
.set_cache(self(idV
))
2303 def _max_random_instance_size():
2304 return 3 # Dimension 9
2307 def random_instance(cls
, **kwargs
):
2309 Return a random instance of this type of algebra.
2311 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2312 return cls(n
, **kwargs
)
2314 class QuaternionMatrixEJA(MatrixEJA
):
2316 # A manual dictionary-cache for the quaternion_extension() method,
2317 # since apparently @classmethods can't also be @cached_methods.
2318 _quaternion_extension
= {}
2321 def quaternion_extension(cls
,field
):
2323 The quaternion field that we embed/unembed, as an extension
2324 of the given ``field``.
2326 if field
in cls
._quaternion
_extension
:
2327 return cls
._quaternion
_extension
[field
]
2329 Q
= QuaternionAlgebra(field
,-1,-1)
2331 cls
._quaternion
_extension
[field
] = Q
2335 def dimension_over_reals():
2339 def real_embed(cls
,M
):
2341 Embed the n-by-n quaternion matrix ``M`` into the space of real
2342 matrices of size 4n-by-4n by first sending each quaternion entry `z
2343 = a + bi + cj + dk` to the block-complex matrix ``[[a + bi,
2344 c+di],[-c + di, a-bi]]`, and then embedding those into a real
2349 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2353 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2354 sage: i,j,k = Q.gens()
2355 sage: x = 1 + 2*i + 3*j + 4*k
2356 sage: M = matrix(Q, 1, [[x]])
2357 sage: QuaternionMatrixEJA.real_embed(M)
2363 Embedding is a homomorphism (isomorphism, in fact)::
2365 sage: set_random_seed()
2366 sage: n = ZZ.random_element(2)
2367 sage: Q = QuaternionAlgebra(QQ,-1,-1)
2368 sage: X = random_matrix(Q, n)
2369 sage: Y = random_matrix(Q, n)
2370 sage: Xe = QuaternionMatrixEJA.real_embed(X)
2371 sage: Ye = QuaternionMatrixEJA.real_embed(Y)
2372 sage: XYe = QuaternionMatrixEJA.real_embed(X*Y)
2377 super().real_embed(M
)
2378 quaternions
= M
.base_ring()
2381 F
= QuadraticField(-1, 'I')
2386 t
= z
.coefficient_tuple()
2391 cplxM
= matrix(F
, 2, [[ a
+ b
*i
, c
+ d
*i
],
2392 [-c
+ d
*i
, a
- b
*i
]])
2393 realM
= ComplexMatrixEJA
.real_embed(cplxM
)
2394 blocks
.append(realM
)
2396 # We should have real entries by now, so use the realest field
2397 # we've got for the return value.
2398 return matrix
.block(quaternions
.base_ring(), n
, blocks
)
2403 def real_unembed(cls
,M
):
2405 The inverse of _embed_quaternion_matrix().
2409 sage: from mjo.eja.eja_algebra import QuaternionMatrixEJA
2413 sage: M = matrix(QQ, [[ 1, 2, 3, 4],
2414 ....: [-2, 1, -4, 3],
2415 ....: [-3, 4, 1, -2],
2416 ....: [-4, -3, 2, 1]])
2417 sage: QuaternionMatrixEJA.real_unembed(M)
2418 [1 + 2*i + 3*j + 4*k]
2422 Unembedding is the inverse of embedding::
2424 sage: set_random_seed()
2425 sage: Q = QuaternionAlgebra(QQ, -1, -1)
2426 sage: M = random_matrix(Q, 3)
2427 sage: Me = QuaternionMatrixEJA.real_embed(M)
2428 sage: QuaternionMatrixEJA.real_unembed(Me) == M
2432 super().real_unembed(M
)
2434 d
= cls
.dimension_over_reals()
2436 # Use the base ring of the matrix to ensure that its entries can be
2437 # multiplied by elements of the quaternion algebra.
2438 Q
= cls
.quaternion_extension(M
.base_ring())
2441 # Go top-left to bottom-right (reading order), converting every
2442 # 4-by-4 block we see to a 2-by-2 complex block, to a 1-by-1
2445 for l
in range(n
/d
):
2446 for m
in range(n
/d
):
2447 submat
= ComplexMatrixEJA
.real_unembed(
2448 M
[d
*l
:d
*l
+d
,d
*m
:d
*m
+d
] )
2449 if submat
[0,0] != submat
[1,1].conjugate():
2450 raise ValueError('bad on-diagonal submatrix')
2451 if submat
[0,1] != -submat
[1,0].conjugate():
2452 raise ValueError('bad off-diagonal submatrix')
2453 z
= submat
[0,0].real()
2454 z
+= submat
[0,0].imag()*i
2455 z
+= submat
[0,1].real()*j
2456 z
+= submat
[0,1].imag()*k
2459 return matrix(Q
, n
/d
, elements
)
2462 class QuaternionHermitianEJA(ConcreteEJA
, QuaternionMatrixEJA
):
2464 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2465 matrices, the usual symmetric Jordan product, and the
2466 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2471 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2475 In theory, our "field" can be any subfield of the reals::
2477 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2478 Euclidean Jordan algebra of dimension 6 over Real Double Field
2479 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2480 Euclidean Jordan algebra of dimension 6 over Real Field with
2481 53 bits of precision
2485 The dimension of this algebra is `2*n^2 - n`::
2487 sage: set_random_seed()
2488 sage: n_max = QuaternionHermitianEJA._max_random_instance_size()
2489 sage: n = ZZ.random_element(1, n_max)
2490 sage: J = QuaternionHermitianEJA(n)
2491 sage: J.dimension() == 2*(n^2) - n
2494 The Jordan multiplication is what we think it is::
2496 sage: set_random_seed()
2497 sage: J = QuaternionHermitianEJA.random_instance()
2498 sage: x,y = J.random_elements(2)
2499 sage: actual = (x*y).to_matrix()
2500 sage: X = x.to_matrix()
2501 sage: Y = y.to_matrix()
2502 sage: expected = (X*Y + Y*X)/2
2503 sage: actual == expected
2505 sage: J(expected) == x*y
2508 We can change the generator prefix::
2510 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2511 (a0, a1, a2, a3, a4, a5)
2513 We can construct the (trivial) algebra of rank zero::
2515 sage: QuaternionHermitianEJA(0)
2516 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2520 def _denormalized_basis(cls
, n
):
2522 Returns a basis for the space of quaternion Hermitian n-by-n matrices.
2524 Why do we embed these? Basically, because all of numerical
2525 linear algebra assumes that you're working with vectors consisting
2526 of `n` entries from a field and scalars from the same field. There's
2527 no way to tell SageMath that (for example) the vectors contain
2528 complex numbers, while the scalar field is real.
2532 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2536 sage: set_random_seed()
2537 sage: n = ZZ.random_element(1,5)
2538 sage: B = QuaternionHermitianEJA._denormalized_basis(n)
2539 sage: all( M.is_symmetric() for M in B )
2544 Q
= QuaternionAlgebra(QQ
,-1,-1)
2547 # This is like the symmetric case, but we need to be careful:
2549 # * We want conjugate-symmetry, not just symmetry.
2550 # * The diagonal will (as a result) be real.
2553 Eij
= matrix
.zero(Q
,n
)
2555 for j
in range(i
+1):
2559 Sij
= cls
.real_embed(Eij
)
2562 # The second, third, and fourth ones have a minus
2563 # because they're conjugated.
2564 # Eij = Eij + Eij.transpose()
2566 Sij_real
= cls
.real_embed(Eij
)
2568 # Eij = I*(Eij - Eij.transpose())
2571 Sij_I
= cls
.real_embed(Eij
)
2573 # Eij = J*(Eij - Eij.transpose())
2576 Sij_J
= cls
.real_embed(Eij
)
2578 # Eij = K*(Eij - Eij.transpose())
2581 Sij_K
= cls
.real_embed(Eij
)
2587 # Since we embedded these, we can drop back to the "field" that we
2588 # started with instead of the quaternion algebra "Q".
2589 return tuple( s
.change_ring(field
) for s
in S
)
2592 def __init__(self
, n
, **kwargs
):
2593 # We know this is a valid EJA, but will double-check
2594 # if the user passes check_axioms=True.
2595 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2601 super().__init
__(self
._denormalized
_basis
(n
),
2602 self
.jordan_product
,
2603 self
.trace_inner_product
,
2604 associative
=associative
,
2607 # TODO: this could be factored out somehow, but is left here
2608 # because the MatrixEJA is not presently a subclass of the
2609 # FDEJA class that defines rank() and one().
2610 self
.rank
.set_cache(n
)
2611 idV
= matrix
.identity(ZZ
, self
.dimension_over_reals()*n
)
2612 self
.one
.set_cache(self(idV
))
2616 def _max_random_instance_size():
2618 The maximum rank of a random QuaternionHermitianEJA.
2620 return 2 # Dimension 6
2623 def random_instance(cls
, **kwargs
):
2625 Return a random instance of this type of algebra.
2627 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2628 return cls(n
, **kwargs
)
2631 class HadamardEJA(ConcreteEJA
):
2633 Return the Euclidean Jordan Algebra corresponding to the set
2634 `R^n` under the Hadamard product.
2636 Note: this is nothing more than the Cartesian product of ``n``
2637 copies of the spin algebra. Once Cartesian product algebras
2638 are implemented, this can go.
2642 sage: from mjo.eja.eja_algebra import HadamardEJA
2646 This multiplication table can be verified by hand::
2648 sage: J = HadamardEJA(3)
2649 sage: e0,e1,e2 = J.gens()
2665 We can change the generator prefix::
2667 sage: HadamardEJA(3, prefix='r').gens()
2671 def __init__(self
, n
, **kwargs
):
2673 jordan_product
= lambda x
,y
: x
2674 inner_product
= lambda x
,y
: x
2676 def jordan_product(x
,y
):
2678 return P( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2680 def inner_product(x
,y
):
2683 # New defaults for keyword arguments. Don't orthonormalize
2684 # because our basis is already orthonormal with respect to our
2685 # inner-product. Don't check the axioms, because we know this
2686 # is a valid EJA... but do double-check if the user passes
2687 # check_axioms=True. Note: we DON'T override the "check_field"
2688 # default here, because the user can pass in a field!
2689 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2690 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2692 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2693 super().__init
__(column_basis
,
2698 self
.rank
.set_cache(n
)
2701 self
.one
.set_cache( self
.zero() )
2703 self
.one
.set_cache( sum(self
.gens()) )
2706 def _max_random_instance_size():
2708 The maximum dimension of a random HadamardEJA.
2713 def random_instance(cls
, **kwargs
):
2715 Return a random instance of this type of algebra.
2717 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2718 return cls(n
, **kwargs
)
2721 class BilinearFormEJA(ConcreteEJA
):
2723 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2724 with the half-trace inner product and jordan product ``x*y =
2725 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2726 a symmetric positive-definite "bilinear form" matrix. Its
2727 dimension is the size of `B`, and it has rank two in dimensions
2728 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2729 the identity matrix of order ``n``.
2731 We insist that the one-by-one upper-left identity block of `B` be
2732 passed in as well so that we can be passed a matrix of size zero
2733 to construct a trivial algebra.
2737 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2738 ....: JordanSpinEJA)
2742 When no bilinear form is specified, the identity matrix is used,
2743 and the resulting algebra is the Jordan spin algebra::
2745 sage: B = matrix.identity(AA,3)
2746 sage: J0 = BilinearFormEJA(B)
2747 sage: J1 = JordanSpinEJA(3)
2748 sage: J0.multiplication_table() == J0.multiplication_table()
2751 An error is raised if the matrix `B` does not correspond to a
2752 positive-definite bilinear form::
2754 sage: B = matrix.random(QQ,2,3)
2755 sage: J = BilinearFormEJA(B)
2756 Traceback (most recent call last):
2758 ValueError: bilinear form is not positive-definite
2759 sage: B = matrix.zero(QQ,3)
2760 sage: J = BilinearFormEJA(B)
2761 Traceback (most recent call last):
2763 ValueError: bilinear form is not positive-definite
2767 We can create a zero-dimensional algebra::
2769 sage: B = matrix.identity(AA,0)
2770 sage: J = BilinearFormEJA(B)
2774 We can check the multiplication condition given in the Jordan, von
2775 Neumann, and Wigner paper (and also discussed on my "On the
2776 symmetry..." paper). Note that this relies heavily on the standard
2777 choice of basis, as does anything utilizing the bilinear form
2778 matrix. We opt not to orthonormalize the basis, because if we
2779 did, we would have to normalize the `s_{i}` in a similar manner::
2781 sage: set_random_seed()
2782 sage: n = ZZ.random_element(5)
2783 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2784 sage: B11 = matrix.identity(QQ,1)
2785 sage: B22 = M.transpose()*M
2786 sage: B = block_matrix(2,2,[ [B11,0 ],
2788 sage: J = BilinearFormEJA(B, orthonormalize=False)
2789 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2790 sage: V = J.vector_space()
2791 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2792 ....: for ei in eis ]
2793 sage: actual = [ sis[i]*sis[j]
2794 ....: for i in range(n-1)
2795 ....: for j in range(n-1) ]
2796 sage: expected = [ J.one() if i == j else J.zero()
2797 ....: for i in range(n-1)
2798 ....: for j in range(n-1) ]
2799 sage: actual == expected
2803 def __init__(self
, B
, **kwargs
):
2804 # The matrix "B" is supplied by the user in most cases,
2805 # so it makes sense to check whether or not its positive-
2806 # definite unless we are specifically asked not to...
2807 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2808 if not B
.is_positive_definite():
2809 raise ValueError("bilinear form is not positive-definite")
2811 # However, all of the other data for this EJA is computed
2812 # by us in manner that guarantees the axioms are
2813 # satisfied. So, again, unless we are specifically asked to
2814 # verify things, we'll skip the rest of the checks.
2815 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2817 def inner_product(x
,y
):
2818 return (y
.T
*B
*x
)[0,0]
2820 def jordan_product(x
,y
):
2826 z0
= inner_product(y
,x
)
2827 zbar
= y0
*xbar
+ x0
*ybar
2828 return P([z0
] + zbar
.list())
2831 column_basis
= tuple( b
.column() for b
in FreeModule(ZZ
, n
).basis() )
2833 # TODO: I haven't actually checked this, but it seems legit.
2838 super().__init
__(column_basis
,
2841 associative
=associative
,
2844 # The rank of this algebra is two, unless we're in a
2845 # one-dimensional ambient space (because the rank is bounded
2846 # by the ambient dimension).
2847 self
.rank
.set_cache(min(n
,2))
2850 self
.one
.set_cache( self
.zero() )
2852 self
.one
.set_cache( self
.monomial(0) )
2855 def _max_random_instance_size():
2857 The maximum dimension of a random BilinearFormEJA.
2862 def random_instance(cls
, **kwargs
):
2864 Return a random instance of this algebra.
2866 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2868 B
= matrix
.identity(ZZ
, n
)
2869 return cls(B
, **kwargs
)
2871 B11
= matrix
.identity(ZZ
, 1)
2872 M
= matrix
.random(ZZ
, n
-1)
2873 I
= matrix
.identity(ZZ
, n
-1)
2875 while alpha
.is_zero():
2876 alpha
= ZZ
.random_element().abs()
2877 B22
= M
.transpose()*M
+ alpha
*I
2879 from sage
.matrix
.special
import block_matrix
2880 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2883 return cls(B
, **kwargs
)
2886 class JordanSpinEJA(BilinearFormEJA
):
2888 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2889 with the usual inner product and jordan product ``x*y =
2890 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2895 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2899 This multiplication table can be verified by hand::
2901 sage: J = JordanSpinEJA(4)
2902 sage: e0,e1,e2,e3 = J.gens()
2918 We can change the generator prefix::
2920 sage: JordanSpinEJA(2, prefix='B').gens()
2925 Ensure that we have the usual inner product on `R^n`::
2927 sage: set_random_seed()
2928 sage: J = JordanSpinEJA.random_instance()
2929 sage: x,y = J.random_elements(2)
2930 sage: actual = x.inner_product(y)
2931 sage: expected = x.to_vector().inner_product(y.to_vector())
2932 sage: actual == expected
2936 def __init__(self
, n
, **kwargs
):
2937 # This is a special case of the BilinearFormEJA with the
2938 # identity matrix as its bilinear form.
2939 B
= matrix
.identity(ZZ
, n
)
2941 # Don't orthonormalize because our basis is already
2942 # orthonormal with respect to our inner-product.
2943 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2945 # But also don't pass check_field=False here, because the user
2946 # can pass in a field!
2947 super().__init
__(B
, **kwargs
)
2950 def _max_random_instance_size():
2952 The maximum dimension of a random JordanSpinEJA.
2957 def random_instance(cls
, **kwargs
):
2959 Return a random instance of this type of algebra.
2961 Needed here to override the implementation for ``BilinearFormEJA``.
2963 n
= ZZ
.random_element(cls
._max
_random
_instance
_size
() + 1)
2964 return cls(n
, **kwargs
)
2967 class TrivialEJA(ConcreteEJA
):
2969 The trivial Euclidean Jordan algebra consisting of only a zero element.
2973 sage: from mjo.eja.eja_algebra import TrivialEJA
2977 sage: J = TrivialEJA()
2984 sage: 7*J.one()*12*J.one()
2986 sage: J.one().inner_product(J.one())
2988 sage: J.one().norm()
2990 sage: J.one().subalgebra_generated_by()
2991 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2996 def __init__(self
, **kwargs
):
2997 jordan_product
= lambda x
,y
: x
2998 inner_product
= lambda x
,y
: 0
3001 # New defaults for keyword arguments
3002 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
3003 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
3005 super().__init
__(basis
,
3011 # The rank is zero using my definition, namely the dimension of the
3012 # largest subalgebra generated by any element.
3013 self
.rank
.set_cache(0)
3014 self
.one
.set_cache( self
.zero() )
3017 def random_instance(cls
, **kwargs
):
3018 # We don't take a "size" argument so the superclass method is
3019 # inappropriate for us.
3020 return cls(**kwargs
)
3023 class CartesianProductEJA(CombinatorialFreeModule_CartesianProduct
,
3024 FiniteDimensionalEJA
):
3026 The external (orthogonal) direct sum of two or more Euclidean
3027 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
3028 orthogonal direct sum of simple Euclidean Jordan algebras which is
3029 then isometric to a Cartesian product, so no generality is lost by
3030 providing only this construction.
3034 sage: from mjo.eja.eja_algebra import (random_eja,
3035 ....: CartesianProductEJA,
3037 ....: JordanSpinEJA,
3038 ....: RealSymmetricEJA)
3042 The Jordan product is inherited from our factors and implemented by
3043 our CombinatorialFreeModule Cartesian product superclass::
3045 sage: set_random_seed()
3046 sage: J1 = HadamardEJA(2)
3047 sage: J2 = RealSymmetricEJA(2)
3048 sage: J = cartesian_product([J1,J2])
3049 sage: x,y = J.random_elements(2)
3053 The ability to retrieve the original factors is implemented by our
3054 CombinatorialFreeModule Cartesian product superclass::
3056 sage: J1 = HadamardEJA(2, field=QQ)
3057 sage: J2 = JordanSpinEJA(3, field=QQ)
3058 sage: J = cartesian_product([J1,J2])
3059 sage: J.cartesian_factors()
3060 (Euclidean Jordan algebra of dimension 2 over Rational Field,
3061 Euclidean Jordan algebra of dimension 3 over Rational Field)
3063 You can provide more than two factors::
3065 sage: J1 = HadamardEJA(2)
3066 sage: J2 = JordanSpinEJA(3)
3067 sage: J3 = RealSymmetricEJA(3)
3068 sage: cartesian_product([J1,J2,J3])
3069 Euclidean Jordan algebra of dimension 2 over Algebraic Real
3070 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
3071 Real Field (+) Euclidean Jordan algebra of dimension 6 over
3072 Algebraic Real Field
3074 Rank is additive on a Cartesian product::
3076 sage: J1 = HadamardEJA(1)
3077 sage: J2 = RealSymmetricEJA(2)
3078 sage: J = cartesian_product([J1,J2])
3079 sage: J1.rank.clear_cache()
3080 sage: J2.rank.clear_cache()
3081 sage: J.rank.clear_cache()
3084 sage: J.rank() == J1.rank() + J2.rank()
3087 The same rank computation works over the rationals, with whatever
3090 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
3091 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
3092 sage: J = cartesian_product([J1,J2])
3093 sage: J1.rank.clear_cache()
3094 sage: J2.rank.clear_cache()
3095 sage: J.rank.clear_cache()
3098 sage: J.rank() == J1.rank() + J2.rank()
3101 The product algebra will be associative if and only if all of its
3102 components are associative::
3104 sage: J1 = HadamardEJA(2)
3105 sage: J1.is_associative()
3107 sage: J2 = HadamardEJA(3)
3108 sage: J2.is_associative()
3110 sage: J3 = RealSymmetricEJA(3)
3111 sage: J3.is_associative()
3113 sage: CP1 = cartesian_product([J1,J2])
3114 sage: CP1.is_associative()
3116 sage: CP2 = cartesian_product([J1,J3])
3117 sage: CP2.is_associative()
3122 All factors must share the same base field::
3124 sage: J1 = HadamardEJA(2, field=QQ)
3125 sage: J2 = RealSymmetricEJA(2)
3126 sage: CartesianProductEJA((J1,J2))
3127 Traceback (most recent call last):
3129 ValueError: all factors must share the same base field
3131 The cached unit element is the same one that would be computed::
3133 sage: set_random_seed() # long time
3134 sage: J1 = random_eja() # long time
3135 sage: J2 = random_eja() # long time
3136 sage: J = cartesian_product([J1,J2]) # long time
3137 sage: actual = J.one() # long time
3138 sage: J.one.clear_cache() # long time
3139 sage: expected = J.one() # long time
3140 sage: actual == expected # long time
3144 Element
= FiniteDimensionalEJAElement
3147 def __init__(self
, algebras
, **kwargs
):
3148 CombinatorialFreeModule_CartesianProduct
.__init
__(self
,
3151 field
= algebras
[0].base_ring()
3152 if not all( J
.base_ring() == field
for J
in algebras
):
3153 raise ValueError("all factors must share the same base field")
3155 associative
= all( m
.is_associative() for m
in algebras
)
3157 # The definition of matrix_space() and self.basis() relies
3158 # only on the stuff in the CFM_CartesianProduct class, which
3159 # we've already initialized.
3160 Js
= self
.cartesian_factors()
3162 MS
= self
.matrix_space()
3164 MS(tuple( self
.cartesian_projection(i
)(b
).to_matrix()
3165 for i
in range(m
) ))
3166 for b
in self
.basis()
3169 # Define jordan/inner products that operate on that matrix_basis.
3170 def jordan_product(x
,y
):
3172 (Js
[i
](x
[i
])*Js
[i
](y
[i
])).to_matrix() for i
in range(m
)
3175 def inner_product(x
, y
):
3177 Js
[i
](x
[i
]).inner_product(Js
[i
](y
[i
])) for i
in range(m
)
3180 # There's no need to check the field since it already came
3181 # from an EJA. Likewise the axioms are guaranteed to be
3182 # satisfied, unless the guy writing this class sucks.
3184 # If you want the basis to be orthonormalized, orthonormalize
3186 FiniteDimensionalEJA
.__init
__(self
,
3191 orthonormalize
=False,
3192 associative
=associative
,
3193 cartesian_product
=True,
3197 ones
= tuple(J
.one() for J
in algebras
)
3198 self
.one
.set_cache(self
._cartesian
_product
_of
_elements
(ones
))
3199 self
.rank
.set_cache(sum(J
.rank() for J
in algebras
))
3201 def matrix_space(self
):
3203 Return the space that our matrix basis lives in as a Cartesian
3208 sage: from mjo.eja.eja_algebra import (HadamardEJA,
3209 ....: RealSymmetricEJA)
3213 sage: J1 = HadamardEJA(1)
3214 sage: J2 = RealSymmetricEJA(2)
3215 sage: J = cartesian_product([J1,J2])
3216 sage: J.matrix_space()
3217 The Cartesian product of (Full MatrixSpace of 1 by 1 dense
3218 matrices over Algebraic Real Field, Full MatrixSpace of 2
3219 by 2 dense matrices over Algebraic Real Field)
3222 from sage
.categories
.cartesian_product
import cartesian_product
3223 return cartesian_product( [J
.matrix_space()
3224 for J
in self
.cartesian_factors()] )
3227 def cartesian_projection(self
, i
):
3231 sage: from mjo.eja.eja_algebra import (random_eja,
3232 ....: JordanSpinEJA,
3234 ....: RealSymmetricEJA,
3235 ....: ComplexHermitianEJA)
3239 The projection morphisms are Euclidean Jordan algebra
3242 sage: J1 = HadamardEJA(2)
3243 sage: J2 = RealSymmetricEJA(2)
3244 sage: J = cartesian_product([J1,J2])
3245 sage: J.cartesian_projection(0)
3246 Linear operator between finite-dimensional Euclidean Jordan
3247 algebras represented by the matrix:
3250 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3251 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3252 Algebraic Real Field
3253 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3255 sage: J.cartesian_projection(1)
3256 Linear operator between finite-dimensional Euclidean Jordan
3257 algebras represented by the matrix:
3261 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3262 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3263 Algebraic Real Field
3264 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3267 The projections work the way you'd expect on the vector
3268 representation of an element::
3270 sage: J1 = JordanSpinEJA(2)
3271 sage: J2 = ComplexHermitianEJA(2)
3272 sage: J = cartesian_product([J1,J2])
3273 sage: pi_left = J.cartesian_projection(0)
3274 sage: pi_right = J.cartesian_projection(1)
3275 sage: pi_left(J.one()).to_vector()
3277 sage: pi_right(J.one()).to_vector()
3279 sage: J.one().to_vector()
3284 The answer never changes::
3286 sage: set_random_seed()
3287 sage: J1 = random_eja()
3288 sage: J2 = random_eja()
3289 sage: J = cartesian_product([J1,J2])
3290 sage: P0 = J.cartesian_projection(0)
3291 sage: P1 = J.cartesian_projection(0)
3296 Ji
= self
.cartesian_factors()[i
]
3297 # Requires the fix on Trac 31421/31422 to work!
3298 Pi
= super().cartesian_projection(i
)
3299 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3302 def cartesian_embedding(self
, i
):
3306 sage: from mjo.eja.eja_algebra import (random_eja,
3307 ....: JordanSpinEJA,
3309 ....: RealSymmetricEJA)
3313 The embedding morphisms are Euclidean Jordan algebra
3316 sage: J1 = HadamardEJA(2)
3317 sage: J2 = RealSymmetricEJA(2)
3318 sage: J = cartesian_product([J1,J2])
3319 sage: J.cartesian_embedding(0)
3320 Linear operator between finite-dimensional Euclidean Jordan
3321 algebras represented by the matrix:
3327 Domain: Euclidean Jordan algebra of dimension 2 over
3328 Algebraic Real Field
3329 Codomain: Euclidean Jordan algebra of dimension 2 over
3330 Algebraic Real Field (+) Euclidean Jordan algebra of
3331 dimension 3 over Algebraic Real Field
3332 sage: J.cartesian_embedding(1)
3333 Linear operator between finite-dimensional Euclidean Jordan
3334 algebras represented by the matrix:
3340 Domain: Euclidean Jordan algebra of dimension 3 over
3341 Algebraic Real Field
3342 Codomain: Euclidean Jordan algebra of dimension 2 over
3343 Algebraic Real Field (+) Euclidean Jordan algebra of
3344 dimension 3 over Algebraic Real Field
3346 The embeddings work the way you'd expect on the vector
3347 representation of an element::
3349 sage: J1 = JordanSpinEJA(3)
3350 sage: J2 = RealSymmetricEJA(2)
3351 sage: J = cartesian_product([J1,J2])
3352 sage: iota_left = J.cartesian_embedding(0)
3353 sage: iota_right = J.cartesian_embedding(1)
3354 sage: iota_left(J1.zero()) == J.zero()
3356 sage: iota_right(J2.zero()) == J.zero()
3358 sage: J1.one().to_vector()
3360 sage: iota_left(J1.one()).to_vector()
3362 sage: J2.one().to_vector()
3364 sage: iota_right(J2.one()).to_vector()
3366 sage: J.one().to_vector()
3371 The answer never changes::
3373 sage: set_random_seed()
3374 sage: J1 = random_eja()
3375 sage: J2 = random_eja()
3376 sage: J = cartesian_product([J1,J2])
3377 sage: E0 = J.cartesian_embedding(0)
3378 sage: E1 = J.cartesian_embedding(0)
3382 Composing a projection with the corresponding inclusion should
3383 produce the identity map, and mismatching them should produce
3386 sage: set_random_seed()
3387 sage: J1 = random_eja()
3388 sage: J2 = random_eja()
3389 sage: J = cartesian_product([J1,J2])
3390 sage: iota_left = J.cartesian_embedding(0)
3391 sage: iota_right = J.cartesian_embedding(1)
3392 sage: pi_left = J.cartesian_projection(0)
3393 sage: pi_right = J.cartesian_projection(1)
3394 sage: pi_left*iota_left == J1.one().operator()
3396 sage: pi_right*iota_right == J2.one().operator()
3398 sage: (pi_left*iota_right).is_zero()
3400 sage: (pi_right*iota_left).is_zero()
3404 Ji
= self
.cartesian_factors()[i
]
3405 # Requires the fix on Trac 31421/31422 to work!
3406 Ei
= super().cartesian_embedding(i
)
3407 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3411 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3413 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3416 A separate class for products of algebras for which we know a
3421 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
3422 ....: RealSymmetricEJA)
3426 This gives us fast characteristic polynomial computations in
3427 product algebras, too::
3430 sage: J1 = JordanSpinEJA(2)
3431 sage: J2 = RealSymmetricEJA(3)
3432 sage: J = cartesian_product([J1,J2])
3433 sage: J.characteristic_polynomial_of().degree()
3439 def __init__(self
, algebras
, **kwargs
):
3440 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3442 self
._rational
_algebra
= None
3443 if self
.vector_space().base_field() is not QQ
:
3444 self
._rational
_algebra
= cartesian_product([
3445 r
._rational
_algebra
for r
in algebras
3449 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3451 random_eja
= ConcreteEJA
.random_instance