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`
35 * :class:`OctonionHermitianEJA`
37 In addition to these, we provide two other example constructions,
39 * :class:`JordanSpinEJA`
40 * :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. The Albert EJA is simply a special
47 case of the :class:`OctonionHermitianEJA` where the matrices are
48 three-by-three and the resulting space has dimension 27. And
49 last/least, the trivial EJA is exactly what you think it is; it could
50 also be obtained by constructing a dimension-zero instance of any of
51 the other algebras. Cartesian products of these are also supported
52 using the usual ``cartesian_product()`` function; as a result, we
53 support (up to isomorphism) all Euclidean Jordan algebras.
55 At a minimum, the following are required to construct a Euclidean
58 * A basis of matrices, column vectors, or MatrixAlgebra elements
59 * A Jordan product defined on the basis
60 * Its inner product defined on the basis
62 The real numbers form a Euclidean Jordan algebra when both the Jordan
63 and inner products are the usual multiplication. We use this as our
64 example, and demonstrate a few ways to construct an EJA.
66 First, we can use one-by-one SageMath matrices with algebraic real
67 entries to represent real numbers. We define the Jordan and inner
68 products to be essentially real-number multiplication, with the only
69 difference being that the Jordan product again returns a one-by-one
70 matrix, whereas the inner product must return a scalar. Our basis for
71 the one-by-one matrices is of course the set consisting of a single
72 matrix with its sole entry non-zero::
74 sage: from mjo.eja.eja_algebra import FiniteDimensionalEJA
75 sage: jp = lambda X,Y: X*Y
76 sage: ip = lambda X,Y: X[0,0]*Y[0,0]
77 sage: b1 = matrix(AA, [[1]])
78 sage: J1 = FiniteDimensionalEJA((b1,), jp, ip)
80 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
82 In fact, any positive scalar multiple of that inner-product would work::
84 sage: ip2 = lambda X,Y: 16*ip(X,Y)
85 sage: J2 = FiniteDimensionalEJA((b1,), jp, ip2)
87 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
89 But beware that your basis will be orthonormalized _with respect to the
90 given inner-product_ unless you pass ``orthonormalize=False`` to the
91 constructor. For example::
93 sage: J3 = FiniteDimensionalEJA((b1,), jp, ip2, orthonormalize=False)
95 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
97 To see the difference, you can take the first and only basis element
98 of the resulting algebra, and ask for it to be converted back into
101 sage: J1.basis()[0].to_matrix()
103 sage: J2.basis()[0].to_matrix()
105 sage: J3.basis()[0].to_matrix()
108 Since square roots are used in that process, the default scalar field
109 that we use is the field of algebraic real numbers, ``AA``. You can
110 also Use rational numbers, but only if you either pass
111 ``orthonormalize=False`` or know that orthonormalizing your basis
112 won't stray beyond the rational numbers. The example above would
113 have worked only because ``sqrt(16) == 4`` is rational.
115 Another option for your basis is to use elemebts of a
116 :class:`MatrixAlgebra`::
118 sage: from mjo.matrix_algebra import MatrixAlgebra
119 sage: A = MatrixAlgebra(1,AA,AA)
120 sage: J4 = FiniteDimensionalEJA(A.gens(), jp, ip)
122 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
123 sage: J4.basis()[0].to_matrix()
128 An easier way to view the entire EJA basis in its original (but
129 perhaps orthonormalized) matrix form is to use the ``matrix_basis``
132 sage: J4.matrix_basis()
137 In particular, a :class:`MatrixAlgebra` is needed to work around the
138 fact that matrices in SageMath must have entries in the same
139 (commutative and associative) ring as its scalars. There are many
140 Euclidean Jordan algebras whose elements are matrices that violate
141 those assumptions. The complex, quaternion, and octonion Hermitian
142 matrices all have entries in a ring (the complex numbers, quaternions,
143 or octonions...) that differs from the algebra's scalar ring (the real
144 numbers). Quaternions are also non-commutative; the octonions are
145 neither commutative nor associative.
149 sage: from mjo.eja.eja_algebra import random_eja
154 Euclidean Jordan algebra of dimension...
157 from sage
.algebras
.quatalg
.quaternion_algebra
import QuaternionAlgebra
158 from sage
.categories
.magmatic_algebras
import MagmaticAlgebras
159 from sage
.categories
.sets_cat
import cartesian_product
160 from sage
.combinat
.free_module
import CombinatorialFreeModule
161 from sage
.matrix
.constructor
import matrix
162 from sage
.matrix
.matrix_space
import MatrixSpace
163 from sage
.misc
.cachefunc
import cached_method
164 from sage
.misc
.table
import table
165 from sage
.modules
.free_module
import FreeModule
, VectorSpace
166 from sage
.rings
.all
import (ZZ
, QQ
, AA
, QQbar
, RR
, RLF
, CLF
,
169 from mjo
.eja
.eja_element
import (CartesianProductEJAElement
,
170 FiniteDimensionalEJAElement
)
171 from mjo
.eja
.eja_operator
import FiniteDimensionalEJAOperator
172 from mjo
.eja
.eja_utils
import _all2list
174 def EuclideanJordanAlgebras(field
):
176 The category of Euclidean Jordan algebras over ``field``, which
177 must be a subfield of the real numbers. For now this is just a
178 convenient wrapper around all of the other category axioms that
181 category
= MagmaticAlgebras(field
).FiniteDimensional()
182 category
= category
.WithBasis().Unital().Commutative()
185 class FiniteDimensionalEJA(CombinatorialFreeModule
):
187 A finite-dimensional Euclidean Jordan algebra.
191 - ``basis`` -- a tuple; a tuple of basis elements in "matrix
192 form," which must be the same form as the arguments to
193 ``jordan_product`` and ``inner_product``. In reality, "matrix
194 form" can be either vectors, matrices, or a Cartesian product
195 (ordered tuple) of vectors or matrices. All of these would
196 ideally be vector spaces in sage with no special-casing
197 needed; but in reality we turn vectors into column-matrices
198 and Cartesian products `(a,b)` into column matrices
199 `(a,b)^{T}` after converting `a` and `b` themselves.
201 - ``jordan_product`` -- a function; afunction of two ``basis``
202 elements (in matrix form) that returns their jordan product,
203 also in matrix form; this will be applied to ``basis`` to
204 compute a multiplication table for the algebra.
206 - ``inner_product`` -- a function; a function of two ``basis``
207 elements (in matrix form) that returns their inner
208 product. This will be applied to ``basis`` to compute an
209 inner-product table (basically a matrix) for this algebra.
211 - ``matrix_space`` -- the space that your matrix basis lives in,
212 or ``None`` (the default). So long as your basis does not have
213 length zero you can omit this. But in trivial algebras, it is
216 - ``field`` -- a subfield of the reals (default: ``AA``); the scalar
217 field for the algebra.
219 - ``orthonormalize`` -- boolean (default: ``True``); whether or
220 not to orthonormalize the basis. Doing so is expensive and
221 generally rules out using the rationals as your ``field``, but
222 is required for spectral decompositions.
226 sage: from mjo.eja.eja_algebra import random_eja
230 We should compute that an element subalgebra is associative even
231 if we circumvent the element method::
233 sage: J = random_eja(field=QQ,orthonormalize=False)
234 sage: x = J.random_element()
235 sage: A = x.subalgebra_generated_by(orthonormalize=False)
236 sage: basis = tuple(b.superalgebra_element() for b in A.basis())
237 sage: J.subalgebra(basis, orthonormalize=False).is_associative()
240 Element
= FiniteDimensionalEJAElement
243 def _check_input_field(field
):
244 if not field
.is_subring(RR
):
245 # Note: this does return true for the real algebraic
246 # field, the rationals, and any quadratic field where
247 # we've specified a real embedding.
248 raise ValueError("scalar field is not real")
251 def _check_input_axioms(basis
, jordan_product
, inner_product
):
252 if not all( jordan_product(bi
,bj
) == jordan_product(bj
,bi
)
255 raise ValueError("Jordan product is not commutative")
257 if not all( inner_product(bi
,bj
) == inner_product(bj
,bi
)
260 raise ValueError("inner-product is not commutative")
277 self
._check
_input
_field
(field
)
280 # Check commutativity of the Jordan and inner-products.
281 # This has to be done before we build the multiplication
282 # and inner-product tables/matrices, because we take
283 # advantage of symmetry in the process.
284 self
._check
_input
_axioms
(basis
, jordan_product
, inner_product
)
287 # All zero- and one-dimensional algebras are just the real
288 # numbers with (some positive multiples of) the usual
289 # multiplication as its Jordan and inner-product.
291 if associative
is None:
292 # We should figure it out. As with check_axioms, we have to do
293 # this without the help of the _jordan_product_is_associative()
294 # method because we need to know the category before we
295 # initialize the algebra.
296 associative
= all( jordan_product(jordan_product(bi
,bj
),bk
)
298 jordan_product(bi
,jordan_product(bj
,bk
))
303 category
= EuclideanJordanAlgebras(field
)
306 # Element subalgebras can take advantage of this.
307 category
= category
.Associative()
309 # Call the superclass constructor so that we can use its from_vector()
310 # method to build our multiplication table.
311 CombinatorialFreeModule
.__init
__(self
,
318 # Now comes all of the hard work. We'll be constructing an
319 # ambient vector space V that our (vectorized) basis lives in,
320 # as well as a subspace W of V spanned by those (vectorized)
321 # basis elements. The W-coordinates are the coefficients that
322 # we see in things like x = 1*b1 + 2*b2.
326 degree
= len(_all2list(basis
[0]))
328 # Build an ambient space that fits our matrix basis when
329 # written out as "long vectors."
330 V
= VectorSpace(field
, degree
)
332 # The matrix that will hold the orthonormal -> unorthonormal
333 # coordinate transformation. Default to an identity matrix of
334 # the appropriate size to avoid special cases for None
336 self
._deortho
_matrix
= matrix
.identity(field
,n
)
339 # Save a copy of the un-orthonormalized basis for later.
340 # Convert it to ambient V (vector) coordinates while we're
341 # at it, because we'd have to do it later anyway.
342 deortho_vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
344 from mjo
.eja
.eja_utils
import gram_schmidt
345 basis
= tuple(gram_schmidt(basis
, inner_product
))
347 # Save the (possibly orthonormalized) matrix basis for
348 # later, as well as the space that its elements live in.
349 # In most cases we can deduce the matrix space, but when
350 # n == 0 (that is, there are no basis elements) we cannot.
351 self
._matrix
_basis
= basis
352 if matrix_space
is None:
353 self
._matrix
_space
= self
._matrix
_basis
[0].parent()
355 self
._matrix
_space
= matrix_space
357 # Now create the vector space for the algebra, which will have
358 # its own set of non-ambient coordinates (in terms of the
360 vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
362 # Save the span of our matrix basis (when written out as long
363 # vectors) because otherwise we'll have to reconstruct it
364 # every time we want to coerce a matrix into the algebra.
365 self
._matrix
_span
= V
.span_of_basis( vector_basis
, check
=check_axioms
)
368 # Now "self._matrix_span" is the vector space of our
369 # algebra coordinates. The variables "X0", "X1",... refer
370 # to the entries of vectors in self._matrix_span. Thus to
371 # convert back and forth between the orthonormal
372 # coordinates and the given ones, we need to stick the
373 # original basis in self._matrix_span.
374 U
= V
.span_of_basis( deortho_vector_basis
, check
=check_axioms
)
375 self
._deortho
_matrix
= matrix
.column( U
.coordinate_vector(q
)
376 for q
in vector_basis
)
379 # Now we actually compute the multiplication and inner-product
380 # tables/matrices using the possibly-orthonormalized basis.
381 self
._inner
_product
_matrix
= matrix
.identity(field
, n
)
383 self
._multiplication
_table
= [ [zed
for j
in range(i
+1)]
386 # Note: the Jordan and inner-products are defined in terms
387 # of the ambient basis. It's important that their arguments
388 # are in ambient coordinates as well.
391 # ortho basis w.r.t. ambient coords
395 # The jordan product returns a matrixy answer, so we
396 # have to convert it to the algebra coordinates.
397 elt
= jordan_product(q_i
, q_j
)
398 elt
= self
._matrix
_span
.coordinate_vector(V(_all2list(elt
)))
399 self
._multiplication
_table
[i
][j
] = self
.from_vector(elt
)
401 if not orthonormalize
:
402 # If we're orthonormalizing the basis with respect
403 # to an inner-product, then the inner-product
404 # matrix with respect to the resulting basis is
405 # just going to be the identity.
406 ip
= inner_product(q_i
, q_j
)
407 self
._inner
_product
_matrix
[i
,j
] = ip
408 self
._inner
_product
_matrix
[j
,i
] = ip
410 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
411 self
._inner
_product
_matrix
.set_immutable()
414 if not self
._is
_jordanian
():
415 raise ValueError("Jordan identity does not hold")
416 if not self
._inner
_product
_is
_associative
():
417 raise ValueError("inner product is not associative")
420 def _coerce_map_from_base_ring(self
):
422 Disable the map from the base ring into the algebra.
424 Performing a nonsense conversion like this automatically
425 is counterpedagogical. The fallback is to try the usual
426 element constructor, which should also fail.
430 sage: from mjo.eja.eja_algebra import random_eja
434 sage: J = random_eja()
436 Traceback (most recent call last):
438 ValueError: not an element of this algebra
444 def product_on_basis(self
, i
, j
):
446 Returns the Jordan product of the `i` and `j`th basis elements.
448 This completely defines the Jordan product on the algebra, and
449 is used direclty by our superclass machinery to implement
454 sage: from mjo.eja.eja_algebra import random_eja
458 sage: J = random_eja()
459 sage: n = J.dimension()
462 sage: bi_bj = J.zero()*J.zero()
464 ....: i = ZZ.random_element(n)
465 ....: j = ZZ.random_element(n)
466 ....: bi = J.monomial(i)
467 ....: bj = J.monomial(j)
468 ....: bi_bj = J.product_on_basis(i,j)
473 # We only stored the lower-triangular portion of the
474 # multiplication table.
476 return self
._multiplication
_table
[i
][j
]
478 return self
._multiplication
_table
[j
][i
]
480 def inner_product(self
, x
, y
):
482 The inner product associated with this Euclidean Jordan algebra.
484 Defaults to the trace inner product, but can be overridden by
485 subclasses if they are sure that the necessary properties are
490 sage: from mjo.eja.eja_algebra import (random_eja,
492 ....: BilinearFormEJA)
496 Our inner product is "associative," which means the following for
497 a symmetric bilinear form::
499 sage: J = random_eja()
500 sage: x,y,z = J.random_elements(3)
501 sage: (x*y).inner_product(z) == y.inner_product(x*z)
506 Ensure that this is the usual inner product for the algebras
509 sage: J = HadamardEJA.random_instance()
510 sage: x,y = J.random_elements(2)
511 sage: actual = x.inner_product(y)
512 sage: expected = x.to_vector().inner_product(y.to_vector())
513 sage: actual == expected
516 Ensure that this is one-half of the trace inner-product in a
517 BilinearFormEJA that isn't just the reals (when ``n`` isn't
518 one). This is in Faraut and Koranyi, and also my "On the
521 sage: J = BilinearFormEJA.random_instance()
522 sage: n = J.dimension()
523 sage: x = J.random_element()
524 sage: y = J.random_element()
525 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
529 B
= self
._inner
_product
_matrix
530 return (B
*x
.to_vector()).inner_product(y
.to_vector())
533 def is_associative(self
):
535 Return whether or not this algebra's Jordan product is associative.
539 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
543 sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False)
544 sage: J.is_associative()
546 sage: x = sum(J.gens())
547 sage: A = x.subalgebra_generated_by(orthonormalize=False)
548 sage: A.is_associative()
552 return "Associative" in self
.category().axioms()
554 def _is_commutative(self
):
556 Whether or not this algebra's multiplication table is commutative.
558 This method should of course always return ``True``, unless
559 this algebra was constructed with ``check_axioms=False`` and
560 passed an invalid multiplication table.
562 return all( x
*y
== y
*x
for x
in self
.gens() for y
in self
.gens() )
564 def _is_jordanian(self
):
566 Whether or not this algebra's multiplication table respects the
567 Jordan identity `(x^{2})(xy) = x(x^{2}y)`.
569 We only check one arrangement of `x` and `y`, so for a
570 ``True`` result to be truly true, you should also check
571 :meth:`_is_commutative`. This method should of course always
572 return ``True``, unless this algebra was constructed with
573 ``check_axioms=False`` and passed an invalid multiplication table.
575 return all( (self
.monomial(i
)**2)*(self
.monomial(i
)*self
.monomial(j
))
577 (self
.monomial(i
))*((self
.monomial(i
)**2)*self
.monomial(j
))
578 for i
in range(self
.dimension())
579 for j
in range(self
.dimension()) )
581 def _jordan_product_is_associative(self
):
583 Return whether or not this algebra's Jordan product is
584 associative; that is, whether or not `x*(y*z) = (x*y)*z`
587 This method should agree with :meth:`is_associative` unless
588 you lied about the value of the ``associative`` parameter
589 when you constructed the algebra.
593 sage: from mjo.eja.eja_algebra import (random_eja,
594 ....: RealSymmetricEJA,
595 ....: ComplexHermitianEJA,
596 ....: QuaternionHermitianEJA)
600 sage: J = RealSymmetricEJA(4, orthonormalize=False)
601 sage: J._jordan_product_is_associative()
603 sage: x = sum(J.gens())
604 sage: A = x.subalgebra_generated_by()
605 sage: A._jordan_product_is_associative()
610 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
611 sage: J._jordan_product_is_associative()
613 sage: x = sum(J.gens())
614 sage: A = x.subalgebra_generated_by(orthonormalize=False)
615 sage: A._jordan_product_is_associative()
620 sage: J = QuaternionHermitianEJA(2)
621 sage: J._jordan_product_is_associative()
623 sage: x = sum(J.gens())
624 sage: A = x.subalgebra_generated_by()
625 sage: A._jordan_product_is_associative()
630 The values we've presupplied to the constructors agree with
633 sage: J = random_eja()
634 sage: J.is_associative() == J._jordan_product_is_associative()
640 # Used to check whether or not something is zero.
643 # I don't know of any examples that make this magnitude
644 # necessary because I don't know how to make an
645 # associative algebra when the element subalgebra
646 # construction is unreliable (as it is over RDF; we can't
647 # find the degree of an element because we can't compute
648 # the rank of a matrix). But even multiplication of floats
649 # is non-associative, so *some* epsilon is needed... let's
650 # just take the one from _inner_product_is_associative?
653 for i
in range(self
.dimension()):
654 for j
in range(self
.dimension()):
655 for k
in range(self
.dimension()):
659 diff
= (x
*y
)*z
- x
*(y
*z
)
661 if diff
.norm() > epsilon
:
666 def _inner_product_is_associative(self
):
668 Return whether or not this algebra's inner product `B` is
669 associative; that is, whether or not `B(xy,z) = B(x,yz)`.
671 This method should of course always return ``True``, unless
672 this algebra was constructed with ``check_axioms=False`` and
673 passed an invalid Jordan or inner-product.
677 # Used to check whether or not something is zero.
680 # This choice is sufficient to allow the construction of
681 # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
684 for i
in range(self
.dimension()):
685 for j
in range(self
.dimension()):
686 for k
in range(self
.dimension()):
690 diff
= (x
*y
).inner_product(z
) - x
.inner_product(y
*z
)
692 if diff
.abs() > epsilon
:
697 def _element_constructor_(self
, elt
):
699 Construct an element of this algebra or a subalgebra from its
700 EJA element, vector, or matrix representation.
702 This gets called only after the parent element _call_ method
703 fails to find a coercion for the argument.
707 sage: from mjo.eja.eja_algebra import (random_eja,
710 ....: RealSymmetricEJA)
714 The identity in `S^n` is converted to the identity in the EJA::
716 sage: J = RealSymmetricEJA(3)
717 sage: I = matrix.identity(QQ,3)
718 sage: J(I) == J.one()
721 This skew-symmetric matrix can't be represented in the EJA::
723 sage: J = RealSymmetricEJA(3)
724 sage: A = matrix(QQ,3, lambda i,j: i-j)
726 Traceback (most recent call last):
728 ValueError: not an element of this algebra
730 Tuples work as well, provided that the matrix basis for the
731 algebra consists of them::
733 sage: J1 = HadamardEJA(3)
734 sage: J2 = RealSymmetricEJA(2)
735 sage: J = cartesian_product([J1,J2])
736 sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
739 Subalgebra elements are embedded into the superalgebra::
741 sage: J = JordanSpinEJA(3)
744 sage: x = sum(J.gens())
745 sage: A = x.subalgebra_generated_by()
751 Ensure that we can convert any element back and forth
752 faithfully between its matrix and algebra representations::
754 sage: J = random_eja()
755 sage: x = J.random_element()
756 sage: J(x.to_matrix()) == x
759 We cannot coerce elements between algebras just because their
760 matrix representations are compatible::
762 sage: J1 = HadamardEJA(3)
763 sage: J2 = JordanSpinEJA(3)
765 Traceback (most recent call last):
767 ValueError: not an element of this algebra
769 Traceback (most recent call last):
771 ValueError: not an element of this algebra
774 msg
= "not an element of this algebra"
775 if elt
in self
.base_ring():
776 # Ensure that no base ring -> algebra coercion is performed
777 # by this method. There's some stupidity in sage that would
778 # otherwise propagate to this method; for example, sage thinks
779 # that the integer 3 belongs to the space of 2-by-2 matrices.
780 raise ValueError(msg
)
782 if hasattr(elt
, 'superalgebra_element'):
783 # Handle subalgebra elements
784 if elt
.parent().superalgebra() == self
:
785 return elt
.superalgebra_element()
787 if hasattr(elt
, 'sparse_vector'):
788 # Convert a vector into a column-matrix. We check for
789 # "sparse_vector" and not "column" because matrices also
790 # have a "column" method.
793 if elt
not in self
.matrix_space():
794 raise ValueError(msg
)
796 # Thanks for nothing! Matrix spaces aren't vector spaces in
797 # Sage, so we have to figure out its matrix-basis coordinates
798 # ourselves. We use the basis space's ring instead of the
799 # element's ring because the basis space might be an algebraic
800 # closure whereas the base ring of the 3-by-3 identity matrix
801 # could be QQ instead of QQbar.
803 # And, we also have to handle Cartesian product bases (when
804 # the matrix basis consists of tuples) here. The "good news"
805 # is that we're already converting everything to long vectors,
806 # and that strategy works for tuples as well.
808 elt
= self
._matrix
_span
.ambient_vector_space()(_all2list(elt
))
811 coords
= self
._matrix
_span
.coordinate_vector(elt
)
812 except ArithmeticError: # vector is not in free module
813 raise ValueError(msg
)
815 return self
.from_vector(coords
)
819 Return a string representation of ``self``.
823 sage: from mjo.eja.eja_algebra import JordanSpinEJA
827 Ensure that it says what we think it says::
829 sage: JordanSpinEJA(2, field=AA)
830 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
831 sage: JordanSpinEJA(3, field=RDF)
832 Euclidean Jordan algebra of dimension 3 over Real Double Field
835 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
836 return fmt
.format(self
.dimension(), self
.base_ring())
840 def characteristic_polynomial_of(self
):
842 Return the algebra's "characteristic polynomial of" function,
843 which is itself a multivariate polynomial that, when evaluated
844 at the coordinates of some algebra element, returns that
845 element's characteristic polynomial.
847 The resulting polynomial has `n+1` variables, where `n` is the
848 dimension of this algebra. The first `n` variables correspond to
849 the coordinates of an algebra element: when evaluated at the
850 coordinates of an algebra element with respect to a certain
851 basis, the result is a univariate polynomial (in the one
852 remaining variable ``t``), namely the characteristic polynomial
857 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
861 The characteristic polynomial in the spin algebra is given in
862 Alizadeh, Example 11.11::
864 sage: J = JordanSpinEJA(3)
865 sage: p = J.characteristic_polynomial_of(); p
866 X0^2 - X1^2 - X2^2 + (-2*t)*X0 + t^2
867 sage: xvec = J.one().to_vector()
871 By definition, the characteristic polynomial is a monic
872 degree-zero polynomial in a rank-zero algebra. Note that
873 Cayley-Hamilton is indeed satisfied since the polynomial
874 ``1`` evaluates to the identity element of the algebra on
877 sage: J = TrivialEJA()
878 sage: J.characteristic_polynomial_of()
885 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
886 a
= self
._charpoly
_coefficients
()
888 # We go to a bit of trouble here to reorder the
889 # indeterminates, so that it's easier to evaluate the
890 # characteristic polynomial at x's coordinates and get back
891 # something in terms of t, which is what we want.
892 S
= PolynomialRing(self
.base_ring(),'t')
896 S
= PolynomialRing(S
, R
.variable_names())
899 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
901 def coordinate_polynomial_ring(self
):
903 The multivariate polynomial ring in which this algebra's
904 :meth:`characteristic_polynomial_of` lives.
908 sage: from mjo.eja.eja_algebra import (HadamardEJA,
909 ....: RealSymmetricEJA)
913 sage: J = HadamardEJA(2)
914 sage: J.coordinate_polynomial_ring()
915 Multivariate Polynomial Ring in X0, X1...
916 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
917 sage: J.coordinate_polynomial_ring()
918 Multivariate Polynomial Ring in X0, X1, X2, X3, X4, X5...
921 var_names
= tuple( "X%d" % z
for z
in range(self
.dimension()) )
922 return PolynomialRing(self
.base_ring(), var_names
)
924 def inner_product(self
, x
, y
):
926 The inner product associated with this Euclidean Jordan algebra.
928 Defaults to the trace inner product, but can be overridden by
929 subclasses if they are sure that the necessary properties are
934 sage: from mjo.eja.eja_algebra import (random_eja,
936 ....: BilinearFormEJA)
940 Our inner product is "associative," which means the following for
941 a symmetric bilinear form::
943 sage: J = random_eja()
944 sage: x,y,z = J.random_elements(3)
945 sage: (x*y).inner_product(z) == y.inner_product(x*z)
950 Ensure that this is the usual inner product for the algebras
953 sage: J = HadamardEJA.random_instance()
954 sage: x,y = J.random_elements(2)
955 sage: actual = x.inner_product(y)
956 sage: expected = x.to_vector().inner_product(y.to_vector())
957 sage: actual == expected
960 Ensure that this is one-half of the trace inner-product in a
961 BilinearFormEJA that isn't just the reals (when ``n`` isn't
962 one). This is in Faraut and Koranyi, and also my "On the
965 sage: J = BilinearFormEJA.random_instance()
966 sage: n = J.dimension()
967 sage: x = J.random_element()
968 sage: y = J.random_element()
969 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
972 B
= self
._inner
_product
_matrix
973 return (B
*x
.to_vector()).inner_product(y
.to_vector())
976 def is_trivial(self
):
978 Return whether or not this algebra is trivial.
980 A trivial algebra contains only the zero element.
984 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
989 sage: J = ComplexHermitianEJA(3)
995 sage: J = TrivialEJA()
1000 return self
.dimension() == 0
1003 def multiplication_table(self
):
1005 Return a visual representation of this algebra's multiplication
1006 table (on basis elements).
1010 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1014 sage: J = JordanSpinEJA(4)
1015 sage: J.multiplication_table()
1016 +----++----+----+----+----+
1017 | * || b0 | b1 | b2 | b3 |
1018 +====++====+====+====+====+
1019 | b0 || b0 | b1 | b2 | b3 |
1020 +----++----+----+----+----+
1021 | b1 || b1 | b0 | 0 | 0 |
1022 +----++----+----+----+----+
1023 | b2 || b2 | 0 | b0 | 0 |
1024 +----++----+----+----+----+
1025 | b3 || b3 | 0 | 0 | b0 |
1026 +----++----+----+----+----+
1029 n
= self
.dimension()
1030 # Prepend the header row.
1031 M
= [["*"] + list(self
.gens())]
1033 # And to each subsequent row, prepend an entry that belongs to
1034 # the left-side "header column."
1035 M
+= [ [self
.monomial(i
)] + [ self
.monomial(i
)*self
.monomial(j
)
1039 return table(M
, header_row
=True, header_column
=True, frame
=True)
1042 def matrix_basis(self
):
1044 Return an (often more natural) representation of this algebras
1045 basis as an ordered tuple of matrices.
1047 Every finite-dimensional Euclidean Jordan Algebra is a, up to
1048 Jordan isomorphism, a direct sum of five simple
1049 algebras---four of which comprise Hermitian matrices. And the
1050 last type of algebra can of course be thought of as `n`-by-`1`
1051 column matrices (ambiguusly called column vectors) to avoid
1052 special cases. As a result, matrices (and column vectors) are
1053 a natural representation format for Euclidean Jordan algebra
1056 But, when we construct an algebra from a basis of matrices,
1057 those matrix representations are lost in favor of coordinate
1058 vectors *with respect to* that basis. We could eventually
1059 convert back if we tried hard enough, but having the original
1060 representations handy is valuable enough that we simply store
1061 them and return them from this method.
1063 Why implement this for non-matrix algebras? Avoiding special
1064 cases for the :class:`BilinearFormEJA` pays with simplicity in
1065 its own right. But mainly, we would like to be able to assume
1066 that elements of a :class:`CartesianProductEJA` can be displayed
1067 nicely, without having to have special classes for direct sums
1068 one of whose components was a matrix algebra.
1072 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1073 ....: RealSymmetricEJA)
1077 sage: J = RealSymmetricEJA(2)
1079 Finite family {0: b0, 1: b1, 2: b2}
1080 sage: J.matrix_basis()
1082 [1 0] [ 0 0.7071067811865475?] [0 0]
1083 [0 0], [0.7071067811865475? 0], [0 1]
1088 sage: J = JordanSpinEJA(2)
1090 Finite family {0: b0, 1: b1}
1091 sage: J.matrix_basis()
1097 return self
._matrix
_basis
1100 def matrix_space(self
):
1102 Return the matrix space in which this algebra's elements live, if
1103 we think of them as matrices (including column vectors of the
1106 "By default" this will be an `n`-by-`1` column-matrix space,
1107 except when the algebra is trivial. There it's `n`-by-`n`
1108 (where `n` is zero), to ensure that two elements of the matrix
1109 space (empty matrices) can be multiplied. For algebras of
1110 matrices, this returns the space in which their
1111 real embeddings live.
1115 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
1116 ....: JordanSpinEJA,
1117 ....: QuaternionHermitianEJA,
1122 By default, the matrix representation is just a column-matrix
1123 equivalent to the vector representation::
1125 sage: J = JordanSpinEJA(3)
1126 sage: J.matrix_space()
1127 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
1130 The matrix representation in the trivial algebra is
1131 zero-by-zero instead of the usual `n`-by-one::
1133 sage: J = TrivialEJA()
1134 sage: J.matrix_space()
1135 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
1138 The matrix space for complex/quaternion Hermitian matrix EJA
1139 is the space in which their real-embeddings live, not the
1140 original complex/quaternion matrix space::
1142 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1143 sage: J.matrix_space()
1144 Module of 2 by 2 matrices with entries in Algebraic Field over
1145 the scalar ring Rational Field
1146 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1147 sage: J.matrix_space()
1148 Module of 1 by 1 matrices with entries in Quaternion
1149 Algebra (-1, -1) with base ring Rational Field over
1150 the scalar ring Rational Field
1153 return self
._matrix
_space
1159 Return the unit element of this algebra.
1163 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1168 We can compute unit element in the Hadamard EJA::
1170 sage: J = HadamardEJA(5)
1172 b0 + b1 + b2 + b3 + b4
1174 The unit element in the Hadamard EJA is inherited in the
1175 subalgebras generated by its elements::
1177 sage: J = HadamardEJA(5)
1179 b0 + b1 + b2 + b3 + b4
1180 sage: x = sum(J.gens())
1181 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1184 sage: A.one().superalgebra_element()
1185 b0 + b1 + b2 + b3 + b4
1189 The identity element acts like the identity, regardless of
1190 whether or not we orthonormalize::
1192 sage: J = random_eja()
1193 sage: x = J.random_element()
1194 sage: J.one()*x == x and x*J.one() == x
1196 sage: A = x.subalgebra_generated_by()
1197 sage: y = A.random_element()
1198 sage: A.one()*y == y and y*A.one() == y
1203 sage: J = random_eja(field=QQ, orthonormalize=False)
1204 sage: x = J.random_element()
1205 sage: J.one()*x == x and x*J.one() == x
1207 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1208 sage: y = A.random_element()
1209 sage: A.one()*y == y and y*A.one() == y
1212 The matrix of the unit element's operator is the identity,
1213 regardless of the base field and whether or not we
1216 sage: J = random_eja()
1217 sage: actual = J.one().operator().matrix()
1218 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1219 sage: actual == expected
1221 sage: x = J.random_element()
1222 sage: A = x.subalgebra_generated_by()
1223 sage: actual = A.one().operator().matrix()
1224 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1225 sage: actual == expected
1230 sage: J = random_eja(field=QQ, orthonormalize=False)
1231 sage: actual = J.one().operator().matrix()
1232 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1233 sage: actual == expected
1235 sage: x = J.random_element()
1236 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1237 sage: actual = A.one().operator().matrix()
1238 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1239 sage: actual == expected
1242 Ensure that the cached unit element (often precomputed by
1243 hand) agrees with the computed one::
1245 sage: J = random_eja()
1246 sage: cached = J.one()
1247 sage: J.one.clear_cache()
1248 sage: J.one() == cached
1253 sage: J = random_eja(field=QQ, orthonormalize=False)
1254 sage: cached = J.one()
1255 sage: J.one.clear_cache()
1256 sage: J.one() == cached
1260 # We can brute-force compute the matrices of the operators
1261 # that correspond to the basis elements of this algebra.
1262 # If some linear combination of those basis elements is the
1263 # algebra identity, then the same linear combination of
1264 # their matrices has to be the identity matrix.
1266 # Of course, matrices aren't vectors in sage, so we have to
1267 # appeal to the "long vectors" isometry.
1269 V
= VectorSpace(self
.base_ring(), self
.dimension()**2)
1270 oper_vecs
= [ V(g
.operator().matrix().list()) for g
in self
.gens() ]
1272 # Now we use basic linear algebra to find the coefficients,
1273 # of the matrices-as-vectors-linear-combination, which should
1274 # work for the original algebra basis too.
1275 A
= matrix(self
.base_ring(), oper_vecs
)
1277 # We used the isometry on the left-hand side already, but we
1278 # still need to do it for the right-hand side. Recall that we
1279 # wanted something that summed to the identity matrix.
1280 b
= V( matrix
.identity(self
.base_ring(), self
.dimension()).list() )
1282 # Now if there's an identity element in the algebra, this
1283 # should work. We solve on the left to avoid having to
1284 # transpose the matrix "A".
1285 return self
.from_vector(A
.solve_left(b
))
1288 def peirce_decomposition(self
, c
):
1290 The Peirce decomposition of this algebra relative to the
1293 In the future, this can be extended to a complete system of
1294 orthogonal idempotents.
1298 - ``c`` -- an idempotent of this algebra.
1302 A triple (J0, J5, J1) containing two subalgebras and one subspace
1305 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1306 corresponding to the eigenvalue zero.
1308 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1309 corresponding to the eigenvalue one-half.
1311 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1312 corresponding to the eigenvalue one.
1314 These are the only possible eigenspaces for that operator, and this
1315 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1316 orthogonal, and are subalgebras of this algebra with the appropriate
1321 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1325 The canonical example comes from the symmetric matrices, which
1326 decompose into diagonal and off-diagonal parts::
1328 sage: J = RealSymmetricEJA(3)
1329 sage: C = matrix(QQ, [ [1,0,0],
1333 sage: J0,J5,J1 = J.peirce_decomposition(c)
1335 Euclidean Jordan algebra of dimension 1...
1337 Vector space of degree 6 and dimension 2...
1339 Euclidean Jordan algebra of dimension 3...
1340 sage: J0.one().to_matrix()
1344 sage: orig_df = AA.options.display_format
1345 sage: AA.options.display_format = 'radical'
1346 sage: J.from_vector(J5.basis()[0]).to_matrix()
1350 sage: J.from_vector(J5.basis()[1]).to_matrix()
1354 sage: AA.options.display_format = orig_df
1355 sage: J1.one().to_matrix()
1362 Every algebra decomposes trivially with respect to its identity
1365 sage: J = random_eja()
1366 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1367 sage: J0.dimension() == 0 and J5.dimension() == 0
1369 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1372 The decomposition is into eigenspaces, and its components are
1373 therefore necessarily orthogonal. Moreover, the identity
1374 elements in the two subalgebras are the projections onto their
1375 respective subspaces of the superalgebra's identity element::
1377 sage: J = random_eja()
1378 sage: x = J.random_element()
1379 sage: if not J.is_trivial():
1380 ....: while x.is_nilpotent():
1381 ....: x = J.random_element()
1382 sage: c = x.subalgebra_idempotent()
1383 sage: J0,J5,J1 = J.peirce_decomposition(c)
1385 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1386 ....: w = w.superalgebra_element()
1387 ....: y = J.from_vector(y)
1388 ....: z = z.superalgebra_element()
1389 ....: ipsum += w.inner_product(y).abs()
1390 ....: ipsum += w.inner_product(z).abs()
1391 ....: ipsum += y.inner_product(z).abs()
1394 sage: J1(c) == J1.one()
1396 sage: J0(J.one() - c) == J0.one()
1400 if not c
.is_idempotent():
1401 raise ValueError("element is not idempotent: %s" % c
)
1403 # Default these to what they should be if they turn out to be
1404 # trivial, because eigenspaces_left() won't return eigenvalues
1405 # corresponding to trivial spaces (e.g. it returns only the
1406 # eigenspace corresponding to lambda=1 if you take the
1407 # decomposition relative to the identity element).
1408 trivial
= self
.subalgebra((), check_axioms
=False)
1409 J0
= trivial
# eigenvalue zero
1410 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1411 J1
= trivial
# eigenvalue one
1413 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1414 if eigval
== ~
(self
.base_ring()(2)):
1417 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1418 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1424 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1429 def random_element(self
, thorough
=False):
1431 Return a random element of this algebra.
1433 Our algebra superclass method only returns a linear
1434 combination of at most two basis elements. We instead
1435 want the vector space "random element" method that
1436 returns a more diverse selection.
1440 - ``thorough`` -- (boolean; default False) whether or not we
1441 should generate irrational coefficients for the random
1442 element when our base ring is irrational; this slows the
1443 algebra operations to a crawl, but any truly random method
1447 # For a general base ring... maybe we can trust this to do the
1448 # right thing? Unlikely, but.
1449 V
= self
.vector_space()
1450 if self
.base_ring() is AA
and not thorough
:
1451 # Now that AA generates actually random random elements
1452 # (post Trac 30875), we only need to de-thorough the
1453 # randomness when asked to.
1454 V
= V
.change_ring(QQ
)
1456 v
= V
.random_element()
1457 return self
.from_vector(V
.coordinate_vector(v
))
1459 def random_elements(self
, count
, thorough
=False):
1461 Return ``count`` random elements as a tuple.
1465 - ``thorough`` -- (boolean; default False) whether or not we
1466 should generate irrational coefficients for the random
1467 elements when our base ring is irrational; this slows the
1468 algebra operations to a crawl, but any truly random method
1473 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1477 sage: J = JordanSpinEJA(3)
1478 sage: x,y,z = J.random_elements(3)
1479 sage: all( [ x in J, y in J, z in J ])
1481 sage: len( J.random_elements(10) ) == 10
1485 return tuple( self
.random_element(thorough
)
1486 for idx
in range(count
) )
1489 def operator_polynomial_matrix(self
):
1491 Return the matrix of polynomials (over this algebra's
1492 :meth:`coordinate_polynomial_ring`) that, when evaluated at
1493 the basis coordinates of an element `x`, produces the basis
1494 representation of `L_{x}`.
1498 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1499 ....: JordanSpinEJA)
1503 sage: J = HadamardEJA(4)
1504 sage: L_x = J.operator_polynomial_matrix()
1511 sage: d = zip(J.coordinate_polynomial_ring().gens(), x.to_vector())
1512 sage: L_x.subs(dict(d))
1520 sage: J = JordanSpinEJA(4)
1521 sage: L_x = J.operator_polynomial_matrix()
1528 sage: d = zip(J.coordinate_polynomial_ring().gens(), x.to_vector())
1529 sage: L_x.subs(dict(d))
1536 R
= self
.coordinate_polynomial_ring()
1539 # From a result in my book, these are the entries of the
1540 # basis representation of L_x.
1541 return sum( v
*self
.monomial(k
).operator().matrix()[i
,j
]
1542 for (k
,v
) in enumerate(R
.gens()) )
1544 n
= self
.dimension()
1545 return matrix(R
, n
, n
, L_x_i_j
)
1548 def _charpoly_coefficients(self
):
1550 The `r` polynomial coefficients of the "characteristic polynomial
1555 sage: from mjo.eja.eja_algebra import random_eja
1559 The theory shows that these are all homogeneous polynomials of
1562 sage: J = random_eja()
1563 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1567 n
= self
.dimension()
1568 R
= self
.coordinate_polynomial_ring()
1569 F
= R
.fraction_field()
1571 L_x
= self
.operator_polynomial_matrix()
1574 if self
.rank
.is_in_cache():
1576 # There's no need to pad the system with redundant
1577 # columns if we *know* they'll be redundant.
1580 # Compute an extra power in case the rank is equal to
1581 # the dimension (otherwise, we would stop at x^(r-1)).
1582 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1583 for k
in range(n
+1) ]
1584 A
= matrix
.column(F
, x_powers
[:n
])
1585 AE
= A
.extended_echelon_form()
1592 # The theory says that only the first "r" coefficients are
1593 # nonzero, and they actually live in the original polynomial
1594 # ring and not the fraction field. We negate them because in
1595 # the actual characteristic polynomial, they get moved to the
1596 # other side where x^r lives. We don't bother to trim A_rref
1597 # down to a square matrix and solve the resulting system,
1598 # because the upper-left r-by-r portion of A_rref is
1599 # guaranteed to be the identity matrix, so e.g.
1601 # A_rref.solve_right(Y)
1603 # would just be returning Y.
1604 return (-E
*b
)[:r
].change_ring(R
)
1609 Return the rank of this EJA.
1611 This is a cached method because we know the rank a priori for
1612 all of the algebras we can construct. Thus we can avoid the
1613 expensive ``_charpoly_coefficients()`` call unless we truly
1614 need to compute the whole characteristic polynomial.
1618 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1619 ....: JordanSpinEJA,
1620 ....: RealSymmetricEJA,
1621 ....: ComplexHermitianEJA,
1622 ....: QuaternionHermitianEJA,
1627 The rank of the Jordan spin algebra is always two::
1629 sage: JordanSpinEJA(2).rank()
1631 sage: JordanSpinEJA(3).rank()
1633 sage: JordanSpinEJA(4).rank()
1636 The rank of the `n`-by-`n` Hermitian real, complex, or
1637 quaternion matrices is `n`::
1639 sage: RealSymmetricEJA(4).rank()
1641 sage: ComplexHermitianEJA(3).rank()
1643 sage: QuaternionHermitianEJA(2).rank()
1648 Ensure that every EJA that we know how to construct has a
1649 positive integer rank, unless the algebra is trivial in
1650 which case its rank will be zero::
1652 sage: J = random_eja()
1656 sage: r > 0 or (r == 0 and J.is_trivial())
1659 Ensure that computing the rank actually works, since the ranks
1660 of all simple algebras are known and will be cached by default::
1662 sage: J = random_eja() # long time
1663 sage: cached = J.rank() # long time
1664 sage: J.rank.clear_cache() # long time
1665 sage: J.rank() == cached # long time
1669 return len(self
._charpoly
_coefficients
())
1672 def subalgebra(self
, basis
, **kwargs
):
1674 Create a subalgebra of this algebra from the given basis.
1676 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1677 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1680 def vector_space(self
):
1682 Return the vector space that underlies this algebra.
1686 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1690 sage: J = RealSymmetricEJA(2)
1691 sage: J.vector_space()
1692 Vector space of dimension 3 over...
1695 return self
.zero().to_vector().parent().ambient_vector_space()
1699 class RationalBasisEJA(FiniteDimensionalEJA
):
1701 Algebras whose supplied basis elements have all rational entries.
1705 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1709 The supplied basis is orthonormalized by default::
1711 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1712 sage: J = BilinearFormEJA(B)
1713 sage: J.matrix_basis()
1730 # Abuse the check_field parameter to check that the entries of
1731 # out basis (in ambient coordinates) are in the field QQ.
1732 # Use _all2list to get the vector coordinates of octonion
1733 # entries and not the octonions themselves (which are not
1735 if not all( all(b_i
in QQ
for b_i
in _all2list(b
))
1737 raise TypeError("basis not rational")
1739 super().__init
__(basis
,
1743 check_field
=check_field
,
1746 self
._rational
_algebra
= None
1748 # There's no point in constructing the extra algebra if this
1749 # one is already rational.
1751 # Note: the same Jordan and inner-products work here,
1752 # because they are necessarily defined with respect to
1753 # ambient coordinates and not any particular basis.
1754 self
._rational
_algebra
= FiniteDimensionalEJA(
1759 matrix_space
=self
.matrix_space(),
1760 associative
=self
.is_associative(),
1761 orthonormalize
=False,
1765 def rational_algebra(self
):
1766 # Using None as a flag here (rather than just assigning "self"
1767 # to self._rational_algebra by default) feels a little bit
1768 # more sane to me in a garbage-collected environment.
1769 if self
._rational
_algebra
is None:
1772 return self
._rational
_algebra
1775 def _charpoly_coefficients(self
):
1779 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1780 ....: JordanSpinEJA)
1784 The base ring of the resulting polynomial coefficients is what
1785 it should be, and not the rationals (unless the algebra was
1786 already over the rationals)::
1788 sage: J = JordanSpinEJA(3)
1789 sage: J._charpoly_coefficients()
1790 (X0^2 - X1^2 - X2^2, -2*X0)
1791 sage: a0 = J._charpoly_coefficients()[0]
1793 Algebraic Real Field
1794 sage: a0.base_ring()
1795 Algebraic Real Field
1798 if self
.rational_algebra() is self
:
1799 # Bypass the hijinks if they won't benefit us.
1800 return super()._charpoly
_coefficients
()
1802 # Do the computation over the rationals. The answer will be
1803 # the same, because all we've done is a change of basis.
1804 # Then, change back from QQ to our real base ring
1805 a
= ( a_i
.change_ring(self
.base_ring())
1806 for a_i
in self
.rational_algebra()._charpoly
_coefficients
() )
1808 # Otherwise, convert the coordinate variables back to the
1809 # deorthonormalized ones.
1810 R
= self
.coordinate_polynomial_ring()
1811 from sage
.modules
.free_module_element
import vector
1812 X
= vector(R
, R
.gens())
1813 BX
= self
._deortho
_matrix
*X
1815 subs_dict
= { X[i]: BX[i] for i in range(len(X)) }
1816 return tuple( a_i
.subs(subs_dict
) for a_i
in a
)
1818 class ConcreteEJA(FiniteDimensionalEJA
):
1820 A class for the Euclidean Jordan algebras that we know by name.
1822 These are the Jordan algebras whose basis, multiplication table,
1823 rank, and so on are known a priori. More to the point, they are
1824 the Euclidean Jordan algebras for which we are able to conjure up
1825 a "random instance."
1829 sage: from mjo.eja.eja_algebra import ConcreteEJA
1833 Our basis is normalized with respect to the algebra's inner
1834 product, unless we specify otherwise::
1836 sage: J = ConcreteEJA.random_instance()
1837 sage: all( b.norm() == 1 for b in J.gens() )
1840 Since our basis is orthonormal with respect to the algebra's inner
1841 product, and since we know that this algebra is an EJA, any
1842 left-multiplication operator's matrix will be symmetric because
1843 natural->EJA basis representation is an isometry and within the
1844 EJA the operator is self-adjoint by the Jordan axiom::
1846 sage: J = ConcreteEJA.random_instance()
1847 sage: x = J.random_element()
1848 sage: x.operator().is_self_adjoint()
1853 def _max_random_instance_dimension():
1855 The maximum dimension of any random instance. Ten dimensions seems
1856 to be about the point where everything takes a turn for the
1857 worse. And dimension ten (but not nine) allows the 4-by-4 real
1858 Hermitian matrices, the 2-by-2 quaternion Hermitian matrices,
1859 and the 2-by-2 octonion Hermitian matrices.
1864 def _max_random_instance_size(max_dimension
):
1866 Return an integer "size" that is an upper bound on the size of
1867 this algebra when it is used in a random test case. This size
1868 (which can be passed to the algebra's constructor) is itself
1869 based on the ``max_dimension`` parameter.
1871 This method must be implemented in each subclass.
1873 raise NotImplementedError
1876 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
1878 Return a random instance of this type of algebra whose dimension
1879 is less than or equal to the lesser of ``max_dimension`` and
1880 the value returned by ``_max_random_instance_dimension()``. If
1881 the dimension bound is omitted, then only the
1882 ``_max_random_instance_dimension()`` is used as a bound.
1884 This method should be implemented in each subclass.
1888 sage: from mjo.eja.eja_algebra import ConcreteEJA
1892 Both the class bound and the ``max_dimension`` argument are upper
1893 bounds on the dimension of the algebra returned::
1895 sage: from sage.misc.prandom import choice
1896 sage: eja_class = choice(ConcreteEJA.__subclasses__())
1897 sage: class_max_d = eja_class._max_random_instance_dimension()
1898 sage: J = eja_class.random_instance(max_dimension=20,
1900 ....: orthonormalize=False)
1901 sage: J.dimension() <= class_max_d
1903 sage: J = eja_class.random_instance(max_dimension=2,
1905 ....: orthonormalize=False)
1906 sage: J.dimension() <= 2
1910 from sage
.misc
.prandom
import choice
1911 eja_class
= choice(cls
.__subclasses
__())
1913 # These all bubble up to the RationalBasisEJA superclass
1914 # constructor, so any (kw)args valid there are also valid
1916 return eja_class
.random_instance(max_dimension
, *args
, **kwargs
)
1919 class HermitianMatrixEJA(FiniteDimensionalEJA
):
1921 def _denormalized_basis(A
):
1923 Returns a basis for the given Hermitian matrix space.
1925 Why do we embed these? Basically, because all of numerical linear
1926 algebra assumes that you're working with vectors consisting of `n`
1927 entries from a field and scalars from the same field. There's no way
1928 to tell SageMath that (for example) the vectors contain complex
1929 numbers, while the scalar field is real.
1933 sage: from mjo.hurwitz import (ComplexMatrixAlgebra,
1934 ....: QuaternionMatrixAlgebra,
1935 ....: OctonionMatrixAlgebra)
1936 sage: from mjo.eja.eja_algebra import HermitianMatrixEJA
1940 sage: n = ZZ.random_element(1,5)
1941 sage: A = MatrixSpace(QQ, n)
1942 sage: B = HermitianMatrixEJA._denormalized_basis(A)
1943 sage: all( M.is_hermitian() for M in B)
1948 sage: n = ZZ.random_element(1,5)
1949 sage: A = ComplexMatrixAlgebra(n, scalars=QQ)
1950 sage: B = HermitianMatrixEJA._denormalized_basis(A)
1951 sage: all( M.is_hermitian() for M in B)
1956 sage: n = ZZ.random_element(1,5)
1957 sage: A = QuaternionMatrixAlgebra(n, scalars=QQ)
1958 sage: B = HermitianMatrixEJA._denormalized_basis(A)
1959 sage: all( M.is_hermitian() for M in B )
1964 sage: n = ZZ.random_element(1,5)
1965 sage: A = OctonionMatrixAlgebra(n, scalars=QQ)
1966 sage: B = HermitianMatrixEJA._denormalized_basis(A)
1967 sage: all( M.is_hermitian() for M in B )
1971 # These work for real MatrixSpace, whose monomials only have
1972 # two coordinates (because the last one would always be "1").
1973 es
= A
.base_ring().gens()
1974 gen
= lambda A
,m
: A
.monomial(m
[:2])
1976 if hasattr(A
, 'entry_algebra_gens'):
1977 # We've got a MatrixAlgebra, and its monomials will have
1978 # three coordinates.
1979 es
= A
.entry_algebra_gens()
1980 gen
= lambda A
,m
: A
.monomial(m
)
1983 for i
in range(A
.nrows()):
1984 for j
in range(i
+1):
1986 E_ii
= gen(A
, (i
,j
,es
[0]))
1990 E_ij
= gen(A
, (i
,j
,e
))
1991 E_ij
+= E_ij
.conjugate_transpose()
1994 return tuple( basis
)
1997 def jordan_product(X
,Y
):
1998 return (X
*Y
+ Y
*X
)/2
2001 def trace_inner_product(X
,Y
):
2003 A trace inner-product for matrices that aren't embedded in the
2004 reals. It takes MATRICES as arguments, not EJA elements.
2008 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
2009 ....: ComplexHermitianEJA,
2010 ....: QuaternionHermitianEJA,
2011 ....: OctonionHermitianEJA)
2015 sage: J = RealSymmetricEJA(2,field=QQ,orthonormalize=False)
2016 sage: I = J.one().to_matrix()
2017 sage: J.trace_inner_product(I, -I)
2022 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
2023 sage: I = J.one().to_matrix()
2024 sage: J.trace_inner_product(I, -I)
2029 sage: J = QuaternionHermitianEJA(2,field=QQ,orthonormalize=False)
2030 sage: I = J.one().to_matrix()
2031 sage: J.trace_inner_product(I, -I)
2036 sage: J = OctonionHermitianEJA(2,field=QQ,orthonormalize=False)
2037 sage: I = J.one().to_matrix()
2038 sage: J.trace_inner_product(I, -I)
2043 if hasattr(tr
, 'coefficient'):
2044 # Works for octonions, and has to come first because they
2045 # also have a "real()" method that doesn't return an
2046 # element of the scalar ring.
2047 return tr
.coefficient(0)
2048 elif hasattr(tr
, 'coefficient_tuple'):
2049 # Works for quaternions.
2050 return tr
.coefficient_tuple()[0]
2052 # Works for real and complex numbers.
2056 def __init__(self
, matrix_space
, **kwargs
):
2057 # We know this is a valid EJA, but will double-check
2058 # if the user passes check_axioms=True.
2059 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2061 super().__init
__(self
._denormalized
_basis
(matrix_space
),
2062 self
.jordan_product
,
2063 self
.trace_inner_product
,
2064 field
=matrix_space
.base_ring(),
2065 matrix_space
=matrix_space
,
2068 self
.rank
.set_cache(matrix_space
.nrows())
2069 self
.one
.set_cache( self(matrix_space
.one()) )
2071 class RealSymmetricEJA(HermitianMatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2073 The rank-n simple EJA consisting of real symmetric n-by-n
2074 matrices, the usual symmetric Jordan product, and the trace inner
2075 product. It has dimension `(n^2 + n)/2` over the reals.
2079 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
2083 sage: J = RealSymmetricEJA(2)
2084 sage: b0, b1, b2 = J.gens()
2092 In theory, our "field" can be any subfield of the reals::
2094 sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
2095 Euclidean Jordan algebra of dimension 3 over Real Double Field
2096 sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
2097 Euclidean Jordan algebra of dimension 3 over Real Field with
2098 53 bits of precision
2102 The dimension of this algebra is `(n^2 + n) / 2`::
2104 sage: d = RealSymmetricEJA._max_random_instance_dimension()
2105 sage: n = RealSymmetricEJA._max_random_instance_size(d)
2106 sage: J = RealSymmetricEJA(n)
2107 sage: J.dimension() == (n^2 + n)/2
2110 The Jordan multiplication is what we think it is::
2112 sage: J = RealSymmetricEJA.random_instance()
2113 sage: x,y = J.random_elements(2)
2114 sage: actual = (x*y).to_matrix()
2115 sage: X = x.to_matrix()
2116 sage: Y = y.to_matrix()
2117 sage: expected = (X*Y + Y*X)/2
2118 sage: actual == expected
2120 sage: J(expected) == x*y
2123 We can change the generator prefix::
2125 sage: RealSymmetricEJA(3, prefix='q').gens()
2126 (q0, q1, q2, q3, q4, q5)
2128 We can construct the (trivial) algebra of rank zero::
2130 sage: RealSymmetricEJA(0)
2131 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2135 def _max_random_instance_size(max_dimension
):
2136 # Obtained by solving d = (n^2 + n)/2.
2137 # The ZZ-int-ZZ thing is just "floor."
2138 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/2 - 1/2))
2141 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2143 Return a random instance of this type of algebra.
2145 class_max_d
= cls
._max
_random
_instance
_dimension
()
2146 if (max_dimension
is None or max_dimension
> class_max_d
):
2147 max_dimension
= class_max_d
2148 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2149 n
= ZZ
.random_element(max_size
+ 1)
2150 return cls(n
, **kwargs
)
2152 def __init__(self
, n
, field
=AA
, **kwargs
):
2153 A
= MatrixSpace(field
, n
)
2154 super().__init
__(A
, **kwargs
)
2156 from mjo
.eja
.eja_cache
import real_symmetric_eja_coeffs
2157 a
= real_symmetric_eja_coeffs(self
)
2159 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2163 class ComplexHermitianEJA(HermitianMatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2165 The rank-n simple EJA consisting of complex Hermitian n-by-n
2166 matrices over the real numbers, the usual symmetric Jordan product,
2167 and the real-part-of-trace inner product. It has dimension `n^2` over
2172 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2176 In theory, our "field" can be any subfield of the reals, but we
2177 can't use inexact real fields at the moment because SageMath
2178 doesn't know how to convert their elements into complex numbers,
2179 or even into algebraic reals::
2182 Traceback (most recent call last):
2184 TypeError: Illegal initializer for algebraic number
2186 Traceback (most recent call last):
2188 TypeError: Illegal initializer for algebraic number
2190 This causes the following error when we try to scale a matrix of
2191 complex numbers by an inexact real number::
2193 sage: ComplexHermitianEJA(2,field=RR)
2194 Traceback (most recent call last):
2196 TypeError: Unable to coerce entries (=(1.00000000000000,
2197 -0.000000000000000)) to coefficients in Algebraic Real Field
2201 The dimension of this algebra is `n^2`::
2203 sage: d = ComplexHermitianEJA._max_random_instance_dimension()
2204 sage: n = ComplexHermitianEJA._max_random_instance_size(d)
2205 sage: J = ComplexHermitianEJA(n)
2206 sage: J.dimension() == n^2
2209 The Jordan multiplication is what we think it is::
2211 sage: J = ComplexHermitianEJA.random_instance()
2212 sage: x,y = J.random_elements(2)
2213 sage: actual = (x*y).to_matrix()
2214 sage: X = x.to_matrix()
2215 sage: Y = y.to_matrix()
2216 sage: expected = (X*Y + Y*X)/2
2217 sage: actual == expected
2219 sage: J(expected) == x*y
2222 We can change the generator prefix::
2224 sage: ComplexHermitianEJA(2, prefix='z').gens()
2227 We can construct the (trivial) algebra of rank zero::
2229 sage: ComplexHermitianEJA(0)
2230 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2233 def __init__(self
, n
, field
=AA
, **kwargs
):
2234 from mjo
.hurwitz
import ComplexMatrixAlgebra
2235 A
= ComplexMatrixAlgebra(n
, scalars
=field
)
2236 super().__init
__(A
, **kwargs
)
2238 from mjo
.eja
.eja_cache
import complex_hermitian_eja_coeffs
2239 a
= complex_hermitian_eja_coeffs(self
)
2241 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2244 def _max_random_instance_size(max_dimension
):
2245 # Obtained by solving d = n^2.
2246 # The ZZ-int-ZZ thing is just "floor."
2247 return ZZ(int(ZZ(max_dimension
).sqrt()))
2250 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2252 Return a random instance of this type of algebra.
2254 class_max_d
= cls
._max
_random
_instance
_dimension
()
2255 if (max_dimension
is None or max_dimension
> class_max_d
):
2256 max_dimension
= class_max_d
2257 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2258 n
= ZZ
.random_element(max_size
+ 1)
2259 return cls(n
, **kwargs
)
2262 class QuaternionHermitianEJA(HermitianMatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2264 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2265 matrices, the usual symmetric Jordan product, and the
2266 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2271 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2275 In theory, our "field" can be any subfield of the reals::
2277 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2278 Euclidean Jordan algebra of dimension 6 over Real Double Field
2279 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2280 Euclidean Jordan algebra of dimension 6 over Real Field with
2281 53 bits of precision
2285 The dimension of this algebra is `2*n^2 - n`::
2287 sage: d = QuaternionHermitianEJA._max_random_instance_dimension()
2288 sage: n = QuaternionHermitianEJA._max_random_instance_size(d)
2289 sage: J = QuaternionHermitianEJA(n)
2290 sage: J.dimension() == 2*(n^2) - n
2293 The Jordan multiplication is what we think it is::
2295 sage: J = QuaternionHermitianEJA.random_instance()
2296 sage: x,y = J.random_elements(2)
2297 sage: actual = (x*y).to_matrix()
2298 sage: X = x.to_matrix()
2299 sage: Y = y.to_matrix()
2300 sage: expected = (X*Y + Y*X)/2
2301 sage: actual == expected
2303 sage: J(expected) == x*y
2306 We can change the generator prefix::
2308 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2309 (a0, a1, a2, a3, a4, a5)
2311 We can construct the (trivial) algebra of rank zero::
2313 sage: QuaternionHermitianEJA(0)
2314 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2317 def __init__(self
, n
, field
=AA
, **kwargs
):
2318 from mjo
.hurwitz
import QuaternionMatrixAlgebra
2319 A
= QuaternionMatrixAlgebra(n
, scalars
=field
)
2320 super().__init
__(A
, **kwargs
)
2322 from mjo
.eja
.eja_cache
import quaternion_hermitian_eja_coeffs
2323 a
= quaternion_hermitian_eja_coeffs(self
)
2325 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2330 def _max_random_instance_size(max_dimension
):
2332 The maximum rank of a random QuaternionHermitianEJA.
2334 # Obtained by solving d = 2n^2 - n.
2335 # The ZZ-int-ZZ thing is just "floor."
2336 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/4 + 1/4))
2339 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2341 Return a random instance of this type of algebra.
2343 class_max_d
= cls
._max
_random
_instance
_dimension
()
2344 if (max_dimension
is None or max_dimension
> class_max_d
):
2345 max_dimension
= class_max_d
2346 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2347 n
= ZZ
.random_element(max_size
+ 1)
2348 return cls(n
, **kwargs
)
2350 class OctonionHermitianEJA(HermitianMatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2354 sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA,
2355 ....: OctonionHermitianEJA)
2356 sage: from mjo.hurwitz import Octonions, OctonionMatrixAlgebra
2360 The 3-by-3 algebra satisfies the axioms of an EJA::
2362 sage: OctonionHermitianEJA(3, # long time
2363 ....: field=QQ, # long time
2364 ....: orthonormalize=False, # long time
2365 ....: check_axioms=True) # long time
2366 Euclidean Jordan algebra of dimension 27 over Rational Field
2368 After a change-of-basis, the 2-by-2 algebra has the same
2369 multiplication table as the ten-dimensional Jordan spin algebra::
2371 sage: A = OctonionMatrixAlgebra(2,Octonions(QQ),QQ)
2372 sage: b = OctonionHermitianEJA._denormalized_basis(A)
2373 sage: basis = (b[0] + b[9],) + b[1:9] + (b[0] - b[9],)
2374 sage: jp = OctonionHermitianEJA.jordan_product
2375 sage: ip = OctonionHermitianEJA.trace_inner_product
2376 sage: J = FiniteDimensionalEJA(basis,
2380 ....: orthonormalize=False)
2381 sage: J.multiplication_table()
2382 +----++----+----+----+----+----+----+----+----+----+----+
2383 | * || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 |
2384 +====++====+====+====+====+====+====+====+====+====+====+
2385 | b0 || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 |
2386 +----++----+----+----+----+----+----+----+----+----+----+
2387 | b1 || b1 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
2388 +----++----+----+----+----+----+----+----+----+----+----+
2389 | b2 || b2 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
2390 +----++----+----+----+----+----+----+----+----+----+----+
2391 | b3 || b3 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 |
2392 +----++----+----+----+----+----+----+----+----+----+----+
2393 | b4 || b4 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 |
2394 +----++----+----+----+----+----+----+----+----+----+----+
2395 | b5 || b5 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 |
2396 +----++----+----+----+----+----+----+----+----+----+----+
2397 | b6 || b6 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 |
2398 +----++----+----+----+----+----+----+----+----+----+----+
2399 | b7 || b7 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 |
2400 +----++----+----+----+----+----+----+----+----+----+----+
2401 | b8 || b8 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 |
2402 +----++----+----+----+----+----+----+----+----+----+----+
2403 | b9 || b9 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 |
2404 +----++----+----+----+----+----+----+----+----+----+----+
2408 We can actually construct the 27-dimensional Albert algebra,
2409 and we get the right unit element if we recompute it::
2411 sage: J = OctonionHermitianEJA(3, # long time
2412 ....: field=QQ, # long time
2413 ....: orthonormalize=False) # long time
2414 sage: J.one.clear_cache() # long time
2415 sage: J.one() # long time
2417 sage: J.one().to_matrix() # long time
2426 The 2-by-2 algebra is isomorphic to the ten-dimensional Jordan
2427 spin algebra, but just to be sure, we recompute its rank::
2429 sage: J = OctonionHermitianEJA(2, # long time
2430 ....: field=QQ, # long time
2431 ....: orthonormalize=False) # long time
2432 sage: J.rank.clear_cache() # long time
2433 sage: J.rank() # long time
2438 def _max_random_instance_size(max_dimension
):
2440 The maximum rank of a random OctonionHermitianEJA.
2442 # There's certainly a formula for this, but with only four
2443 # cases to worry about, I'm not that motivated to derive it.
2444 if max_dimension
>= 27:
2446 elif max_dimension
>= 10:
2448 elif max_dimension
>= 1:
2454 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2456 Return a random instance of this type of algebra.
2458 class_max_d
= cls
._max
_random
_instance
_dimension
()
2459 if (max_dimension
is None or max_dimension
> class_max_d
):
2460 max_dimension
= class_max_d
2461 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2462 n
= ZZ
.random_element(max_size
+ 1)
2463 return cls(n
, **kwargs
)
2465 def __init__(self
, n
, field
=AA
, **kwargs
):
2467 # Otherwise we don't get an EJA.
2468 raise ValueError("n cannot exceed 3")
2470 # We know this is a valid EJA, but will double-check
2471 # if the user passes check_axioms=True.
2472 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2474 from mjo
.hurwitz
import OctonionMatrixAlgebra
2475 A
= OctonionMatrixAlgebra(n
, scalars
=field
)
2476 super().__init
__(A
, **kwargs
)
2478 from mjo
.eja
.eja_cache
import octonion_hermitian_eja_coeffs
2479 a
= octonion_hermitian_eja_coeffs(self
)
2481 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2484 class AlbertEJA(OctonionHermitianEJA
):
2486 The Albert algebra is the algebra of three-by-three Hermitian
2487 matrices whose entries are octonions.
2491 sage: from mjo.eja.eja_algebra import AlbertEJA
2495 sage: AlbertEJA(field=QQ, orthonormalize=False)
2496 Euclidean Jordan algebra of dimension 27 over Rational Field
2497 sage: AlbertEJA() # long time
2498 Euclidean Jordan algebra of dimension 27 over Algebraic Real Field
2501 def __init__(self
, *args
, **kwargs
):
2502 super().__init
__(3, *args
, **kwargs
)
2505 class HadamardEJA(RationalBasisEJA
, ConcreteEJA
):
2507 Return the Euclidean Jordan algebra on `R^n` with the Hadamard
2508 (pointwise real-number multiplication) Jordan product and the
2509 usual inner-product.
2511 This is nothing more than the Cartesian product of ``n`` copies of
2512 the one-dimensional Jordan spin algebra, and is the most common
2513 example of a non-simple Euclidean Jordan algebra.
2517 sage: from mjo.eja.eja_algebra import HadamardEJA
2521 This multiplication table can be verified by hand::
2523 sage: J = HadamardEJA(3)
2524 sage: b0,b1,b2 = J.gens()
2540 We can change the generator prefix::
2542 sage: HadamardEJA(3, prefix='r').gens()
2545 def __init__(self
, n
, field
=AA
, **kwargs
):
2546 MS
= MatrixSpace(field
, n
, 1)
2549 jordan_product
= lambda x
,y
: x
2550 inner_product
= lambda x
,y
: x
2552 def jordan_product(x
,y
):
2553 return MS( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2555 def inner_product(x
,y
):
2558 # New defaults for keyword arguments. Don't orthonormalize
2559 # because our basis is already orthonormal with respect to our
2560 # inner-product. Don't check the axioms, because we know this
2561 # is a valid EJA... but do double-check if the user passes
2562 # check_axioms=True. Note: we DON'T override the "check_field"
2563 # default here, because the user can pass in a field!
2564 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2565 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2567 column_basis
= tuple( MS(b
) for b
in FreeModule(field
, n
).basis() )
2568 super().__init
__(column_basis
,
2575 self
.rank
.set_cache(n
)
2577 self
.one
.set_cache( self
.sum(self
.gens()) )
2580 def _max_random_instance_dimension():
2582 There's no reason to go higher than five here. That's
2583 enough to get the point across.
2588 def _max_random_instance_size(max_dimension
):
2590 The maximum size (=dimension) of a random HadamardEJA.
2592 return max_dimension
2595 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2597 Return a random instance of this type of algebra.
2599 class_max_d
= cls
._max
_random
_instance
_dimension
()
2600 if (max_dimension
is None or max_dimension
> class_max_d
):
2601 max_dimension
= class_max_d
2602 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2603 n
= ZZ
.random_element(max_size
+ 1)
2604 return cls(n
, **kwargs
)
2607 class BilinearFormEJA(RationalBasisEJA
, ConcreteEJA
):
2609 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2610 with the half-trace inner product and jordan product ``x*y =
2611 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2612 a symmetric positive-definite "bilinear form" matrix. Its
2613 dimension is the size of `B`, and it has rank two in dimensions
2614 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2615 the identity matrix of order ``n``.
2617 We insist that the one-by-one upper-left identity block of `B` be
2618 passed in as well so that we can be passed a matrix of size zero
2619 to construct a trivial algebra.
2623 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2624 ....: JordanSpinEJA)
2628 When no bilinear form is specified, the identity matrix is used,
2629 and the resulting algebra is the Jordan spin algebra::
2631 sage: B = matrix.identity(AA,3)
2632 sage: J0 = BilinearFormEJA(B)
2633 sage: J1 = JordanSpinEJA(3)
2634 sage: J0.multiplication_table() == J0.multiplication_table()
2637 An error is raised if the matrix `B` does not correspond to a
2638 positive-definite bilinear form::
2640 sage: B = matrix.random(QQ,2,3)
2641 sage: J = BilinearFormEJA(B)
2642 Traceback (most recent call last):
2644 ValueError: bilinear form is not positive-definite
2645 sage: B = matrix.zero(QQ,3)
2646 sage: J = BilinearFormEJA(B)
2647 Traceback (most recent call last):
2649 ValueError: bilinear form is not positive-definite
2653 We can create a zero-dimensional algebra::
2655 sage: B = matrix.identity(AA,0)
2656 sage: J = BilinearFormEJA(B)
2660 We can check the multiplication condition given in the Jordan, von
2661 Neumann, and Wigner paper (and also discussed on my "On the
2662 symmetry..." paper). Note that this relies heavily on the standard
2663 choice of basis, as does anything utilizing the bilinear form
2664 matrix. We opt not to orthonormalize the basis, because if we
2665 did, we would have to normalize the `s_{i}` in a similar manner::
2667 sage: n = ZZ.random_element(5)
2668 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2669 sage: B11 = matrix.identity(QQ,1)
2670 sage: B22 = M.transpose()*M
2671 sage: B = block_matrix(2,2,[ [B11,0 ],
2673 sage: J = BilinearFormEJA(B, orthonormalize=False)
2674 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2675 sage: V = J.vector_space()
2676 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2677 ....: for ei in eis ]
2678 sage: actual = [ sis[i]*sis[j]
2679 ....: for i in range(n-1)
2680 ....: for j in range(n-1) ]
2681 sage: expected = [ J.one() if i == j else J.zero()
2682 ....: for i in range(n-1)
2683 ....: for j in range(n-1) ]
2684 sage: actual == expected
2688 def __init__(self
, B
, field
=AA
, **kwargs
):
2689 # The matrix "B" is supplied by the user in most cases,
2690 # so it makes sense to check whether or not its positive-
2691 # definite unless we are specifically asked not to...
2692 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2693 if not B
.is_positive_definite():
2694 raise ValueError("bilinear form is not positive-definite")
2696 # However, all of the other data for this EJA is computed
2697 # by us in manner that guarantees the axioms are
2698 # satisfied. So, again, unless we are specifically asked to
2699 # verify things, we'll skip the rest of the checks.
2700 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2703 MS
= MatrixSpace(field
, n
, 1)
2705 def inner_product(x
,y
):
2706 return (y
.T
*B
*x
)[0,0]
2708 def jordan_product(x
,y
):
2713 z0
= inner_product(y
,x
)
2714 zbar
= y0
*xbar
+ x0
*ybar
2715 return MS([z0
] + zbar
.list())
2717 column_basis
= tuple( MS(b
) for b
in FreeModule(field
, n
).basis() )
2719 # TODO: I haven't actually checked this, but it seems legit.
2724 super().__init
__(column_basis
,
2729 associative
=associative
,
2732 # The rank of this algebra is two, unless we're in a
2733 # one-dimensional ambient space (because the rank is bounded
2734 # by the ambient dimension).
2735 self
.rank
.set_cache(min(n
,2))
2737 self
.one
.set_cache( self
.zero() )
2739 self
.one
.set_cache( self
.monomial(0) )
2742 def _max_random_instance_dimension():
2744 There's no reason to go higher than five here. That's
2745 enough to get the point across.
2750 def _max_random_instance_size(max_dimension
):
2752 The maximum size (=dimension) of a random BilinearFormEJA.
2754 return max_dimension
2757 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2759 Return a random instance of this algebra.
2761 class_max_d
= cls
._max
_random
_instance
_dimension
()
2762 if (max_dimension
is None or max_dimension
> class_max_d
):
2763 max_dimension
= class_max_d
2764 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2765 n
= ZZ
.random_element(max_size
+ 1)
2768 B
= matrix
.identity(ZZ
, n
)
2769 return cls(B
, **kwargs
)
2771 B11
= matrix
.identity(ZZ
, 1)
2772 M
= matrix
.random(ZZ
, n
-1)
2773 I
= matrix
.identity(ZZ
, n
-1)
2775 while alpha
.is_zero():
2776 alpha
= ZZ
.random_element().abs()
2778 B22
= M
.transpose()*M
+ alpha
*I
2780 from sage
.matrix
.special
import block_matrix
2781 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2784 return cls(B
, **kwargs
)
2787 class JordanSpinEJA(BilinearFormEJA
):
2789 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2790 with the usual inner product and jordan product ``x*y =
2791 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2796 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2800 This multiplication table can be verified by hand::
2802 sage: J = JordanSpinEJA(4)
2803 sage: b0,b1,b2,b3 = J.gens()
2819 We can change the generator prefix::
2821 sage: JordanSpinEJA(2, prefix='B').gens()
2826 Ensure that we have the usual inner product on `R^n`::
2828 sage: J = JordanSpinEJA.random_instance()
2829 sage: x,y = J.random_elements(2)
2830 sage: actual = x.inner_product(y)
2831 sage: expected = x.to_vector().inner_product(y.to_vector())
2832 sage: actual == expected
2836 def __init__(self
, n
, *args
, **kwargs
):
2837 # This is a special case of the BilinearFormEJA with the
2838 # identity matrix as its bilinear form.
2839 B
= matrix
.identity(ZZ
, n
)
2841 # Don't orthonormalize because our basis is already
2842 # orthonormal with respect to our inner-product.
2843 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2845 # But also don't pass check_field=False here, because the user
2846 # can pass in a field!
2847 super().__init
__(B
, *args
, **kwargs
)
2850 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2852 Return a random instance of this type of algebra.
2854 Needed here to override the implementation for ``BilinearFormEJA``.
2856 class_max_d
= cls
._max
_random
_instance
_dimension
()
2857 if (max_dimension
is None or max_dimension
> class_max_d
):
2858 max_dimension
= class_max_d
2859 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2860 n
= ZZ
.random_element(max_size
+ 1)
2861 return cls(n
, **kwargs
)
2864 class TrivialEJA(RationalBasisEJA
, ConcreteEJA
):
2866 The trivial Euclidean Jordan algebra consisting of only a zero element.
2870 sage: from mjo.eja.eja_algebra import TrivialEJA
2874 sage: J = TrivialEJA()
2881 sage: 7*J.one()*12*J.one()
2883 sage: J.one().inner_product(J.one())
2885 sage: J.one().norm()
2887 sage: J.one().subalgebra_generated_by()
2888 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2893 def __init__(self
, field
=AA
, **kwargs
):
2894 jordan_product
= lambda x
,y
: x
2895 inner_product
= lambda x
,y
: field
.zero()
2897 MS
= MatrixSpace(field
,0)
2899 # New defaults for keyword arguments
2900 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2901 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2903 super().__init
__(basis
,
2911 # The rank is zero using my definition, namely the dimension of the
2912 # largest subalgebra generated by any element.
2913 self
.rank
.set_cache(0)
2914 self
.one
.set_cache( self
.zero() )
2917 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2918 # We don't take a "size" argument so the superclass method is
2919 # inappropriate for us. The ``max_dimension`` argument is
2920 # included so that if this method is called generically with a
2921 # ``max_dimension=<whatever>`` argument, we don't try to pass
2922 # it on to the algebra constructor.
2923 return cls(**kwargs
)
2926 class CartesianProductEJA(FiniteDimensionalEJA
):
2928 The external (orthogonal) direct sum of two or more Euclidean
2929 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
2930 orthogonal direct sum of simple Euclidean Jordan algebras which is
2931 then isometric to a Cartesian product, so no generality is lost by
2932 providing only this construction.
2936 sage: from mjo.eja.eja_algebra import (random_eja,
2937 ....: CartesianProductEJA,
2938 ....: ComplexHermitianEJA,
2940 ....: JordanSpinEJA,
2941 ....: RealSymmetricEJA)
2945 The Jordan product is inherited from our factors and implemented by
2946 our CombinatorialFreeModule Cartesian product superclass::
2948 sage: J1 = HadamardEJA(2)
2949 sage: J2 = RealSymmetricEJA(2)
2950 sage: J = cartesian_product([J1,J2])
2951 sage: x,y = J.random_elements(2)
2955 The ability to retrieve the original factors is implemented by our
2956 CombinatorialFreeModule Cartesian product superclass::
2958 sage: J1 = HadamardEJA(2, field=QQ)
2959 sage: J2 = JordanSpinEJA(3, field=QQ)
2960 sage: J = cartesian_product([J1,J2])
2961 sage: J.cartesian_factors()
2962 (Euclidean Jordan algebra of dimension 2 over Rational Field,
2963 Euclidean Jordan algebra of dimension 3 over Rational Field)
2965 You can provide more than two factors::
2967 sage: J1 = HadamardEJA(2)
2968 sage: J2 = JordanSpinEJA(3)
2969 sage: J3 = RealSymmetricEJA(3)
2970 sage: cartesian_product([J1,J2,J3])
2971 Euclidean Jordan algebra of dimension 2 over Algebraic Real
2972 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
2973 Real Field (+) Euclidean Jordan algebra of dimension 6 over
2974 Algebraic Real Field
2976 Rank is additive on a Cartesian product::
2978 sage: J1 = HadamardEJA(1)
2979 sage: J2 = RealSymmetricEJA(2)
2980 sage: J = cartesian_product([J1,J2])
2981 sage: J1.rank.clear_cache()
2982 sage: J2.rank.clear_cache()
2983 sage: J.rank.clear_cache()
2986 sage: J.rank() == J1.rank() + J2.rank()
2989 The same rank computation works over the rationals, with whatever
2992 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
2993 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
2994 sage: J = cartesian_product([J1,J2])
2995 sage: J1.rank.clear_cache()
2996 sage: J2.rank.clear_cache()
2997 sage: J.rank.clear_cache()
3000 sage: J.rank() == J1.rank() + J2.rank()
3003 The product algebra will be associative if and only if all of its
3004 components are associative::
3006 sage: J1 = HadamardEJA(2)
3007 sage: J1.is_associative()
3009 sage: J2 = HadamardEJA(3)
3010 sage: J2.is_associative()
3012 sage: J3 = RealSymmetricEJA(3)
3013 sage: J3.is_associative()
3015 sage: CP1 = cartesian_product([J1,J2])
3016 sage: CP1.is_associative()
3018 sage: CP2 = cartesian_product([J1,J3])
3019 sage: CP2.is_associative()
3022 Cartesian products of Cartesian products work::
3024 sage: J1 = JordanSpinEJA(1)
3025 sage: J2 = JordanSpinEJA(1)
3026 sage: J3 = JordanSpinEJA(1)
3027 sage: J = cartesian_product([J1,cartesian_product([J2,J3])])
3028 sage: J.multiplication_table()
3029 +----++----+----+----+
3030 | * || b0 | b1 | b2 |
3031 +====++====+====+====+
3032 | b0 || b0 | 0 | 0 |
3033 +----++----+----+----+
3034 | b1 || 0 | b1 | 0 |
3035 +----++----+----+----+
3036 | b2 || 0 | 0 | b2 |
3037 +----++----+----+----+
3038 sage: HadamardEJA(3).multiplication_table()
3039 +----++----+----+----+
3040 | * || b0 | b1 | b2 |
3041 +====++====+====+====+
3042 | b0 || b0 | 0 | 0 |
3043 +----++----+----+----+
3044 | b1 || 0 | b1 | 0 |
3045 +----++----+----+----+
3046 | b2 || 0 | 0 | b2 |
3047 +----++----+----+----+
3049 The "matrix space" of a Cartesian product always consists of
3050 ordered pairs (or triples, or...) whose components are the
3051 matrix spaces of its factors::
3053 sage: J1 = HadamardEJA(2)
3054 sage: J2 = ComplexHermitianEJA(2)
3055 sage: J = cartesian_product([J1,J2])
3056 sage: J.matrix_space()
3057 The Cartesian product of (Full MatrixSpace of 2 by 1 dense
3058 matrices over Algebraic Real Field, Module of 2 by 2 matrices
3059 with entries in Algebraic Field over the scalar ring Algebraic
3061 sage: J.one().to_matrix()[0]
3064 sage: J.one().to_matrix()[1]
3073 All factors must share the same base field::
3075 sage: J1 = HadamardEJA(2, field=QQ)
3076 sage: J2 = RealSymmetricEJA(2)
3077 sage: CartesianProductEJA((J1,J2))
3078 Traceback (most recent call last):
3080 ValueError: all factors must share the same base field
3082 The cached unit element is the same one that would be computed::
3084 sage: J1 = random_eja() # long time
3085 sage: J2 = random_eja() # long time
3086 sage: J = cartesian_product([J1,J2]) # long time
3087 sage: actual = J.one() # long time
3088 sage: J.one.clear_cache() # long time
3089 sage: expected = J.one() # long time
3090 sage: actual == expected # long time
3093 Element
= CartesianProductEJAElement
3094 def __init__(self
, factors
, **kwargs
):
3099 self
._sets
= factors
3101 field
= factors
[0].base_ring()
3102 if not all( J
.base_ring() == field
for J
in factors
):
3103 raise ValueError("all factors must share the same base field")
3105 # Figure out the category to use.
3106 associative
= all( f
.is_associative() for f
in factors
)
3107 category
= EuclideanJordanAlgebras(field
)
3108 if associative
: category
= category
.Associative()
3109 category
= category
.join([category
, category
.CartesianProducts()])
3111 # Compute my matrix space. We don't simply use the
3112 # ``cartesian_product()`` functor here because it acts
3113 # differently on SageMath MatrixSpaces and our custom
3114 # MatrixAlgebras, which are CombinatorialFreeModules. We
3115 # always want the result to be represented (and indexed) as an
3116 # ordered tuple. This category isn't perfect, but is good
3117 # enough for what we need to do.
3118 MS_cat
= MagmaticAlgebras(field
).FiniteDimensional().WithBasis()
3119 MS_cat
= MS_cat
.Unital().CartesianProducts()
3120 MS_factors
= tuple( J
.matrix_space() for J
in factors
)
3121 from sage
.sets
.cartesian_product
import CartesianProduct
3122 self
._matrix
_space
= CartesianProduct(MS_factors
, MS_cat
)
3124 self
._matrix
_basis
= []
3125 zero
= self
._matrix
_space
.zero()
3127 for b
in factors
[i
].matrix_basis():
3130 self
._matrix
_basis
.append(z
)
3132 self
._matrix
_basis
= tuple( self
._matrix
_space
(b
)
3133 for b
in self
._matrix
_basis
)
3134 n
= len(self
._matrix
_basis
)
3136 # We already have what we need for the super-superclass constructor.
3137 CombinatorialFreeModule
.__init
__(self
,
3144 # Now create the vector space for the algebra, which will have
3145 # its own set of non-ambient coordinates (in terms of the
3147 degree
= sum( f
._matrix
_span
.ambient_vector_space().degree()
3149 V
= VectorSpace(field
, degree
)
3150 vector_basis
= tuple( V(_all2list(b
)) for b
in self
._matrix
_basis
)
3152 # Save the span of our matrix basis (when written out as long
3153 # vectors) because otherwise we'll have to reconstruct it
3154 # every time we want to coerce a matrix into the algebra.
3155 self
._matrix
_span
= V
.span_of_basis( vector_basis
, check
=False)
3157 # Since we don't (re)orthonormalize the basis, the FDEJA
3158 # constructor is going to set self._deortho_matrix to the
3159 # identity matrix. Here we set it to the correct value using
3160 # the deortho matrices from our factors.
3161 self
._deortho
_matrix
= matrix
.block_diagonal(
3162 [J
._deortho
_matrix
for J
in factors
]
3165 self
._inner
_product
_matrix
= matrix
.block_diagonal(
3166 [J
._inner
_product
_matrix
for J
in factors
]
3168 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
3169 self
._inner
_product
_matrix
.set_immutable()
3171 # Building the multiplication table is a bit more tricky
3172 # because we have to embed the entries of the factors'
3173 # multiplication tables into the product EJA.
3175 self
._multiplication
_table
= [ [zed
for j
in range(i
+1)]
3178 # Keep track of an offset that tallies the dimensions of all
3179 # previous factors. If the second factor is dim=2 and if the
3180 # first one is dim=3, then we want to skip the first 3x3 block
3181 # when copying the multiplication table for the second factor.
3184 phi_f
= self
.cartesian_embedding(f
)
3185 factor_dim
= factors
[f
].dimension()
3186 for i
in range(factor_dim
):
3187 for j
in range(i
+1):
3188 f_ij
= factors
[f
]._multiplication
_table
[i
][j
]
3190 self
._multiplication
_table
[offset
+i
][offset
+j
] = e
3191 offset
+= factor_dim
3193 self
.rank
.set_cache(sum(J
.rank() for J
in factors
))
3194 ones
= tuple(J
.one().to_matrix() for J
in factors
)
3195 self
.one
.set_cache(self(ones
))
3197 def _sets_keys(self
):
3202 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
3203 ....: RealSymmetricEJA)
3207 The superclass uses ``_sets_keys()`` to implement its
3208 ``cartesian_factors()`` method::
3210 sage: J1 = RealSymmetricEJA(2,
3212 ....: orthonormalize=False,
3214 sage: J2 = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
3215 sage: J = cartesian_product([J1,J2])
3216 sage: x = sum(i*J.gens()[i] for i in range(len(J.gens())))
3217 sage: x.cartesian_factors()
3218 (a1 + 2*a2, 3*b0 + 4*b1 + 5*b2 + 6*b3)
3221 # Copy/pasted from CombinatorialFreeModule_CartesianProduct,
3222 # but returning a tuple instead of a list.
3223 return tuple(range(len(self
.cartesian_factors())))
3225 def cartesian_factors(self
):
3226 # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
3229 def cartesian_factor(self
, i
):
3231 Return the ``i``th factor of this algebra.
3233 return self
._sets
[i
]
3236 # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
3237 from sage
.categories
.cartesian_product
import cartesian_product
3238 return cartesian_product
.symbol
.join("%s" % factor
3239 for factor
in self
._sets
)
3243 def cartesian_projection(self
, i
):
3247 sage: from mjo.eja.eja_algebra import (random_eja,
3248 ....: JordanSpinEJA,
3250 ....: RealSymmetricEJA,
3251 ....: ComplexHermitianEJA)
3255 The projection morphisms are Euclidean Jordan algebra
3258 sage: J1 = HadamardEJA(2)
3259 sage: J2 = RealSymmetricEJA(2)
3260 sage: J = cartesian_product([J1,J2])
3261 sage: J.cartesian_projection(0)
3262 Linear operator between finite-dimensional Euclidean Jordan
3263 algebras represented by the matrix:
3266 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3267 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3268 Algebraic Real Field
3269 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3271 sage: J.cartesian_projection(1)
3272 Linear operator between finite-dimensional Euclidean Jordan
3273 algebras represented by the matrix:
3277 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3278 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3279 Algebraic Real Field
3280 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3283 The projections work the way you'd expect on the vector
3284 representation of an element::
3286 sage: J1 = JordanSpinEJA(2)
3287 sage: J2 = ComplexHermitianEJA(2)
3288 sage: J = cartesian_product([J1,J2])
3289 sage: pi_left = J.cartesian_projection(0)
3290 sage: pi_right = J.cartesian_projection(1)
3291 sage: pi_left(J.one()).to_vector()
3293 sage: pi_right(J.one()).to_vector()
3295 sage: J.one().to_vector()
3300 The answer never changes::
3302 sage: J1 = random_eja()
3303 sage: J2 = random_eja()
3304 sage: J = cartesian_product([J1,J2])
3305 sage: P0 = J.cartesian_projection(0)
3306 sage: P1 = J.cartesian_projection(0)
3311 offset
= sum( self
.cartesian_factor(k
).dimension()
3313 Ji
= self
.cartesian_factor(i
)
3314 Pi
= self
._module
_morphism
(lambda j
: Ji
.monomial(j
- offset
),
3317 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3320 def cartesian_embedding(self
, i
):
3324 sage: from mjo.eja.eja_algebra import (random_eja,
3325 ....: JordanSpinEJA,
3327 ....: RealSymmetricEJA)
3331 The embedding morphisms are Euclidean Jordan algebra
3334 sage: J1 = HadamardEJA(2)
3335 sage: J2 = RealSymmetricEJA(2)
3336 sage: J = cartesian_product([J1,J2])
3337 sage: J.cartesian_embedding(0)
3338 Linear operator between finite-dimensional Euclidean Jordan
3339 algebras represented by the matrix:
3345 Domain: Euclidean Jordan algebra of dimension 2 over
3346 Algebraic Real Field
3347 Codomain: Euclidean Jordan algebra of dimension 2 over
3348 Algebraic Real Field (+) Euclidean Jordan algebra of
3349 dimension 3 over Algebraic Real Field
3350 sage: J.cartesian_embedding(1)
3351 Linear operator between finite-dimensional Euclidean Jordan
3352 algebras represented by the matrix:
3358 Domain: Euclidean Jordan algebra of dimension 3 over
3359 Algebraic Real Field
3360 Codomain: Euclidean Jordan algebra of dimension 2 over
3361 Algebraic Real Field (+) Euclidean Jordan algebra of
3362 dimension 3 over Algebraic Real Field
3364 The embeddings work the way you'd expect on the vector
3365 representation of an element::
3367 sage: J1 = JordanSpinEJA(3)
3368 sage: J2 = RealSymmetricEJA(2)
3369 sage: J = cartesian_product([J1,J2])
3370 sage: iota_left = J.cartesian_embedding(0)
3371 sage: iota_right = J.cartesian_embedding(1)
3372 sage: iota_left(J1.zero()) == J.zero()
3374 sage: iota_right(J2.zero()) == J.zero()
3376 sage: J1.one().to_vector()
3378 sage: iota_left(J1.one()).to_vector()
3380 sage: J2.one().to_vector()
3382 sage: iota_right(J2.one()).to_vector()
3384 sage: J.one().to_vector()
3389 The answer never changes::
3391 sage: J1 = random_eja()
3392 sage: J2 = random_eja()
3393 sage: J = cartesian_product([J1,J2])
3394 sage: E0 = J.cartesian_embedding(0)
3395 sage: E1 = J.cartesian_embedding(0)
3399 Composing a projection with the corresponding inclusion should
3400 produce the identity map, and mismatching them should produce
3403 sage: J1 = random_eja()
3404 sage: J2 = random_eja()
3405 sage: J = cartesian_product([J1,J2])
3406 sage: iota_left = J.cartesian_embedding(0)
3407 sage: iota_right = J.cartesian_embedding(1)
3408 sage: pi_left = J.cartesian_projection(0)
3409 sage: pi_right = J.cartesian_projection(1)
3410 sage: pi_left*iota_left == J1.one().operator()
3412 sage: pi_right*iota_right == J2.one().operator()
3414 sage: (pi_left*iota_right).is_zero()
3416 sage: (pi_right*iota_left).is_zero()
3420 offset
= sum( self
.cartesian_factor(k
).dimension()
3422 Ji
= self
.cartesian_factor(i
)
3423 Ei
= Ji
._module
_morphism
(lambda j
: self
.monomial(j
+ offset
),
3425 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3429 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3431 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3434 A separate class for products of algebras for which we know a
3439 sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA,
3441 ....: JordanSpinEJA,
3442 ....: RealSymmetricEJA)
3446 This gives us fast characteristic polynomial computations in
3447 product algebras, too::
3450 sage: J1 = JordanSpinEJA(2)
3451 sage: J2 = RealSymmetricEJA(3)
3452 sage: J = cartesian_product([J1,J2])
3453 sage: J.characteristic_polynomial_of().degree()
3460 The ``cartesian_product()`` function only uses the first factor to
3461 decide where the result will live; thus we have to be careful to
3462 check that all factors do indeed have a ``rational_algebra()`` method
3463 before we construct an algebra that claims to have a rational basis::
3465 sage: J1 = HadamardEJA(2)
3466 sage: jp = lambda X,Y: X*Y
3467 sage: ip = lambda X,Y: X[0,0]*Y[0,0]
3468 sage: b1 = matrix(QQ, [[1]])
3469 sage: J2 = FiniteDimensionalEJA((b1,), jp, ip)
3470 sage: cartesian_product([J2,J1]) # factor one not RationalBasisEJA
3471 Euclidean Jordan algebra of dimension 1 over Algebraic Real
3472 Field (+) Euclidean Jordan algebra of dimension 2 over Algebraic
3474 sage: cartesian_product([J1,J2]) # factor one is RationalBasisEJA
3475 Traceback (most recent call last):
3477 ValueError: factor not a RationalBasisEJA
3480 def __init__(self
, algebras
, **kwargs
):
3481 if not all( hasattr(r
, "rational_algebra") for r
in algebras
):
3482 raise ValueError("factor not a RationalBasisEJA")
3484 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3487 def rational_algebra(self
):
3488 if self
.base_ring() is QQ
:
3491 return cartesian_product([
3492 r
.rational_algebra() for r
in self
.cartesian_factors()
3496 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3498 def random_eja(max_dimension
=None, *args
, **kwargs
):
3503 sage: from mjo.eja.eja_algebra import random_eja
3507 sage: n = ZZ.random_element(1,5)
3508 sage: J = random_eja(max_dimension=n, field=QQ, orthonormalize=False)
3509 sage: J.dimension() <= n
3513 # Use the ConcreteEJA default as the total upper bound (regardless
3514 # of any whether or not any individual factors set a lower limit).
3515 if max_dimension
is None:
3516 max_dimension
= ConcreteEJA
._max
_random
_instance
_dimension
()
3517 J1
= ConcreteEJA
.random_instance(max_dimension
, *args
, **kwargs
)
3520 # Roll the dice to see if we attempt a Cartesian product.
3521 dice_roll
= ZZ
.random_element(len(ConcreteEJA
.__subclasses
__()) + 1)
3522 new_max_dimension
= max_dimension
- J1
.dimension()
3523 if new_max_dimension
== 0 or dice_roll
!= 0:
3524 # If it's already as big as we're willing to tolerate, just
3525 # return it and don't worry about Cartesian products.
3528 # Use random_eja() again so we can get more than two factors
3529 # if the sub-call also Decides on a cartesian product.
3530 J2
= random_eja(new_max_dimension
, *args
, **kwargs
)
3531 return cartesian_product([J1
,J2
])
3534 class ComplexSkewSymmetricEJA(RationalBasisEJA
, ConcreteEJA
):
3536 The skew-symmetric EJA of order `2n` described in Faraut and
3537 Koranyi's Exercise III.1.b. It has dimension `2n^2 - n`.
3539 It is (not obviously) isomorphic to the QuaternionHermitianEJA of
3540 order `n`, as can be inferred by comparing rank/dimension or
3541 explicitly from their "characteristic polynomial of" functions,
3542 which just so happen to align nicely.
3546 sage: from mjo.eja.eja_algebra import (ComplexSkewSymmetricEJA,
3547 ....: QuaternionHermitianEJA)
3548 sage: from mjo.eja.eja_operator import FiniteDimensionalEJAOperator
3552 This EJA is isomorphic to the quaternions::
3554 sage: J = ComplexSkewSymmetricEJA(2, field=QQ, orthonormalize=False)
3555 sage: K = QuaternionHermitianEJA(2, field=QQ, orthonormalize=False)
3556 sage: jordan_isom_matrix = matrix.diagonal(QQ,[-1,1,1,1,1,-1])
3557 sage: phi = FiniteDimensionalEJAOperator(J,K,jordan_isom_matrix)
3558 sage: all( phi(x*y) == phi(x)*phi(y)
3559 ....: for x in J.gens()
3560 ....: for y in J.gens() )
3562 sage: x,y = J.random_elements(2)
3563 sage: phi(x*y) == phi(x)*phi(y)
3568 Random elements should satisfy the same conditions that the basis
3571 sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ,
3572 ....: orthonormalize=False)
3573 sage: x,y = K.random_elements(2)
3575 sage: x = x.to_matrix()
3576 sage: y = y.to_matrix()
3577 sage: z = z.to_matrix()
3578 sage: all( e.is_skew_symmetric() for e in (x,y,z) )
3580 sage: J = -K.one().to_matrix()
3581 sage: all( e*J == J*e.conjugate() for e in (x,y,z) )
3584 The power law in Faraut & Koranyi's II.7.a is satisfied.
3585 We're in a subalgebra of theirs, but powers are still
3588 sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ,
3589 ....: orthonormalize=False)
3590 sage: x = K.random_element()
3591 sage: k = ZZ.random_element(5)
3593 sage: J = -K.one().to_matrix()
3594 sage: expected = K(-J*(J*x.to_matrix())^k)
3595 sage: actual == expected
3600 def _max_random_instance_size(max_dimension
):
3601 # Obtained by solving d = 2n^2 - n, which comes from noticing
3602 # that, in 2x2 block form, any element of this algebra has a
3603 # free skew-symmetric top-left block, a Hermitian top-right
3604 # block, and two bottom blocks that are determined by the top.
3605 # The ZZ-int-ZZ thing is just "floor."
3606 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/4 + 1/4))
3609 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
3611 Return a random instance of this type of algebra.
3613 class_max_d
= cls
._max
_random
_instance
_dimension
()
3614 if (max_dimension
is None or max_dimension
> class_max_d
):
3615 max_dimension
= class_max_d
3616 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
3617 n
= ZZ
.random_element(max_size
+ 1)
3618 return cls(n
, **kwargs
)
3621 def _denormalized_basis(A
):
3625 sage: from mjo.hurwitz import ComplexMatrixAlgebra
3626 sage: from mjo.eja.eja_algebra import ComplexSkewSymmetricEJA
3630 The basis elements are all skew-Hermitian::
3632 sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension()
3633 sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max)
3634 sage: n = ZZ.random_element(n_max + 1)
3635 sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ)
3636 sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A)
3637 sage: all( M.is_skew_symmetric() for M in B)
3640 The basis elements ``b`` all satisfy ``b*J == J*b.conjugate()``,
3641 as in the definition of the algebra::
3643 sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension()
3644 sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max)
3645 sage: n = ZZ.random_element(n_max + 1)
3646 sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ)
3647 sage: I_n = matrix.identity(ZZ, n)
3648 sage: J = matrix.block(ZZ, 2, 2, (0, I_n, -I_n, 0), subdivide=False)
3649 sage: J = A.from_list(J.rows())
3650 sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A)
3651 sage: all( b*J == J*b.conjugate() for b in B )
3655 es
= A
.entry_algebra_gens()
3656 gen
= lambda A
,m
: A
.monomial(m
)
3660 # The size of the blocks. We're going to treat these thing as
3661 # 2x2 block matrices,
3664 # [ -x2-conj x1-conj ]
3666 # where x1 is skew-symmetric and x2 is Hermitian.
3670 # We only loop through the top half of the matrix, because the
3671 # bottom can be constructed from the top.
3673 # First do the top-left block, which is skew-symmetric.
3674 # We can compute the bottom-right block in the process.
3675 for j
in range(i
+1):
3677 # Skew-symmetry implies zeros for (i == j).
3679 # Top-left block's entry.
3680 E_ij
= gen(A
, (i
,j
,e
))
3681 E_ij
-= gen(A
, (j
,i
,e
))
3683 # Bottom-right block's entry.
3684 F_ij
= gen(A
, (i
+m
,j
+m
,e
)).conjugate()
3685 F_ij
-= gen(A
, (j
+m
,i
+m
,e
)).conjugate()
3687 basis
.append(E_ij
+ F_ij
)
3689 # Now do the top-right block, which is Hermitian, and compute
3690 # the bottom-left block along the way.
3691 for j
in range(m
,i
+m
+1):
3693 # Hermitian matrices have real diagonal entries.
3694 # Top-right block's entry.
3695 E_ii
= gen(A
, (i
,j
,es
[0]))
3697 # Bottom-left block's entry. Don't conjugate
3699 E_ii
-= gen(A
, (i
+m
,j
-m
,es
[0]))
3703 # Top-right block's entry. BEWARE! We're not
3704 # reflecting across the main diagonal as in
3705 # (i,j)~(j,i). We're only reflecting across
3706 # the diagonal for the top-right block.
3707 E_ij
= gen(A
, (i
,j
,e
))
3709 # Shift it back to non-offset coords, transpose,
3710 # conjugate, and put it back:
3712 # (i,j) -> (i,j-m) -> (j-m, i) -> (j-m, i+m)
3713 E_ij
+= gen(A
, (j
-m
,i
+m
,e
)).conjugate()
3715 # Bottom-left's block's below-diagonal entry.
3716 # Just shift the top-right coords down m and
3718 F_ij
= -gen(A
, (i
+m
,j
-m
,e
)).conjugate()
3719 F_ij
+= -gen(A
, (j
,i
,e
)) # double-conjugate cancels
3721 basis
.append(E_ij
+ F_ij
)
3723 return tuple( basis
)
3727 def _J_matrix(matrix_space
):
3728 n
= matrix_space
.nrows() // 2
3729 F
= matrix_space
.base_ring()
3730 I_n
= matrix
.identity(F
, n
)
3731 J
= matrix
.block(F
, 2, 2, (0, I_n
, -I_n
, 0), subdivide
=False)
3732 return matrix_space
.from_list(J
.rows())
3735 return ComplexSkewSymmetricEJA
._J
_matrix
(self
.matrix_space())
3737 def __init__(self
, n
, field
=AA
, **kwargs
):
3738 # New code; always check the axioms.
3739 #if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
3741 from mjo
.hurwitz
import ComplexMatrixAlgebra
3742 A
= ComplexMatrixAlgebra(2*n
, scalars
=field
)
3743 J
= ComplexSkewSymmetricEJA
._J
_matrix
(A
)
3745 def jordan_product(X
,Y
):
3746 return (X
*J
*Y
+ Y
*J
*X
)/2
3748 def inner_product(X
,Y
):
3749 return (X
.conjugate_transpose()*Y
).trace().real()
3751 super().__init
__(self
._denormalized
_basis
(A
),
3758 # This algebra is conjectured (by me) to be isomorphic to
3759 # the quaternion Hermitian EJA of size n, and the rank
3760 # would follow from that.
3761 #self.rank.set_cache(n)
3762 self
.one
.set_cache( self(-J
) )