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 a few other example constructions,
39 * :class:`JordanSpinEJA`
40 * :class:`HadamardEJA`
43 * :class:`ComplexSkewSymmetricEJA`
45 The Jordan spin algebra is a bilinear form algebra where the bilinear
46 form is the identity. The Hadamard EJA is simply a Cartesian product
47 of one-dimensional spin algebras. The Albert EJA is simply a special
48 case of the :class:`OctonionHermitianEJA` where the matrices are
49 three-by-three and the resulting space has dimension 27. And
50 last/least, the trivial EJA is exactly what you think it is; it could
51 also be obtained by constructing a dimension-zero instance of any of
52 the other algebras. Cartesian products of these are also supported
53 using the usual ``cartesian_product()`` function; as a result, we
54 support (up to isomorphism) all Euclidean Jordan algebras.
56 At a minimum, the following are required to construct a Euclidean
59 * A basis of matrices, column vectors, or MatrixAlgebra elements
60 * A Jordan product defined on the basis
61 * Its inner product defined on the basis
63 The real numbers form a Euclidean Jordan algebra when both the Jordan
64 and inner products are the usual multiplication. We use this as our
65 example, and demonstrate a few ways to construct an EJA.
67 First, we can use one-by-one SageMath matrices with algebraic real
68 entries to represent real numbers. We define the Jordan and inner
69 products to be essentially real-number multiplication, with the only
70 difference being that the Jordan product again returns a one-by-one
71 matrix, whereas the inner product must return a scalar. Our basis for
72 the one-by-one matrices is of course the set consisting of a single
73 matrix with its sole entry non-zero::
75 sage: from mjo.eja.eja_algebra import FiniteDimensionalEJA
76 sage: jp = lambda X,Y: X*Y
77 sage: ip = lambda X,Y: X[0,0]*Y[0,0]
78 sage: b1 = matrix(AA, [[1]])
79 sage: J1 = FiniteDimensionalEJA((b1,), jp, ip)
81 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
83 In fact, any positive scalar multiple of that inner-product would work::
85 sage: ip2 = lambda X,Y: 16*ip(X,Y)
86 sage: J2 = FiniteDimensionalEJA((b1,), jp, ip2)
88 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
90 But beware that your basis will be orthonormalized _with respect to the
91 given inner-product_ unless you pass ``orthonormalize=False`` to the
92 constructor. For example::
94 sage: J3 = FiniteDimensionalEJA((b1,), jp, ip2, orthonormalize=False)
96 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
98 To see the difference, you can take the first and only basis element
99 of the resulting algebra, and ask for it to be converted back into
102 sage: J1.basis()[0].to_matrix()
104 sage: J2.basis()[0].to_matrix()
106 sage: J3.basis()[0].to_matrix()
109 Since square roots are used in that process, the default scalar field
110 that we use is the field of algebraic real numbers, ``AA``. You can
111 also Use rational numbers, but only if you either pass
112 ``orthonormalize=False`` or know that orthonormalizing your basis
113 won't stray beyond the rational numbers. The example above would
114 have worked only because ``sqrt(16) == 4`` is rational.
116 Another option for your basis is to use elemebts of a
117 :class:`MatrixAlgebra`::
119 sage: from mjo.matrix_algebra import MatrixAlgebra
120 sage: A = MatrixAlgebra(1,AA,AA)
121 sage: J4 = FiniteDimensionalEJA(A.gens(), jp, ip)
123 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
124 sage: J4.basis()[0].to_matrix()
129 An easier way to view the entire EJA basis in its original (but
130 perhaps orthonormalized) matrix form is to use the ``matrix_basis``
133 sage: J4.matrix_basis()
138 In particular, a :class:`MatrixAlgebra` is needed to work around the
139 fact that matrices in SageMath must have entries in the same
140 (commutative and associative) ring as its scalars. There are many
141 Euclidean Jordan algebras whose elements are matrices that violate
142 those assumptions. The complex, quaternion, and octonion Hermitian
143 matrices all have entries in a ring (the complex numbers, quaternions,
144 or octonions...) that differs from the algebra's scalar ring (the real
145 numbers). Quaternions are also non-commutative; the octonions are
146 neither commutative nor associative.
150 sage: from mjo.eja.eja_algebra import random_eja
155 Euclidean Jordan algebra of dimension...
158 from sage
.algebras
.quatalg
.quaternion_algebra
import QuaternionAlgebra
159 from sage
.categories
.magmatic_algebras
import MagmaticAlgebras
160 from sage
.categories
.sets_cat
import cartesian_product
161 from sage
.combinat
.free_module
import CombinatorialFreeModule
162 from sage
.matrix
.constructor
import matrix
163 from sage
.matrix
.matrix_space
import MatrixSpace
164 from sage
.misc
.cachefunc
import cached_method
165 from sage
.misc
.table
import table
166 from sage
.modules
.free_module
import FreeModule
, VectorSpace
167 from sage
.rings
.all
import (ZZ
, QQ
, AA
, QQbar
, RR
, RLF
, CLF
,
170 from mjo
.eja
.eja_element
import (CartesianProductEJAElement
,
171 FiniteDimensionalEJAElement
)
172 from mjo
.eja
.eja_operator
import FiniteDimensionalEJAOperator
173 from mjo
.eja
.eja_utils
import _all2list
175 def EuclideanJordanAlgebras(field
):
177 The category of Euclidean Jordan algebras over ``field``, which
178 must be a subfield of the real numbers. For now this is just a
179 convenient wrapper around all of the other category axioms that
182 category
= MagmaticAlgebras(field
).FiniteDimensional()
183 category
= category
.WithBasis().Unital().Commutative()
186 class FiniteDimensionalEJA(CombinatorialFreeModule
):
188 A finite-dimensional Euclidean Jordan algebra.
192 - ``basis`` -- a tuple; a tuple of basis elements in "matrix
193 form," which must be the same form as the arguments to
194 ``jordan_product`` and ``inner_product``. In reality, "matrix
195 form" can be either vectors, matrices, or a Cartesian product
196 (ordered tuple) of vectors or matrices. All of these would
197 ideally be vector spaces in sage with no special-casing
198 needed; but in reality we turn vectors into column-matrices
199 and Cartesian products `(a,b)` into column matrices
200 `(a,b)^{T}` after converting `a` and `b` themselves.
202 - ``jordan_product`` -- a function; afunction of two ``basis``
203 elements (in matrix form) that returns their jordan product,
204 also in matrix form; this will be applied to ``basis`` to
205 compute a multiplication table for the algebra.
207 - ``inner_product`` -- a function; a function of two ``basis``
208 elements (in matrix form) that returns their inner
209 product. This will be applied to ``basis`` to compute an
210 inner-product table (basically a matrix) for this algebra.
212 - ``matrix_space`` -- the space that your matrix basis lives in,
213 or ``None`` (the default). So long as your basis does not have
214 length zero you can omit this. But in trivial algebras, it is
217 - ``field`` -- a subfield of the reals (default: ``AA``); the scalar
218 field for the algebra.
220 - ``orthonormalize`` -- boolean (default: ``True``); whether or
221 not to orthonormalize the basis. Doing so is expensive and
222 generally rules out using the rationals as your ``field``, but
223 is required for spectral decompositions.
227 sage: from mjo.eja.eja_algebra import random_eja
231 We should compute that an element subalgebra is associative even
232 if we circumvent the element method::
234 sage: J = random_eja(field=QQ,orthonormalize=False)
235 sage: x = J.random_element()
236 sage: A = x.subalgebra_generated_by(orthonormalize=False)
237 sage: basis = tuple(b.superalgebra_element() for b in A.basis())
238 sage: J.subalgebra(basis, orthonormalize=False).is_associative()
241 Element
= FiniteDimensionalEJAElement
244 def _check_input_field(field
):
245 if not field
.is_subring(RR
):
246 # Note: this does return true for the real algebraic
247 # field, the rationals, and any quadratic field where
248 # we've specified a real embedding.
249 raise ValueError("scalar field is not real")
252 def _check_input_axioms(basis
, jordan_product
, inner_product
):
253 if not all( jordan_product(bi
,bj
) == jordan_product(bj
,bi
)
256 raise ValueError("Jordan product is not commutative")
258 if not all( inner_product(bi
,bj
) == inner_product(bj
,bi
)
261 raise ValueError("inner-product is not commutative")
278 self
._check
_input
_field
(field
)
281 # Check commutativity of the Jordan and inner-products.
282 # This has to be done before we build the multiplication
283 # and inner-product tables/matrices, because we take
284 # advantage of symmetry in the process.
285 self
._check
_input
_axioms
(basis
, jordan_product
, inner_product
)
288 # All zero- and one-dimensional algebras are just the real
289 # numbers with (some positive multiples of) the usual
290 # multiplication as its Jordan and inner-product.
292 if associative
is None:
293 # We should figure it out. As with check_axioms, we have to do
294 # this without the help of the _jordan_product_is_associative()
295 # method because we need to know the category before we
296 # initialize the algebra.
297 associative
= all( jordan_product(jordan_product(bi
,bj
),bk
)
299 jordan_product(bi
,jordan_product(bj
,bk
))
304 category
= EuclideanJordanAlgebras(field
)
307 # Element subalgebras can take advantage of this.
308 category
= category
.Associative()
310 # Call the superclass constructor so that we can use its from_vector()
311 # method to build our multiplication table.
312 CombinatorialFreeModule
.__init
__(self
,
319 # Now comes all of the hard work. We'll be constructing an
320 # ambient vector space V that our (vectorized) basis lives in,
321 # as well as a subspace W of V spanned by those (vectorized)
322 # basis elements. The W-coordinates are the coefficients that
323 # we see in things like x = 1*b1 + 2*b2.
327 degree
= len(_all2list(basis
[0]))
329 # Build an ambient space that fits our matrix basis when
330 # written out as "long vectors."
331 V
= VectorSpace(field
, degree
)
333 # The matrix that will hold the orthonormal -> unorthonormal
334 # coordinate transformation. Default to an identity matrix of
335 # the appropriate size to avoid special cases for None
337 self
._deortho
_matrix
= matrix
.identity(field
,n
)
340 # Save a copy of the un-orthonormalized basis for later.
341 # Convert it to ambient V (vector) coordinates while we're
342 # at it, because we'd have to do it later anyway.
343 deortho_vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
345 from mjo
.eja
.eja_utils
import gram_schmidt
346 basis
= tuple(gram_schmidt(basis
, inner_product
))
348 # Save the (possibly orthonormalized) matrix basis for
349 # later, as well as the space that its elements live in.
350 # In most cases we can deduce the matrix space, but when
351 # n == 0 (that is, there are no basis elements) we cannot.
352 self
._matrix
_basis
= basis
353 if matrix_space
is None:
354 self
._matrix
_space
= self
._matrix
_basis
[0].parent()
356 self
._matrix
_space
= matrix_space
358 # Now create the vector space for the algebra, which will have
359 # its own set of non-ambient coordinates (in terms of the
361 vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
363 # Save the span of our matrix basis (when written out as long
364 # vectors) because otherwise we'll have to reconstruct it
365 # every time we want to coerce a matrix into the algebra.
366 self
._matrix
_span
= V
.span_of_basis( vector_basis
, check
=check_axioms
)
369 # Now "self._matrix_span" is the vector space of our
370 # algebra coordinates. The variables "X0", "X1",... refer
371 # to the entries of vectors in self._matrix_span. Thus to
372 # convert back and forth between the orthonormal
373 # coordinates and the given ones, we need to stick the
374 # original basis in self._matrix_span.
375 U
= V
.span_of_basis( deortho_vector_basis
, check
=check_axioms
)
376 self
._deortho
_matrix
= matrix
.column( U
.coordinate_vector(q
)
377 for q
in vector_basis
)
380 # Now we actually compute the multiplication and inner-product
381 # tables/matrices using the possibly-orthonormalized basis.
382 self
._inner
_product
_matrix
= matrix
.identity(field
, n
)
384 self
._multiplication
_table
= [ [zed
for j
in range(i
+1)]
387 # Note: the Jordan and inner-products are defined in terms
388 # of the ambient basis. It's important that their arguments
389 # are in ambient coordinates as well.
392 # ortho basis w.r.t. ambient coords
396 # The jordan product returns a matrixy answer, so we
397 # have to convert it to the algebra coordinates.
398 elt
= jordan_product(q_i
, q_j
)
399 elt
= self
._matrix
_span
.coordinate_vector(V(_all2list(elt
)))
400 self
._multiplication
_table
[i
][j
] = self
.from_vector(elt
)
402 if not orthonormalize
:
403 # If we're orthonormalizing the basis with respect
404 # to an inner-product, then the inner-product
405 # matrix with respect to the resulting basis is
406 # just going to be the identity.
407 ip
= inner_product(q_i
, q_j
)
408 self
._inner
_product
_matrix
[i
,j
] = ip
409 self
._inner
_product
_matrix
[j
,i
] = ip
411 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
412 self
._inner
_product
_matrix
.set_immutable()
415 if not self
._is
_jordanian
():
416 raise ValueError("Jordan identity does not hold")
417 if not self
._inner
_product
_is
_associative
():
418 raise ValueError("inner product is not associative")
421 def _coerce_map_from_base_ring(self
):
423 Disable the map from the base ring into the algebra.
425 Performing a nonsense conversion like this automatically
426 is counterpedagogical. The fallback is to try the usual
427 element constructor, which should also fail.
431 sage: from mjo.eja.eja_algebra import random_eja
435 sage: J = random_eja()
437 Traceback (most recent call last):
439 ValueError: not an element of this algebra
445 def product_on_basis(self
, i
, j
):
447 Returns the Jordan product of the `i` and `j`th basis elements.
449 This completely defines the Jordan product on the algebra, and
450 is used direclty by our superclass machinery to implement
455 sage: from mjo.eja.eja_algebra import random_eja
459 sage: J = random_eja()
460 sage: n = J.dimension()
463 sage: bi_bj = J.zero()*J.zero()
465 ....: i = ZZ.random_element(n)
466 ....: j = ZZ.random_element(n)
467 ....: bi = J.monomial(i)
468 ....: bj = J.monomial(j)
469 ....: bi_bj = J.product_on_basis(i,j)
474 # We only stored the lower-triangular portion of the
475 # multiplication table.
477 return self
._multiplication
_table
[i
][j
]
479 return self
._multiplication
_table
[j
][i
]
481 def inner_product(self
, x
, y
):
483 The inner product associated with this Euclidean Jordan algebra.
485 Defaults to the trace inner product, but can be overridden by
486 subclasses if they are sure that the necessary properties are
491 sage: from mjo.eja.eja_algebra import (random_eja,
493 ....: BilinearFormEJA)
497 Our inner product is "associative," which means the following for
498 a symmetric bilinear form::
500 sage: J = random_eja()
501 sage: x,y,z = J.random_elements(3)
502 sage: (x*y).inner_product(z) == y.inner_product(x*z)
507 Ensure that this is the usual inner product for the algebras
510 sage: J = HadamardEJA.random_instance()
511 sage: x,y = J.random_elements(2)
512 sage: actual = x.inner_product(y)
513 sage: expected = x.to_vector().inner_product(y.to_vector())
514 sage: actual == expected
517 Ensure that this is one-half of the trace inner-product in a
518 BilinearFormEJA that isn't just the reals (when ``n`` isn't
519 one). This is in Faraut and Koranyi, and also my "On the
522 sage: J = BilinearFormEJA.random_instance()
523 sage: n = J.dimension()
524 sage: x = J.random_element()
525 sage: y = J.random_element()
526 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
530 B
= self
._inner
_product
_matrix
531 return (B
*x
.to_vector()).inner_product(y
.to_vector())
534 def is_associative(self
):
536 Return whether or not this algebra's Jordan product is associative.
540 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
544 sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False)
545 sage: J.is_associative()
547 sage: x = sum(J.gens())
548 sage: A = x.subalgebra_generated_by(orthonormalize=False)
549 sage: A.is_associative()
553 return "Associative" in self
.category().axioms()
555 def _is_commutative(self
):
557 Whether or not this algebra's multiplication table is commutative.
559 This method should of course always return ``True``, unless
560 this algebra was constructed with ``check_axioms=False`` and
561 passed an invalid multiplication table.
563 return all( x
*y
== y
*x
for x
in self
.gens() for y
in self
.gens() )
565 def _is_jordanian(self
):
567 Whether or not this algebra's multiplication table respects the
568 Jordan identity `(x^{2})(xy) = x(x^{2}y)`.
570 We only check one arrangement of `x` and `y`, so for a
571 ``True`` result to be truly true, you should also check
572 :meth:`_is_commutative`. This method should of course always
573 return ``True``, unless this algebra was constructed with
574 ``check_axioms=False`` and passed an invalid multiplication table.
576 return all( (self
.monomial(i
)**2)*(self
.monomial(i
)*self
.monomial(j
))
578 (self
.monomial(i
))*((self
.monomial(i
)**2)*self
.monomial(j
))
579 for i
in range(self
.dimension())
580 for j
in range(self
.dimension()) )
582 def _jordan_product_is_associative(self
):
584 Return whether or not this algebra's Jordan product is
585 associative; that is, whether or not `x*(y*z) = (x*y)*z`
588 This method should agree with :meth:`is_associative` unless
589 you lied about the value of the ``associative`` parameter
590 when you constructed the algebra.
594 sage: from mjo.eja.eja_algebra import (random_eja,
595 ....: RealSymmetricEJA,
596 ....: ComplexHermitianEJA,
597 ....: QuaternionHermitianEJA)
601 sage: J = RealSymmetricEJA(4, orthonormalize=False)
602 sage: J._jordan_product_is_associative()
604 sage: x = sum(J.gens())
605 sage: A = x.subalgebra_generated_by()
606 sage: A._jordan_product_is_associative()
611 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
612 sage: J._jordan_product_is_associative()
614 sage: x = sum(J.gens())
615 sage: A = x.subalgebra_generated_by(orthonormalize=False)
616 sage: A._jordan_product_is_associative()
621 sage: J = QuaternionHermitianEJA(2)
622 sage: J._jordan_product_is_associative()
624 sage: x = sum(J.gens())
625 sage: A = x.subalgebra_generated_by()
626 sage: A._jordan_product_is_associative()
631 The values we've presupplied to the constructors agree with
634 sage: J = random_eja()
635 sage: J.is_associative() == J._jordan_product_is_associative()
641 # Used to check whether or not something is zero.
644 # I don't know of any examples that make this magnitude
645 # necessary because I don't know how to make an
646 # associative algebra when the element subalgebra
647 # construction is unreliable (as it is over RDF; we can't
648 # find the degree of an element because we can't compute
649 # the rank of a matrix). But even multiplication of floats
650 # is non-associative, so *some* epsilon is needed... let's
651 # just take the one from _inner_product_is_associative?
654 for i
in range(self
.dimension()):
655 for j
in range(self
.dimension()):
656 for k
in range(self
.dimension()):
660 diff
= (x
*y
)*z
- x
*(y
*z
)
662 if diff
.norm() > epsilon
:
667 def _inner_product_is_associative(self
):
669 Return whether or not this algebra's inner product `B` is
670 associative; that is, whether or not `B(xy,z) = B(x,yz)`.
672 This method should of course always return ``True``, unless
673 this algebra was constructed with ``check_axioms=False`` and
674 passed an invalid Jordan or inner-product.
678 # Used to check whether or not something is zero.
681 # This choice is sufficient to allow the construction of
682 # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
685 for i
in range(self
.dimension()):
686 for j
in range(self
.dimension()):
687 for k
in range(self
.dimension()):
691 diff
= (x
*y
).inner_product(z
) - x
.inner_product(y
*z
)
693 if diff
.abs() > epsilon
:
698 def _element_constructor_(self
, elt
):
700 Construct an element of this algebra or a subalgebra from its
701 EJA element, vector, or matrix representation.
703 This gets called only after the parent element _call_ method
704 fails to find a coercion for the argument.
708 sage: from mjo.eja.eja_algebra import (random_eja,
711 ....: RealSymmetricEJA)
715 The identity in `S^n` is converted to the identity in the EJA::
717 sage: J = RealSymmetricEJA(3)
718 sage: I = matrix.identity(QQ,3)
719 sage: J(I) == J.one()
722 This skew-symmetric matrix can't be represented in the EJA::
724 sage: J = RealSymmetricEJA(3)
725 sage: A = matrix(QQ,3, lambda i,j: i-j)
727 Traceback (most recent call last):
729 ValueError: not an element of this algebra
731 Tuples work as well, provided that the matrix basis for the
732 algebra consists of them::
734 sage: J1 = HadamardEJA(3)
735 sage: J2 = RealSymmetricEJA(2)
736 sage: J = cartesian_product([J1,J2])
737 sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
740 Subalgebra elements are embedded into the superalgebra::
742 sage: J = JordanSpinEJA(3)
745 sage: x = sum(J.gens())
746 sage: A = x.subalgebra_generated_by()
752 Ensure that we can convert any element back and forth
753 faithfully between its matrix and algebra representations::
755 sage: J = random_eja()
756 sage: x = J.random_element()
757 sage: J(x.to_matrix()) == x
760 We cannot coerce elements between algebras just because their
761 matrix representations are compatible::
763 sage: J1 = HadamardEJA(3)
764 sage: J2 = JordanSpinEJA(3)
766 Traceback (most recent call last):
768 ValueError: not an element of this algebra
770 Traceback (most recent call last):
772 ValueError: not an element of this algebra
775 msg
= "not an element of this algebra"
776 if elt
in self
.base_ring():
777 # Ensure that no base ring -> algebra coercion is performed
778 # by this method. There's some stupidity in sage that would
779 # otherwise propagate to this method; for example, sage thinks
780 # that the integer 3 belongs to the space of 2-by-2 matrices.
781 raise ValueError(msg
)
783 if hasattr(elt
, 'superalgebra_element'):
784 # Handle subalgebra elements
785 if elt
.parent().superalgebra() == self
:
786 return elt
.superalgebra_element()
788 if hasattr(elt
, 'sparse_vector'):
789 # Convert a vector into a column-matrix. We check for
790 # "sparse_vector" and not "column" because matrices also
791 # have a "column" method.
794 if elt
not in self
.matrix_space():
795 raise ValueError(msg
)
797 # Thanks for nothing! Matrix spaces aren't vector spaces in
798 # Sage, so we have to figure out its matrix-basis coordinates
799 # ourselves. We use the basis space's ring instead of the
800 # element's ring because the basis space might be an algebraic
801 # closure whereas the base ring of the 3-by-3 identity matrix
802 # could be QQ instead of QQbar.
804 # And, we also have to handle Cartesian product bases (when
805 # the matrix basis consists of tuples) here. The "good news"
806 # is that we're already converting everything to long vectors,
807 # and that strategy works for tuples as well.
809 elt
= self
._matrix
_span
.ambient_vector_space()(_all2list(elt
))
812 coords
= self
._matrix
_span
.coordinate_vector(elt
)
813 except ArithmeticError: # vector is not in free module
814 raise ValueError(msg
)
816 return self
.from_vector(coords
)
820 Return a string representation of ``self``.
824 sage: from mjo.eja.eja_algebra import JordanSpinEJA
828 Ensure that it says what we think it says::
830 sage: JordanSpinEJA(2, field=AA)
831 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
832 sage: JordanSpinEJA(3, field=RDF)
833 Euclidean Jordan algebra of dimension 3 over Real Double Field
836 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
837 return fmt
.format(self
.dimension(), self
.base_ring())
841 def characteristic_polynomial_of(self
):
843 Return the algebra's "characteristic polynomial of" function,
844 which is itself a multivariate polynomial that, when evaluated
845 at the coordinates of some algebra element, returns that
846 element's characteristic polynomial.
848 The resulting polynomial has `n+1` variables, where `n` is the
849 dimension of this algebra. The first `n` variables correspond to
850 the coordinates of an algebra element: when evaluated at the
851 coordinates of an algebra element with respect to a certain
852 basis, the result is a univariate polynomial (in the one
853 remaining variable ``t``), namely the characteristic polynomial
858 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
862 The characteristic polynomial in the spin algebra is given in
863 Alizadeh, Example 11.11::
865 sage: J = JordanSpinEJA(3)
866 sage: p = J.characteristic_polynomial_of(); p
867 X0^2 - X1^2 - X2^2 + (-2*t)*X0 + t^2
868 sage: xvec = J.one().to_vector()
872 By definition, the characteristic polynomial is a monic
873 degree-zero polynomial in a rank-zero algebra. Note that
874 Cayley-Hamilton is indeed satisfied since the polynomial
875 ``1`` evaluates to the identity element of the algebra on
878 sage: J = TrivialEJA()
879 sage: J.characteristic_polynomial_of()
886 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
887 a
= self
._charpoly
_coefficients
()
889 # We go to a bit of trouble here to reorder the
890 # indeterminates, so that it's easier to evaluate the
891 # characteristic polynomial at x's coordinates and get back
892 # something in terms of t, which is what we want.
893 S
= PolynomialRing(self
.base_ring(),'t')
897 S
= PolynomialRing(S
, R
.variable_names())
900 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
902 def coordinate_polynomial_ring(self
):
904 The multivariate polynomial ring in which this algebra's
905 :meth:`characteristic_polynomial_of` lives.
909 sage: from mjo.eja.eja_algebra import (HadamardEJA,
910 ....: RealSymmetricEJA)
914 sage: J = HadamardEJA(2)
915 sage: J.coordinate_polynomial_ring()
916 Multivariate Polynomial Ring in X0, X1...
917 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
918 sage: J.coordinate_polynomial_ring()
919 Multivariate Polynomial Ring in X0, X1, X2, X3, X4, X5...
922 var_names
= tuple( "X%d" % z
for z
in range(self
.dimension()) )
923 return PolynomialRing(self
.base_ring(), var_names
)
925 def inner_product(self
, x
, y
):
927 The inner product associated with this Euclidean Jordan algebra.
929 Defaults to the trace inner product, but can be overridden by
930 subclasses if they are sure that the necessary properties are
935 sage: from mjo.eja.eja_algebra import (random_eja,
937 ....: BilinearFormEJA)
941 Our inner product is "associative," which means the following for
942 a symmetric bilinear form::
944 sage: J = random_eja()
945 sage: x,y,z = J.random_elements(3)
946 sage: (x*y).inner_product(z) == y.inner_product(x*z)
951 Ensure that this is the usual inner product for the algebras
954 sage: J = HadamardEJA.random_instance()
955 sage: x,y = J.random_elements(2)
956 sage: actual = x.inner_product(y)
957 sage: expected = x.to_vector().inner_product(y.to_vector())
958 sage: actual == expected
961 Ensure that this is one-half of the trace inner-product in a
962 BilinearFormEJA that isn't just the reals (when ``n`` isn't
963 one). This is in Faraut and Koranyi, and also my "On the
966 sage: J = BilinearFormEJA.random_instance()
967 sage: n = J.dimension()
968 sage: x = J.random_element()
969 sage: y = J.random_element()
970 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
973 B
= self
._inner
_product
_matrix
974 return (B
*x
.to_vector()).inner_product(y
.to_vector())
977 def is_trivial(self
):
979 Return whether or not this algebra is trivial.
981 A trivial algebra contains only the zero element.
985 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
990 sage: J = ComplexHermitianEJA(3)
996 sage: J = TrivialEJA()
1001 return self
.dimension() == 0
1004 def multiplication_table(self
):
1006 Return a visual representation of this algebra's multiplication
1007 table (on basis elements).
1011 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1015 sage: J = JordanSpinEJA(4)
1016 sage: J.multiplication_table()
1017 +----++----+----+----+----+
1018 | * || b0 | b1 | b2 | b3 |
1019 +====++====+====+====+====+
1020 | b0 || b0 | b1 | b2 | b3 |
1021 +----++----+----+----+----+
1022 | b1 || b1 | b0 | 0 | 0 |
1023 +----++----+----+----+----+
1024 | b2 || b2 | 0 | b0 | 0 |
1025 +----++----+----+----+----+
1026 | b3 || b3 | 0 | 0 | b0 |
1027 +----++----+----+----+----+
1030 n
= self
.dimension()
1031 # Prepend the header row.
1032 M
= [["*"] + list(self
.gens())]
1034 # And to each subsequent row, prepend an entry that belongs to
1035 # the left-side "header column."
1036 M
+= [ [self
.monomial(i
)] + [ self
.monomial(i
)*self
.monomial(j
)
1040 return table(M
, header_row
=True, header_column
=True, frame
=True)
1043 def matrix_basis(self
):
1045 Return an (often more natural) representation of this algebras
1046 basis as an ordered tuple of matrices.
1048 Every finite-dimensional Euclidean Jordan Algebra is a, up to
1049 Jordan isomorphism, a direct sum of five simple
1050 algebras---four of which comprise Hermitian matrices. And the
1051 last type of algebra can of course be thought of as `n`-by-`1`
1052 column matrices (ambiguusly called column vectors) to avoid
1053 special cases. As a result, matrices (and column vectors) are
1054 a natural representation format for Euclidean Jordan algebra
1057 But, when we construct an algebra from a basis of matrices,
1058 those matrix representations are lost in favor of coordinate
1059 vectors *with respect to* that basis. We could eventually
1060 convert back if we tried hard enough, but having the original
1061 representations handy is valuable enough that we simply store
1062 them and return them from this method.
1064 Why implement this for non-matrix algebras? Avoiding special
1065 cases for the :class:`BilinearFormEJA` pays with simplicity in
1066 its own right. But mainly, we would like to be able to assume
1067 that elements of a :class:`CartesianProductEJA` can be displayed
1068 nicely, without having to have special classes for direct sums
1069 one of whose components was a matrix algebra.
1073 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1074 ....: RealSymmetricEJA)
1078 sage: J = RealSymmetricEJA(2)
1080 Finite family {0: b0, 1: b1, 2: b2}
1081 sage: J.matrix_basis()
1083 [1 0] [ 0 0.7071067811865475?] [0 0]
1084 [0 0], [0.7071067811865475? 0], [0 1]
1089 sage: J = JordanSpinEJA(2)
1091 Finite family {0: b0, 1: b1}
1092 sage: J.matrix_basis()
1098 return self
._matrix
_basis
1101 def matrix_space(self
):
1103 Return the matrix space in which this algebra's elements live, if
1104 we think of them as matrices (including column vectors of the
1107 "By default" this will be an `n`-by-`1` column-matrix space,
1108 except when the algebra is trivial. There it's `n`-by-`n`
1109 (where `n` is zero), to ensure that two elements of the matrix
1110 space (empty matrices) can be multiplied. For algebras of
1111 matrices, this returns the space in which their
1112 real embeddings live.
1116 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
1117 ....: JordanSpinEJA,
1118 ....: QuaternionHermitianEJA,
1123 By default, the matrix representation is just a column-matrix
1124 equivalent to the vector representation::
1126 sage: J = JordanSpinEJA(3)
1127 sage: J.matrix_space()
1128 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
1131 The matrix representation in the trivial algebra is
1132 zero-by-zero instead of the usual `n`-by-one::
1134 sage: J = TrivialEJA()
1135 sage: J.matrix_space()
1136 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
1139 The matrix space for complex/quaternion Hermitian matrix EJA
1140 is the space in which their real-embeddings live, not the
1141 original complex/quaternion matrix space::
1143 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1144 sage: J.matrix_space()
1145 Module of 2 by 2 matrices with entries in Algebraic Field over
1146 the scalar ring Rational Field
1147 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1148 sage: J.matrix_space()
1149 Module of 1 by 1 matrices with entries in Quaternion
1150 Algebra (-1, -1) with base ring Rational Field over
1151 the scalar ring Rational Field
1154 return self
._matrix
_space
1160 Return the unit element of this algebra.
1164 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1169 We can compute unit element in the Hadamard EJA::
1171 sage: J = HadamardEJA(5)
1173 b0 + b1 + b2 + b3 + b4
1175 The unit element in the Hadamard EJA is inherited in the
1176 subalgebras generated by its elements::
1178 sage: J = HadamardEJA(5)
1180 b0 + b1 + b2 + b3 + b4
1181 sage: x = sum(J.gens())
1182 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1185 sage: A.one().superalgebra_element()
1186 b0 + b1 + b2 + b3 + b4
1190 The identity element acts like the identity, regardless of
1191 whether or not we orthonormalize::
1193 sage: J = random_eja()
1194 sage: x = J.random_element()
1195 sage: J.one()*x == x and x*J.one() == x
1197 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1198 sage: y = A.random_element()
1199 sage: A.one()*y == y and y*A.one() == y
1204 sage: J = random_eja(field=QQ, orthonormalize=False)
1205 sage: x = J.random_element()
1206 sage: J.one()*x == x and x*J.one() == x
1208 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1209 sage: y = A.random_element()
1210 sage: A.one()*y == y and y*A.one() == y
1213 The matrix of the unit element's operator is the identity,
1214 regardless of the base field and whether or not we
1217 sage: J = random_eja()
1218 sage: actual = J.one().operator().matrix()
1219 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1220 sage: actual == expected
1222 sage: x = J.random_element()
1223 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1224 sage: actual = A.one().operator().matrix()
1225 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1226 sage: actual == expected
1231 sage: J = random_eja(field=QQ, orthonormalize=False)
1232 sage: actual = J.one().operator().matrix()
1233 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1234 sage: actual == expected
1236 sage: x = J.random_element()
1237 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1238 sage: actual = A.one().operator().matrix()
1239 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1240 sage: actual == expected
1243 Ensure that the cached unit element (often precomputed by
1244 hand) agrees with the computed one::
1246 sage: J = random_eja()
1247 sage: cached = J.one()
1248 sage: J.one.clear_cache()
1249 sage: J.one() == cached
1254 sage: J = random_eja(field=QQ, orthonormalize=False)
1255 sage: cached = J.one()
1256 sage: J.one.clear_cache()
1257 sage: J.one() == cached
1261 # We can brute-force compute the matrices of the operators
1262 # that correspond to the basis elements of this algebra.
1263 # If some linear combination of those basis elements is the
1264 # algebra identity, then the same linear combination of
1265 # their matrices has to be the identity matrix.
1267 # Of course, matrices aren't vectors in sage, so we have to
1268 # appeal to the "long vectors" isometry.
1270 V
= VectorSpace(self
.base_ring(), self
.dimension()**2)
1271 oper_vecs
= [ V(g
.operator().matrix().list()) for g
in self
.gens() ]
1273 # Now we use basic linear algebra to find the coefficients,
1274 # of the matrices-as-vectors-linear-combination, which should
1275 # work for the original algebra basis too.
1276 A
= matrix(self
.base_ring(), oper_vecs
)
1278 # We used the isometry on the left-hand side already, but we
1279 # still need to do it for the right-hand side. Recall that we
1280 # wanted something that summed to the identity matrix.
1281 b
= V( matrix
.identity(self
.base_ring(), self
.dimension()).list() )
1283 # Now if there's an identity element in the algebra, this
1284 # should work. We solve on the left to avoid having to
1285 # transpose the matrix "A".
1286 return self
.from_vector(A
.solve_left(b
))
1289 def peirce_decomposition(self
, c
):
1291 The Peirce decomposition of this algebra relative to the
1294 In the future, this can be extended to a complete system of
1295 orthogonal idempotents.
1299 - ``c`` -- an idempotent of this algebra.
1303 A triple (J0, J5, J1) containing two subalgebras and one subspace
1306 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1307 corresponding to the eigenvalue zero.
1309 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1310 corresponding to the eigenvalue one-half.
1312 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1313 corresponding to the eigenvalue one.
1315 These are the only possible eigenspaces for that operator, and this
1316 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1317 orthogonal, and are subalgebras of this algebra with the appropriate
1322 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1326 The canonical example comes from the symmetric matrices, which
1327 decompose into diagonal and off-diagonal parts::
1329 sage: J = RealSymmetricEJA(3)
1330 sage: C = matrix(QQ, [ [1,0,0],
1334 sage: J0,J5,J1 = J.peirce_decomposition(c)
1336 Euclidean Jordan algebra of dimension 1...
1338 Vector space of degree 6 and dimension 2...
1340 Euclidean Jordan algebra of dimension 3...
1341 sage: J0.one().to_matrix()
1345 sage: orig_df = AA.options.display_format
1346 sage: AA.options.display_format = 'radical'
1347 sage: J.from_vector(J5.basis()[0]).to_matrix()
1351 sage: J.from_vector(J5.basis()[1]).to_matrix()
1355 sage: AA.options.display_format = orig_df
1356 sage: J1.one().to_matrix()
1363 Every algebra decomposes trivially with respect to its identity
1366 sage: J = random_eja()
1367 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1368 sage: J0.dimension() == 0 and J5.dimension() == 0
1370 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1373 The decomposition is into eigenspaces, and its components are
1374 therefore necessarily orthogonal. Moreover, the identity
1375 elements in the two subalgebras are the projections onto their
1376 respective subspaces of the superalgebra's identity element::
1378 sage: J = random_eja()
1379 sage: x = J.random_element()
1380 sage: if not J.is_trivial():
1381 ....: while x.is_nilpotent():
1382 ....: x = J.random_element()
1383 sage: c = x.subalgebra_idempotent()
1384 sage: J0,J5,J1 = J.peirce_decomposition(c)
1386 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1387 ....: w = w.superalgebra_element()
1388 ....: y = J.from_vector(y)
1389 ....: z = z.superalgebra_element()
1390 ....: ipsum += w.inner_product(y).abs()
1391 ....: ipsum += w.inner_product(z).abs()
1392 ....: ipsum += y.inner_product(z).abs()
1395 sage: J1(c) == J1.one()
1397 sage: J0(J.one() - c) == J0.one()
1401 if not c
.is_idempotent():
1402 raise ValueError("element is not idempotent: %s" % c
)
1404 # Default these to what they should be if they turn out to be
1405 # trivial, because eigenspaces_left() won't return eigenvalues
1406 # corresponding to trivial spaces (e.g. it returns only the
1407 # eigenspace corresponding to lambda=1 if you take the
1408 # decomposition relative to the identity element).
1409 trivial
= self
.subalgebra((), check_axioms
=False)
1410 J0
= trivial
# eigenvalue zero
1411 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1412 J1
= trivial
# eigenvalue one
1414 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1415 if eigval
== ~
(self
.base_ring()(2)):
1418 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1419 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1425 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1430 def random_element(self
, thorough
=False):
1432 Return a random element of this algebra.
1434 Our algebra superclass method only returns a linear
1435 combination of at most two basis elements. We instead
1436 want the vector space "random element" method that
1437 returns a more diverse selection.
1441 - ``thorough`` -- (boolean; default False) whether or not we
1442 should generate irrational coefficients for the random
1443 element when our base ring is irrational; this slows the
1444 algebra operations to a crawl, but any truly random method
1448 # For a general base ring... maybe we can trust this to do the
1449 # right thing? Unlikely, but.
1450 V
= self
.vector_space()
1451 if self
.base_ring() is AA
and not thorough
:
1452 # Now that AA generates actually random random elements
1453 # (post Trac 30875), we only need to de-thorough the
1454 # randomness when asked to.
1455 V
= V
.change_ring(QQ
)
1457 v
= V
.random_element()
1458 return self
.from_vector(V
.coordinate_vector(v
))
1460 def random_elements(self
, count
, thorough
=False):
1462 Return ``count`` random elements as a tuple.
1466 - ``thorough`` -- (boolean; default False) whether or not we
1467 should generate irrational coefficients for the random
1468 elements when our base ring is irrational; this slows the
1469 algebra operations to a crawl, but any truly random method
1474 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1478 sage: J = JordanSpinEJA(3)
1479 sage: x,y,z = J.random_elements(3)
1480 sage: all( [ x in J, y in J, z in J ])
1482 sage: len( J.random_elements(10) ) == 10
1486 return tuple( self
.random_element(thorough
)
1487 for idx
in range(count
) )
1490 def operator_polynomial_matrix(self
):
1492 Return the matrix of polynomials (over this algebra's
1493 :meth:`coordinate_polynomial_ring`) that, when evaluated at
1494 the basis coordinates of an element `x`, produces the basis
1495 representation of `L_{x}`.
1499 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1500 ....: JordanSpinEJA)
1504 sage: J = HadamardEJA(4)
1505 sage: L_x = J.operator_polynomial_matrix()
1512 sage: d = zip(J.coordinate_polynomial_ring().gens(), x.to_vector())
1513 sage: L_x.subs(dict(d))
1521 sage: J = JordanSpinEJA(4)
1522 sage: L_x = J.operator_polynomial_matrix()
1529 sage: d = zip(J.coordinate_polynomial_ring().gens(), x.to_vector())
1530 sage: L_x.subs(dict(d))
1537 R
= self
.coordinate_polynomial_ring()
1540 # From a result in my book, these are the entries of the
1541 # basis representation of L_x.
1542 return sum( v
*self
.monomial(k
).operator().matrix()[i
,j
]
1543 for (k
,v
) in enumerate(R
.gens()) )
1545 n
= self
.dimension()
1546 return matrix(R
, n
, n
, L_x_i_j
)
1549 def _charpoly_coefficients(self
):
1551 The `r` polynomial coefficients of the "characteristic polynomial
1556 sage: from mjo.eja.eja_algebra import random_eja
1560 The theory shows that these are all homogeneous polynomials of
1563 sage: J = random_eja()
1564 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1568 n
= self
.dimension()
1569 R
= self
.coordinate_polynomial_ring()
1570 F
= R
.fraction_field()
1572 L_x
= self
.operator_polynomial_matrix()
1575 if self
.rank
.is_in_cache():
1577 # There's no need to pad the system with redundant
1578 # columns if we *know* they'll be redundant.
1581 # Compute an extra power in case the rank is equal to
1582 # the dimension (otherwise, we would stop at x^(r-1)).
1583 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1584 for k
in range(n
+1) ]
1585 A
= matrix
.column(F
, x_powers
[:n
])
1586 AE
= A
.extended_echelon_form()
1593 # The theory says that only the first "r" coefficients are
1594 # nonzero, and they actually live in the original polynomial
1595 # ring and not the fraction field. We negate them because in
1596 # the actual characteristic polynomial, they get moved to the
1597 # other side where x^r lives. We don't bother to trim A_rref
1598 # down to a square matrix and solve the resulting system,
1599 # because the upper-left r-by-r portion of A_rref is
1600 # guaranteed to be the identity matrix, so e.g.
1602 # A_rref.solve_right(Y)
1604 # would just be returning Y.
1605 return (-E
*b
)[:r
].change_ring(R
)
1610 Return the rank of this EJA.
1612 This is a cached method because we know the rank a priori for
1613 all of the algebras we can construct. Thus we can avoid the
1614 expensive ``_charpoly_coefficients()`` call unless we truly
1615 need to compute the whole characteristic polynomial.
1619 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1620 ....: JordanSpinEJA,
1621 ....: RealSymmetricEJA,
1622 ....: ComplexHermitianEJA,
1623 ....: QuaternionHermitianEJA,
1628 The rank of the Jordan spin algebra is always two::
1630 sage: JordanSpinEJA(2).rank()
1632 sage: JordanSpinEJA(3).rank()
1634 sage: JordanSpinEJA(4).rank()
1637 The rank of the `n`-by-`n` Hermitian real, complex, or
1638 quaternion matrices is `n`::
1640 sage: RealSymmetricEJA(4).rank()
1642 sage: ComplexHermitianEJA(3).rank()
1644 sage: QuaternionHermitianEJA(2).rank()
1649 Ensure that every EJA that we know how to construct has a
1650 positive integer rank, unless the algebra is trivial in
1651 which case its rank will be zero::
1653 sage: J = random_eja()
1657 sage: r > 0 or (r == 0 and J.is_trivial())
1660 Ensure that computing the rank actually works, since the ranks
1661 of all simple algebras are known and will be cached by default::
1663 sage: J = random_eja() # long time
1664 sage: cached = J.rank() # long time
1665 sage: J.rank.clear_cache() # long time
1666 sage: J.rank() == cached # long time
1670 return len(self
._charpoly
_coefficients
())
1673 def subalgebra(self
, basis
, **kwargs
):
1675 Create a subalgebra of this algebra from the given basis.
1677 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1678 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1681 def vector_space(self
):
1683 Return the vector space that underlies this algebra.
1687 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1691 sage: J = RealSymmetricEJA(2)
1692 sage: J.vector_space()
1693 Vector space of dimension 3 over...
1696 return self
.zero().to_vector().parent().ambient_vector_space()
1700 class RationalBasisEJA(FiniteDimensionalEJA
):
1702 Algebras whose supplied basis elements have all rational entries.
1706 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1710 The supplied basis is orthonormalized by default::
1712 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1713 sage: J = BilinearFormEJA(B)
1714 sage: J.matrix_basis()
1731 # Abuse the check_field parameter to check that the entries of
1732 # out basis (in ambient coordinates) are in the field QQ.
1733 # Use _all2list to get the vector coordinates of octonion
1734 # entries and not the octonions themselves (which are not
1736 if not all( all(b_i
in QQ
for b_i
in _all2list(b
))
1738 raise TypeError("basis not rational")
1740 super().__init
__(basis
,
1744 check_field
=check_field
,
1747 self
._rational
_algebra
= None
1749 # There's no point in constructing the extra algebra if this
1750 # one is already rational.
1752 # Note: the same Jordan and inner-products work here,
1753 # because they are necessarily defined with respect to
1754 # ambient coordinates and not any particular basis.
1755 self
._rational
_algebra
= FiniteDimensionalEJA(
1760 matrix_space
=self
.matrix_space(),
1761 associative
=self
.is_associative(),
1762 orthonormalize
=False,
1766 def rational_algebra(self
):
1767 # Using None as a flag here (rather than just assigning "self"
1768 # to self._rational_algebra by default) feels a little bit
1769 # more sane to me in a garbage-collected environment.
1770 if self
._rational
_algebra
is None:
1773 return self
._rational
_algebra
1776 def _charpoly_coefficients(self
):
1780 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1781 ....: JordanSpinEJA)
1785 The base ring of the resulting polynomial coefficients is what
1786 it should be, and not the rationals (unless the algebra was
1787 already over the rationals)::
1789 sage: J = JordanSpinEJA(3)
1790 sage: J._charpoly_coefficients()
1791 (X0^2 - X1^2 - X2^2, -2*X0)
1792 sage: a0 = J._charpoly_coefficients()[0]
1794 Algebraic Real Field
1795 sage: a0.base_ring()
1796 Algebraic Real Field
1799 if self
.rational_algebra() is self
:
1800 # Bypass the hijinks if they won't benefit us.
1801 return super()._charpoly
_coefficients
()
1803 # Do the computation over the rationals.
1804 a
= ( a_i
.change_ring(self
.base_ring())
1805 for a_i
in self
.rational_algebra()._charpoly
_coefficients
() )
1807 # Convert our coordinate variables into deorthonormalized ones
1808 # and substitute them into the deorthonormalized charpoly
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
2192 The dimension of this algebra is `n^2`::
2194 sage: d = ComplexHermitianEJA._max_random_instance_dimension()
2195 sage: n = ComplexHermitianEJA._max_random_instance_size(d)
2196 sage: J = ComplexHermitianEJA(n)
2197 sage: J.dimension() == n^2
2200 The Jordan multiplication is what we think it is::
2202 sage: J = ComplexHermitianEJA.random_instance()
2203 sage: x,y = J.random_elements(2)
2204 sage: actual = (x*y).to_matrix()
2205 sage: X = x.to_matrix()
2206 sage: Y = y.to_matrix()
2207 sage: expected = (X*Y + Y*X)/2
2208 sage: actual == expected
2210 sage: J(expected) == x*y
2213 We can change the generator prefix::
2215 sage: ComplexHermitianEJA(2, prefix='z').gens()
2218 We can construct the (trivial) algebra of rank zero::
2220 sage: ComplexHermitianEJA(0)
2221 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2224 def __init__(self
, n
, field
=AA
, **kwargs
):
2225 from mjo
.hurwitz
import ComplexMatrixAlgebra
2226 A
= ComplexMatrixAlgebra(n
, scalars
=field
)
2227 super().__init
__(A
, **kwargs
)
2229 from mjo
.eja
.eja_cache
import complex_hermitian_eja_coeffs
2230 a
= complex_hermitian_eja_coeffs(self
)
2232 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2235 def _max_random_instance_size(max_dimension
):
2236 # Obtained by solving d = n^2.
2237 # The ZZ-int-ZZ thing is just "floor."
2238 return ZZ(int(ZZ(max_dimension
).sqrt()))
2241 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2243 Return a random instance of this type of algebra.
2245 class_max_d
= cls
._max
_random
_instance
_dimension
()
2246 if (max_dimension
is None or max_dimension
> class_max_d
):
2247 max_dimension
= class_max_d
2248 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2249 n
= ZZ
.random_element(max_size
+ 1)
2250 return cls(n
, **kwargs
)
2253 class QuaternionHermitianEJA(HermitianMatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2255 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2256 matrices, the usual symmetric Jordan product, and the
2257 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2262 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2266 In theory, our "field" can be any subfield of the reals::
2268 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2269 Euclidean Jordan algebra of dimension 6 over Real Double Field
2270 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2271 Euclidean Jordan algebra of dimension 6 over Real Field with
2272 53 bits of precision
2276 The dimension of this algebra is `2*n^2 - n`::
2278 sage: d = QuaternionHermitianEJA._max_random_instance_dimension()
2279 sage: n = QuaternionHermitianEJA._max_random_instance_size(d)
2280 sage: J = QuaternionHermitianEJA(n)
2281 sage: J.dimension() == 2*(n^2) - n
2284 The Jordan multiplication is what we think it is::
2286 sage: J = QuaternionHermitianEJA.random_instance()
2287 sage: x,y = J.random_elements(2)
2288 sage: actual = (x*y).to_matrix()
2289 sage: X = x.to_matrix()
2290 sage: Y = y.to_matrix()
2291 sage: expected = (X*Y + Y*X)/2
2292 sage: actual == expected
2294 sage: J(expected) == x*y
2297 We can change the generator prefix::
2299 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2300 (a0, a1, a2, a3, a4, a5)
2302 We can construct the (trivial) algebra of rank zero::
2304 sage: QuaternionHermitianEJA(0)
2305 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2308 def __init__(self
, n
, field
=AA
, **kwargs
):
2309 from mjo
.hurwitz
import QuaternionMatrixAlgebra
2310 A
= QuaternionMatrixAlgebra(n
, scalars
=field
)
2311 super().__init
__(A
, **kwargs
)
2313 from mjo
.eja
.eja_cache
import quaternion_hermitian_eja_coeffs
2314 a
= quaternion_hermitian_eja_coeffs(self
)
2316 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2321 def _max_random_instance_size(max_dimension
):
2323 The maximum rank of a random QuaternionHermitianEJA.
2325 # Obtained by solving d = 2n^2 - n.
2326 # The ZZ-int-ZZ thing is just "floor."
2327 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/4 + 1/4))
2330 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2332 Return a random instance of this type of algebra.
2334 class_max_d
= cls
._max
_random
_instance
_dimension
()
2335 if (max_dimension
is None or max_dimension
> class_max_d
):
2336 max_dimension
= class_max_d
2337 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2338 n
= ZZ
.random_element(max_size
+ 1)
2339 return cls(n
, **kwargs
)
2341 class OctonionHermitianEJA(HermitianMatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2345 sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA,
2346 ....: OctonionHermitianEJA)
2347 sage: from mjo.hurwitz import Octonions, OctonionMatrixAlgebra
2351 The 3-by-3 algebra satisfies the axioms of an EJA::
2353 sage: OctonionHermitianEJA(3, # long time
2354 ....: field=QQ, # long time
2355 ....: orthonormalize=False, # long time
2356 ....: check_axioms=True) # long time
2357 Euclidean Jordan algebra of dimension 27 over Rational Field
2359 After a change-of-basis, the 2-by-2 algebra has the same
2360 multiplication table as the ten-dimensional Jordan spin algebra::
2362 sage: A = OctonionMatrixAlgebra(2,Octonions(QQ),QQ)
2363 sage: b = OctonionHermitianEJA._denormalized_basis(A)
2364 sage: basis = (b[0] + b[9],) + b[1:9] + (b[0] - b[9],)
2365 sage: jp = OctonionHermitianEJA.jordan_product
2366 sage: ip = OctonionHermitianEJA.trace_inner_product
2367 sage: J = FiniteDimensionalEJA(basis,
2371 ....: orthonormalize=False)
2372 sage: J.multiplication_table()
2373 +----++----+----+----+----+----+----+----+----+----+----+
2374 | * || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 |
2375 +====++====+====+====+====+====+====+====+====+====+====+
2376 | b0 || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 |
2377 +----++----+----+----+----+----+----+----+----+----+----+
2378 | b1 || b1 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
2379 +----++----+----+----+----+----+----+----+----+----+----+
2380 | b2 || b2 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
2381 +----++----+----+----+----+----+----+----+----+----+----+
2382 | b3 || b3 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 |
2383 +----++----+----+----+----+----+----+----+----+----+----+
2384 | b4 || b4 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 |
2385 +----++----+----+----+----+----+----+----+----+----+----+
2386 | b5 || b5 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 |
2387 +----++----+----+----+----+----+----+----+----+----+----+
2388 | b6 || b6 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 |
2389 +----++----+----+----+----+----+----+----+----+----+----+
2390 | b7 || b7 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 |
2391 +----++----+----+----+----+----+----+----+----+----+----+
2392 | b8 || b8 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 |
2393 +----++----+----+----+----+----+----+----+----+----+----+
2394 | b9 || b9 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 |
2395 +----++----+----+----+----+----+----+----+----+----+----+
2399 We can actually construct the 27-dimensional Albert algebra,
2400 and we get the right unit element if we recompute it::
2402 sage: J = OctonionHermitianEJA(3, # long time
2403 ....: field=QQ, # long time
2404 ....: orthonormalize=False) # long time
2405 sage: J.one.clear_cache() # long time
2406 sage: J.one() # long time
2408 sage: J.one().to_matrix() # long time
2417 The 2-by-2 algebra is isomorphic to the ten-dimensional Jordan
2418 spin algebra, but just to be sure, we recompute its rank::
2420 sage: J = OctonionHermitianEJA(2, # long time
2421 ....: field=QQ, # long time
2422 ....: orthonormalize=False) # long time
2423 sage: J.rank.clear_cache() # long time
2424 sage: J.rank() # long time
2429 def _max_random_instance_size(max_dimension
):
2431 The maximum rank of a random OctonionHermitianEJA.
2433 # There's certainly a formula for this, but with only four
2434 # cases to worry about, I'm not that motivated to derive it.
2435 if max_dimension
>= 27:
2437 elif max_dimension
>= 10:
2439 elif max_dimension
>= 1:
2445 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2447 Return a random instance of this type of algebra.
2449 class_max_d
= cls
._max
_random
_instance
_dimension
()
2450 if (max_dimension
is None or max_dimension
> class_max_d
):
2451 max_dimension
= class_max_d
2452 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2453 n
= ZZ
.random_element(max_size
+ 1)
2454 return cls(n
, **kwargs
)
2456 def __init__(self
, n
, field
=AA
, **kwargs
):
2458 # Otherwise we don't get an EJA.
2459 raise ValueError("n cannot exceed 3")
2461 # We know this is a valid EJA, but will double-check
2462 # if the user passes check_axioms=True.
2463 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2465 from mjo
.hurwitz
import OctonionMatrixAlgebra
2466 A
= OctonionMatrixAlgebra(n
, scalars
=field
)
2467 super().__init
__(A
, **kwargs
)
2469 from mjo
.eja
.eja_cache
import octonion_hermitian_eja_coeffs
2470 a
= octonion_hermitian_eja_coeffs(self
)
2472 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2475 class AlbertEJA(OctonionHermitianEJA
):
2477 The Albert algebra is the algebra of three-by-three Hermitian
2478 matrices whose entries are octonions.
2482 sage: from mjo.eja.eja_algebra import AlbertEJA
2486 sage: AlbertEJA(field=QQ, orthonormalize=False)
2487 Euclidean Jordan algebra of dimension 27 over Rational Field
2488 sage: AlbertEJA() # long time
2489 Euclidean Jordan algebra of dimension 27 over Algebraic Real Field
2492 def __init__(self
, *args
, **kwargs
):
2493 super().__init
__(3, *args
, **kwargs
)
2496 class HadamardEJA(RationalBasisEJA
, ConcreteEJA
):
2498 Return the Euclidean Jordan algebra on `R^n` with the Hadamard
2499 (pointwise real-number multiplication) Jordan product and the
2500 usual inner-product.
2502 This is nothing more than the Cartesian product of ``n`` copies of
2503 the one-dimensional Jordan spin algebra, and is the most common
2504 example of a non-simple Euclidean Jordan algebra.
2508 sage: from mjo.eja.eja_algebra import HadamardEJA
2512 This multiplication table can be verified by hand::
2514 sage: J = HadamardEJA(3)
2515 sage: b0,b1,b2 = J.gens()
2531 We can change the generator prefix::
2533 sage: HadamardEJA(3, prefix='r').gens()
2536 def __init__(self
, n
, field
=AA
, **kwargs
):
2537 MS
= MatrixSpace(field
, n
, 1)
2540 jordan_product
= lambda x
,y
: x
2541 inner_product
= lambda x
,y
: x
2543 def jordan_product(x
,y
):
2544 return MS( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2546 def inner_product(x
,y
):
2549 # New defaults for keyword arguments. Don't orthonormalize
2550 # because our basis is already orthonormal with respect to our
2551 # inner-product. Don't check the axioms, because we know this
2552 # is a valid EJA... but do double-check if the user passes
2553 # check_axioms=True. Note: we DON'T override the "check_field"
2554 # default here, because the user can pass in a field!
2555 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2556 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2558 column_basis
= tuple( MS(b
) for b
in FreeModule(field
, n
).basis() )
2559 super().__init
__(column_basis
,
2566 self
.rank
.set_cache(n
)
2568 self
.one
.set_cache( self
.sum(self
.gens()) )
2571 def _max_random_instance_dimension():
2573 There's no reason to go higher than five here. That's
2574 enough to get the point across.
2579 def _max_random_instance_size(max_dimension
):
2581 The maximum size (=dimension) of a random HadamardEJA.
2583 return max_dimension
2586 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2588 Return a random instance of this type of algebra.
2590 class_max_d
= cls
._max
_random
_instance
_dimension
()
2591 if (max_dimension
is None or max_dimension
> class_max_d
):
2592 max_dimension
= class_max_d
2593 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2594 n
= ZZ
.random_element(max_size
+ 1)
2595 return cls(n
, **kwargs
)
2598 class BilinearFormEJA(RationalBasisEJA
, ConcreteEJA
):
2600 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2601 with the half-trace inner product and jordan product ``x*y =
2602 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2603 a symmetric positive-definite "bilinear form" matrix. Its
2604 dimension is the size of `B`, and it has rank two in dimensions
2605 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2606 the identity matrix of order ``n``.
2608 We insist that the one-by-one upper-left identity block of `B` be
2609 passed in as well so that we can be passed a matrix of size zero
2610 to construct a trivial algebra.
2614 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2615 ....: JordanSpinEJA)
2619 When no bilinear form is specified, the identity matrix is used,
2620 and the resulting algebra is the Jordan spin algebra::
2622 sage: B = matrix.identity(AA,3)
2623 sage: J0 = BilinearFormEJA(B)
2624 sage: J1 = JordanSpinEJA(3)
2625 sage: J0.multiplication_table() == J0.multiplication_table()
2628 An error is raised if the matrix `B` does not correspond to a
2629 positive-definite bilinear form::
2631 sage: B = matrix.random(QQ,2,3)
2632 sage: J = BilinearFormEJA(B)
2633 Traceback (most recent call last):
2635 ValueError: bilinear form is not positive-definite
2636 sage: B = matrix.zero(QQ,3)
2637 sage: J = BilinearFormEJA(B)
2638 Traceback (most recent call last):
2640 ValueError: bilinear form is not positive-definite
2644 We can create a zero-dimensional algebra::
2646 sage: B = matrix.identity(AA,0)
2647 sage: J = BilinearFormEJA(B)
2651 We can check the multiplication condition given in the Jordan, von
2652 Neumann, and Wigner paper (and also discussed on my "On the
2653 symmetry..." paper). Note that this relies heavily on the standard
2654 choice of basis, as does anything utilizing the bilinear form
2655 matrix. We opt not to orthonormalize the basis, because if we
2656 did, we would have to normalize the `s_{i}` in a similar manner::
2658 sage: n = ZZ.random_element(5)
2659 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2660 sage: B11 = matrix.identity(QQ,1)
2661 sage: B22 = M.transpose()*M
2662 sage: B = block_matrix(2,2,[ [B11,0 ],
2664 sage: J = BilinearFormEJA(B, orthonormalize=False)
2665 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2666 sage: V = J.vector_space()
2667 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2668 ....: for ei in eis ]
2669 sage: actual = [ sis[i]*sis[j]
2670 ....: for i in range(n-1)
2671 ....: for j in range(n-1) ]
2672 sage: expected = [ J.one() if i == j else J.zero()
2673 ....: for i in range(n-1)
2674 ....: for j in range(n-1) ]
2675 sage: actual == expected
2679 def __init__(self
, B
, field
=AA
, **kwargs
):
2680 # The matrix "B" is supplied by the user in most cases,
2681 # so it makes sense to check whether or not its positive-
2682 # definite unless we are specifically asked not to...
2683 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2684 if not B
.is_positive_definite():
2685 raise ValueError("bilinear form is not positive-definite")
2687 # However, all of the other data for this EJA is computed
2688 # by us in manner that guarantees the axioms are
2689 # satisfied. So, again, unless we are specifically asked to
2690 # verify things, we'll skip the rest of the checks.
2691 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2694 MS
= MatrixSpace(field
, n
, 1)
2696 def inner_product(x
,y
):
2697 return (y
.T
*B
*x
)[0,0]
2699 def jordan_product(x
,y
):
2704 z0
= inner_product(y
,x
)
2705 zbar
= y0
*xbar
+ x0
*ybar
2706 return MS([z0
] + zbar
.list())
2708 column_basis
= tuple( MS(b
) for b
in FreeModule(field
, n
).basis() )
2710 # TODO: I haven't actually checked this, but it seems legit.
2715 super().__init
__(column_basis
,
2720 associative
=associative
,
2723 # The rank of this algebra is two, unless we're in a
2724 # one-dimensional ambient space (because the rank is bounded
2725 # by the ambient dimension).
2726 self
.rank
.set_cache(min(n
,2))
2728 self
.one
.set_cache( self
.zero() )
2730 self
.one
.set_cache( self
.monomial(0) )
2733 def _max_random_instance_dimension():
2735 There's no reason to go higher than five here. That's
2736 enough to get the point across.
2741 def _max_random_instance_size(max_dimension
):
2743 The maximum size (=dimension) of a random BilinearFormEJA.
2745 return max_dimension
2748 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2750 Return a random instance of this algebra.
2752 class_max_d
= cls
._max
_random
_instance
_dimension
()
2753 if (max_dimension
is None or max_dimension
> class_max_d
):
2754 max_dimension
= class_max_d
2755 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2756 n
= ZZ
.random_element(max_size
+ 1)
2759 B
= matrix
.identity(ZZ
, n
)
2760 return cls(B
, **kwargs
)
2762 B11
= matrix
.identity(ZZ
, 1)
2763 M
= matrix
.random(ZZ
, n
-1)
2764 I
= matrix
.identity(ZZ
, n
-1)
2766 while alpha
.is_zero():
2767 alpha
= ZZ
.random_element().abs()
2769 B22
= M
.transpose()*M
+ alpha
*I
2771 from sage
.matrix
.special
import block_matrix
2772 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2775 return cls(B
, **kwargs
)
2778 class JordanSpinEJA(BilinearFormEJA
):
2780 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2781 with the usual inner product and jordan product ``x*y =
2782 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2787 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2791 This multiplication table can be verified by hand::
2793 sage: J = JordanSpinEJA(4)
2794 sage: b0,b1,b2,b3 = J.gens()
2810 We can change the generator prefix::
2812 sage: JordanSpinEJA(2, prefix='B').gens()
2817 Ensure that we have the usual inner product on `R^n`::
2819 sage: J = JordanSpinEJA.random_instance()
2820 sage: x,y = J.random_elements(2)
2821 sage: actual = x.inner_product(y)
2822 sage: expected = x.to_vector().inner_product(y.to_vector())
2823 sage: actual == expected
2827 def __init__(self
, n
, *args
, **kwargs
):
2828 # This is a special case of the BilinearFormEJA with the
2829 # identity matrix as its bilinear form.
2830 B
= matrix
.identity(ZZ
, n
)
2832 # Don't orthonormalize because our basis is already
2833 # orthonormal with respect to our inner-product.
2834 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2836 # But also don't pass check_field=False here, because the user
2837 # can pass in a field!
2838 super().__init
__(B
, *args
, **kwargs
)
2841 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2843 Return a random instance of this type of algebra.
2845 Needed here to override the implementation for ``BilinearFormEJA``.
2847 class_max_d
= cls
._max
_random
_instance
_dimension
()
2848 if (max_dimension
is None or max_dimension
> class_max_d
):
2849 max_dimension
= class_max_d
2850 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2851 n
= ZZ
.random_element(max_size
+ 1)
2852 return cls(n
, **kwargs
)
2855 class TrivialEJA(RationalBasisEJA
, ConcreteEJA
):
2857 The trivial Euclidean Jordan algebra consisting of only a zero element.
2861 sage: from mjo.eja.eja_algebra import TrivialEJA
2865 sage: J = TrivialEJA()
2872 sage: 7*J.one()*12*J.one()
2874 sage: J.one().inner_product(J.one())
2876 sage: J.one().norm()
2878 sage: J.one().subalgebra_generated_by()
2879 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2884 def __init__(self
, field
=AA
, **kwargs
):
2885 jordan_product
= lambda x
,y
: x
2886 inner_product
= lambda x
,y
: field
.zero()
2888 MS
= MatrixSpace(field
,0)
2890 # New defaults for keyword arguments
2891 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2892 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2894 super().__init
__(basis
,
2902 # The rank is zero using my definition, namely the dimension of the
2903 # largest subalgebra generated by any element.
2904 self
.rank
.set_cache(0)
2905 self
.one
.set_cache( self
.zero() )
2908 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2909 # We don't take a "size" argument so the superclass method is
2910 # inappropriate for us. The ``max_dimension`` argument is
2911 # included so that if this method is called generically with a
2912 # ``max_dimension=<whatever>`` argument, we don't try to pass
2913 # it on to the algebra constructor.
2914 return cls(**kwargs
)
2917 class CartesianProductEJA(FiniteDimensionalEJA
):
2919 The external (orthogonal) direct sum of two or more Euclidean
2920 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
2921 orthogonal direct sum of simple Euclidean Jordan algebras which is
2922 then isometric to a Cartesian product, so no generality is lost by
2923 providing only this construction.
2927 sage: from mjo.eja.eja_algebra import (random_eja,
2928 ....: CartesianProductEJA,
2929 ....: ComplexHermitianEJA,
2931 ....: JordanSpinEJA,
2932 ....: RealSymmetricEJA)
2936 The Jordan product is inherited from our factors and implemented by
2937 our CombinatorialFreeModule Cartesian product superclass::
2939 sage: J1 = HadamardEJA(2)
2940 sage: J2 = RealSymmetricEJA(2)
2941 sage: J = cartesian_product([J1,J2])
2942 sage: x,y = J.random_elements(2)
2946 The ability to retrieve the original factors is implemented by our
2947 CombinatorialFreeModule Cartesian product superclass::
2949 sage: J1 = HadamardEJA(2, field=QQ)
2950 sage: J2 = JordanSpinEJA(3, field=QQ)
2951 sage: J = cartesian_product([J1,J2])
2952 sage: J.cartesian_factors()
2953 (Euclidean Jordan algebra of dimension 2 over Rational Field,
2954 Euclidean Jordan algebra of dimension 3 over Rational Field)
2956 You can provide more than two factors::
2958 sage: J1 = HadamardEJA(2)
2959 sage: J2 = JordanSpinEJA(3)
2960 sage: J3 = RealSymmetricEJA(3)
2961 sage: cartesian_product([J1,J2,J3])
2962 Euclidean Jordan algebra of dimension 2 over Algebraic Real
2963 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
2964 Real Field (+) Euclidean Jordan algebra of dimension 6 over
2965 Algebraic Real Field
2967 Rank is additive on a Cartesian product::
2969 sage: J1 = HadamardEJA(1)
2970 sage: J2 = RealSymmetricEJA(2)
2971 sage: J = cartesian_product([J1,J2])
2972 sage: J1.rank.clear_cache()
2973 sage: J2.rank.clear_cache()
2974 sage: J.rank.clear_cache()
2977 sage: J.rank() == J1.rank() + J2.rank()
2980 The same rank computation works over the rationals, with whatever
2983 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
2984 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
2985 sage: J = cartesian_product([J1,J2])
2986 sage: J1.rank.clear_cache()
2987 sage: J2.rank.clear_cache()
2988 sage: J.rank.clear_cache()
2991 sage: J.rank() == J1.rank() + J2.rank()
2994 The product algebra will be associative if and only if all of its
2995 components are associative::
2997 sage: J1 = HadamardEJA(2)
2998 sage: J1.is_associative()
3000 sage: J2 = HadamardEJA(3)
3001 sage: J2.is_associative()
3003 sage: J3 = RealSymmetricEJA(3)
3004 sage: J3.is_associative()
3006 sage: CP1 = cartesian_product([J1,J2])
3007 sage: CP1.is_associative()
3009 sage: CP2 = cartesian_product([J1,J3])
3010 sage: CP2.is_associative()
3013 Cartesian products of Cartesian products work::
3015 sage: J1 = JordanSpinEJA(1)
3016 sage: J2 = JordanSpinEJA(1)
3017 sage: J3 = JordanSpinEJA(1)
3018 sage: J = cartesian_product([J1,cartesian_product([J2,J3])])
3019 sage: J.multiplication_table()
3020 +----++----+----+----+
3021 | * || b0 | b1 | b2 |
3022 +====++====+====+====+
3023 | b0 || b0 | 0 | 0 |
3024 +----++----+----+----+
3025 | b1 || 0 | b1 | 0 |
3026 +----++----+----+----+
3027 | b2 || 0 | 0 | b2 |
3028 +----++----+----+----+
3029 sage: HadamardEJA(3).multiplication_table()
3030 +----++----+----+----+
3031 | * || b0 | b1 | b2 |
3032 +====++====+====+====+
3033 | b0 || b0 | 0 | 0 |
3034 +----++----+----+----+
3035 | b1 || 0 | b1 | 0 |
3036 +----++----+----+----+
3037 | b2 || 0 | 0 | b2 |
3038 +----++----+----+----+
3040 The "matrix space" of a Cartesian product always consists of
3041 ordered pairs (or triples, or...) whose components are the
3042 matrix spaces of its factors::
3044 sage: J1 = HadamardEJA(2)
3045 sage: J2 = ComplexHermitianEJA(2)
3046 sage: J = cartesian_product([J1,J2])
3047 sage: J.matrix_space()
3048 The Cartesian product of (Full MatrixSpace of 2 by 1 dense
3049 matrices over Algebraic Real Field, Module of 2 by 2 matrices
3050 with entries in Algebraic Field over the scalar ring Algebraic
3052 sage: J.one().to_matrix()[0]
3055 sage: J.one().to_matrix()[1]
3064 All factors must share the same base field::
3066 sage: J1 = HadamardEJA(2, field=QQ)
3067 sage: J2 = RealSymmetricEJA(2)
3068 sage: CartesianProductEJA((J1,J2))
3069 Traceback (most recent call last):
3071 ValueError: all factors must share the same base field
3073 The cached unit element is the same one that would be computed::
3075 sage: J1 = random_eja() # long time
3076 sage: J2 = random_eja() # long time
3077 sage: J = cartesian_product([J1,J2]) # long time
3078 sage: actual = J.one() # long time
3079 sage: J.one.clear_cache() # long time
3080 sage: expected = J.one() # long time
3081 sage: actual == expected # long time
3084 Element
= CartesianProductEJAElement
3085 def __init__(self
, factors
, **kwargs
):
3090 self
._sets
= factors
3092 field
= factors
[0].base_ring()
3093 if not all( J
.base_ring() == field
for J
in factors
):
3094 raise ValueError("all factors must share the same base field")
3096 # Figure out the category to use.
3097 associative
= all( f
.is_associative() for f
in factors
)
3098 category
= EuclideanJordanAlgebras(field
)
3099 if associative
: category
= category
.Associative()
3100 category
= category
.join([category
, category
.CartesianProducts()])
3102 # Compute my matrix space. We don't simply use the
3103 # ``cartesian_product()`` functor here because it acts
3104 # differently on SageMath MatrixSpaces and our custom
3105 # MatrixAlgebras, which are CombinatorialFreeModules. We
3106 # always want the result to be represented (and indexed) as an
3107 # ordered tuple. This category isn't perfect, but is good
3108 # enough for what we need to do.
3109 MS_cat
= MagmaticAlgebras(field
).FiniteDimensional().WithBasis()
3110 MS_cat
= MS_cat
.Unital().CartesianProducts()
3111 MS_factors
= tuple( J
.matrix_space() for J
in factors
)
3112 from sage
.sets
.cartesian_product
import CartesianProduct
3113 self
._matrix
_space
= CartesianProduct(MS_factors
, MS_cat
)
3115 self
._matrix
_basis
= []
3116 zero
= self
._matrix
_space
.zero()
3118 for b
in factors
[i
].matrix_basis():
3121 self
._matrix
_basis
.append(z
)
3123 self
._matrix
_basis
= tuple( self
._matrix
_space
(b
)
3124 for b
in self
._matrix
_basis
)
3125 n
= len(self
._matrix
_basis
)
3127 # We already have what we need for the super-superclass constructor.
3128 CombinatorialFreeModule
.__init
__(self
,
3135 # Now create the vector space for the algebra, which will have
3136 # its own set of non-ambient coordinates (in terms of the
3138 degree
= sum( f
._matrix
_span
.ambient_vector_space().degree()
3140 V
= VectorSpace(field
, degree
)
3141 vector_basis
= tuple( V(_all2list(b
)) for b
in self
._matrix
_basis
)
3143 # Save the span of our matrix basis (when written out as long
3144 # vectors) because otherwise we'll have to reconstruct it
3145 # every time we want to coerce a matrix into the algebra.
3146 self
._matrix
_span
= V
.span_of_basis( vector_basis
, check
=False)
3148 # Since we don't (re)orthonormalize the basis, the FDEJA
3149 # constructor is going to set self._deortho_matrix to the
3150 # identity matrix. Here we set it to the correct value using
3151 # the deortho matrices from our factors.
3152 self
._deortho
_matrix
= matrix
.block_diagonal(
3153 [J
._deortho
_matrix
for J
in factors
]
3156 self
._inner
_product
_matrix
= matrix
.block_diagonal(
3157 [J
._inner
_product
_matrix
for J
in factors
]
3159 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
3160 self
._inner
_product
_matrix
.set_immutable()
3162 # Building the multiplication table is a bit more tricky
3163 # because we have to embed the entries of the factors'
3164 # multiplication tables into the product EJA.
3166 self
._multiplication
_table
= [ [zed
for j
in range(i
+1)]
3169 # Keep track of an offset that tallies the dimensions of all
3170 # previous factors. If the second factor is dim=2 and if the
3171 # first one is dim=3, then we want to skip the first 3x3 block
3172 # when copying the multiplication table for the second factor.
3175 phi_f
= self
.cartesian_embedding(f
)
3176 factor_dim
= factors
[f
].dimension()
3177 for i
in range(factor_dim
):
3178 for j
in range(i
+1):
3179 f_ij
= factors
[f
]._multiplication
_table
[i
][j
]
3181 self
._multiplication
_table
[offset
+i
][offset
+j
] = e
3182 offset
+= factor_dim
3184 self
.rank
.set_cache(sum(J
.rank() for J
in factors
))
3185 ones
= tuple(J
.one().to_matrix() for J
in factors
)
3186 self
.one
.set_cache(self(ones
))
3188 def _sets_keys(self
):
3193 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
3194 ....: RealSymmetricEJA)
3198 The superclass uses ``_sets_keys()`` to implement its
3199 ``cartesian_factors()`` method::
3201 sage: J1 = RealSymmetricEJA(2,
3203 ....: orthonormalize=False,
3205 sage: J2 = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
3206 sage: J = cartesian_product([J1,J2])
3207 sage: x = sum(i*J.gens()[i] for i in range(len(J.gens())))
3208 sage: x.cartesian_factors()
3209 (a1 + 2*a2, 3*b0 + 4*b1 + 5*b2 + 6*b3)
3212 # Copy/pasted from CombinatorialFreeModule_CartesianProduct,
3213 # but returning a tuple instead of a list.
3214 return tuple(range(len(self
.cartesian_factors())))
3216 def cartesian_factors(self
):
3217 # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
3220 def cartesian_factor(self
, i
):
3222 Return the ``i``th factor of this algebra.
3224 return self
._sets
[i
]
3227 # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
3228 from sage
.categories
.cartesian_product
import cartesian_product
3229 return cartesian_product
.symbol
.join("%s" % factor
3230 for factor
in self
._sets
)
3234 def cartesian_projection(self
, i
):
3238 sage: from mjo.eja.eja_algebra import (random_eja,
3239 ....: JordanSpinEJA,
3241 ....: RealSymmetricEJA,
3242 ....: ComplexHermitianEJA)
3246 The projection morphisms are Euclidean Jordan algebra
3249 sage: J1 = HadamardEJA(2)
3250 sage: J2 = RealSymmetricEJA(2)
3251 sage: J = cartesian_product([J1,J2])
3252 sage: J.cartesian_projection(0)
3253 Linear operator between finite-dimensional Euclidean Jordan
3254 algebras represented by the matrix:
3257 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3258 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3259 Algebraic Real Field
3260 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3262 sage: J.cartesian_projection(1)
3263 Linear operator between finite-dimensional Euclidean Jordan
3264 algebras represented by the matrix:
3268 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3269 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3270 Algebraic Real Field
3271 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3274 The projections work the way you'd expect on the vector
3275 representation of an element::
3277 sage: J1 = JordanSpinEJA(2)
3278 sage: J2 = ComplexHermitianEJA(2)
3279 sage: J = cartesian_product([J1,J2])
3280 sage: pi_left = J.cartesian_projection(0)
3281 sage: pi_right = J.cartesian_projection(1)
3282 sage: pi_left(J.one()).to_vector()
3284 sage: pi_right(J.one()).to_vector()
3286 sage: J.one().to_vector()
3291 The answer never changes::
3293 sage: J1 = random_eja()
3294 sage: J2 = random_eja()
3295 sage: J = cartesian_product([J1,J2])
3296 sage: P0 = J.cartesian_projection(0)
3297 sage: P1 = J.cartesian_projection(0)
3302 offset
= sum( self
.cartesian_factor(k
).dimension()
3304 Ji
= self
.cartesian_factor(i
)
3305 Pi
= self
._module
_morphism
(lambda j
: Ji
.monomial(j
- offset
),
3308 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3311 def cartesian_embedding(self
, i
):
3315 sage: from mjo.eja.eja_algebra import (random_eja,
3316 ....: JordanSpinEJA,
3318 ....: RealSymmetricEJA)
3322 The embedding morphisms are Euclidean Jordan algebra
3325 sage: J1 = HadamardEJA(2)
3326 sage: J2 = RealSymmetricEJA(2)
3327 sage: J = cartesian_product([J1,J2])
3328 sage: J.cartesian_embedding(0)
3329 Linear operator between finite-dimensional Euclidean Jordan
3330 algebras represented by the matrix:
3336 Domain: Euclidean Jordan algebra of dimension 2 over
3337 Algebraic Real Field
3338 Codomain: Euclidean Jordan algebra of dimension 2 over
3339 Algebraic Real Field (+) Euclidean Jordan algebra of
3340 dimension 3 over Algebraic Real Field
3341 sage: J.cartesian_embedding(1)
3342 Linear operator between finite-dimensional Euclidean Jordan
3343 algebras represented by the matrix:
3349 Domain: Euclidean Jordan algebra of dimension 3 over
3350 Algebraic Real Field
3351 Codomain: Euclidean Jordan algebra of dimension 2 over
3352 Algebraic Real Field (+) Euclidean Jordan algebra of
3353 dimension 3 over Algebraic Real Field
3355 The embeddings work the way you'd expect on the vector
3356 representation of an element::
3358 sage: J1 = JordanSpinEJA(3)
3359 sage: J2 = RealSymmetricEJA(2)
3360 sage: J = cartesian_product([J1,J2])
3361 sage: iota_left = J.cartesian_embedding(0)
3362 sage: iota_right = J.cartesian_embedding(1)
3363 sage: iota_left(J1.zero()) == J.zero()
3365 sage: iota_right(J2.zero()) == J.zero()
3367 sage: J1.one().to_vector()
3369 sage: iota_left(J1.one()).to_vector()
3371 sage: J2.one().to_vector()
3373 sage: iota_right(J2.one()).to_vector()
3375 sage: J.one().to_vector()
3380 The answer never changes::
3382 sage: J1 = random_eja()
3383 sage: J2 = random_eja()
3384 sage: J = cartesian_product([J1,J2])
3385 sage: E0 = J.cartesian_embedding(0)
3386 sage: E1 = J.cartesian_embedding(0)
3390 Composing a projection with the corresponding inclusion should
3391 produce the identity map, and mismatching them should produce
3394 sage: J1 = random_eja()
3395 sage: J2 = random_eja()
3396 sage: J = cartesian_product([J1,J2])
3397 sage: iota_left = J.cartesian_embedding(0)
3398 sage: iota_right = J.cartesian_embedding(1)
3399 sage: pi_left = J.cartesian_projection(0)
3400 sage: pi_right = J.cartesian_projection(1)
3401 sage: pi_left*iota_left == J1.one().operator()
3403 sage: pi_right*iota_right == J2.one().operator()
3405 sage: (pi_left*iota_right).is_zero()
3407 sage: (pi_right*iota_left).is_zero()
3411 offset
= sum( self
.cartesian_factor(k
).dimension()
3413 Ji
= self
.cartesian_factor(i
)
3414 Ei
= Ji
._module
_morphism
(lambda j
: self
.monomial(j
+ offset
),
3416 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3419 def subalgebra(self
, basis
, **kwargs
):
3421 Create a subalgebra of this algebra from the given basis.
3423 Only overridden to allow us to use a special Cartesian product
3428 sage: from mjo.eja.eja_algebra import (HadamardEJA,
3429 ....: QuaternionHermitianEJA)
3433 Subalgebras of Cartesian product EJAs have a different class
3434 than those of non-Cartesian-product EJAs::
3436 sage: J1 = HadamardEJA(2,field=QQ,orthonormalize=False)
3437 sage: J2 = QuaternionHermitianEJA(0,field=QQ,orthonormalize=False)
3438 sage: J = cartesian_product([J1,J2])
3439 sage: K1 = J1.subalgebra((J1.one(),), orthonormalize=False)
3440 sage: K = J.subalgebra((J.one(),), orthonormalize=False)
3441 sage: K1.__class__ is K.__class__
3445 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalCartesianProductEJASubalgebra
3446 return FiniteDimensionalCartesianProductEJASubalgebra(self
, basis
, **kwargs
)
3448 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3450 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3453 A separate class for products of algebras for which we know a
3458 sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA,
3460 ....: JordanSpinEJA,
3461 ....: RealSymmetricEJA)
3465 This gives us fast characteristic polynomial computations in
3466 product algebras, too::
3469 sage: J1 = JordanSpinEJA(2)
3470 sage: J2 = RealSymmetricEJA(3)
3471 sage: J = cartesian_product([J1,J2])
3472 sage: J.characteristic_polynomial_of().degree()
3479 The ``cartesian_product()`` function only uses the first factor to
3480 decide where the result will live; thus we have to be careful to
3481 check that all factors do indeed have a ``rational_algebra()`` method
3482 before we construct an algebra that claims to have a rational basis::
3484 sage: J1 = HadamardEJA(2)
3485 sage: jp = lambda X,Y: X*Y
3486 sage: ip = lambda X,Y: X[0,0]*Y[0,0]
3487 sage: b1 = matrix(QQ, [[1]])
3488 sage: J2 = FiniteDimensionalEJA((b1,), jp, ip)
3489 sage: cartesian_product([J2,J1]) # factor one not RationalBasisEJA
3490 Euclidean Jordan algebra of dimension 1 over Algebraic Real
3491 Field (+) Euclidean Jordan algebra of dimension 2 over Algebraic
3493 sage: cartesian_product([J1,J2]) # factor one is RationalBasisEJA
3494 Traceback (most recent call last):
3496 ValueError: factor not a RationalBasisEJA
3499 def __init__(self
, algebras
, **kwargs
):
3500 if not all( hasattr(r
, "rational_algebra") for r
in algebras
):
3501 raise ValueError("factor not a RationalBasisEJA")
3503 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3506 def rational_algebra(self
):
3507 if self
.base_ring() is QQ
:
3510 return cartesian_product([
3511 r
.rational_algebra() for r
in self
.cartesian_factors()
3515 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3517 def random_eja(max_dimension
=None, *args
, **kwargs
):
3522 sage: from mjo.eja.eja_algebra import random_eja
3526 sage: n = ZZ.random_element(1,5)
3527 sage: J = random_eja(max_dimension=n, field=QQ, orthonormalize=False)
3528 sage: J.dimension() <= n
3532 # Use the ConcreteEJA default as the total upper bound (regardless
3533 # of any whether or not any individual factors set a lower limit).
3534 if max_dimension
is None:
3535 max_dimension
= ConcreteEJA
._max
_random
_instance
_dimension
()
3536 J1
= ConcreteEJA
.random_instance(max_dimension
, *args
, **kwargs
)
3539 # Roll the dice to see if we attempt a Cartesian product.
3540 dice_roll
= ZZ
.random_element(len(ConcreteEJA
.__subclasses
__()) + 1)
3541 new_max_dimension
= max_dimension
- J1
.dimension()
3542 if new_max_dimension
== 0 or dice_roll
!= 0:
3543 # If it's already as big as we're willing to tolerate, just
3544 # return it and don't worry about Cartesian products.
3547 # Use random_eja() again so we can get more than two factors
3548 # if the sub-call also Decides on a cartesian product.
3549 J2
= random_eja(new_max_dimension
, *args
, **kwargs
)
3550 return cartesian_product([J1
,J2
])
3553 class ComplexSkewSymmetricEJA(RationalBasisEJA
, ConcreteEJA
):
3555 The skew-symmetric EJA of order `2n` described in Faraut and
3556 Koranyi's Exercise III.1.b. It has dimension `2n^2 - n`.
3558 It is (not obviously) isomorphic to the QuaternionHermitianEJA of
3559 order `n`, as can be inferred by comparing rank/dimension or
3560 explicitly from their "characteristic polynomial of" functions,
3561 which just so happen to align nicely.
3565 sage: from mjo.eja.eja_algebra import (ComplexSkewSymmetricEJA,
3566 ....: QuaternionHermitianEJA)
3567 sage: from mjo.eja.eja_operator import FiniteDimensionalEJAOperator
3571 This EJA is isomorphic to the quaternions::
3573 sage: J = ComplexSkewSymmetricEJA(2, field=QQ, orthonormalize=False)
3574 sage: K = QuaternionHermitianEJA(2, field=QQ, orthonormalize=False)
3575 sage: jordan_isom_matrix = matrix.diagonal(QQ,[-1,1,1,1,1,-1])
3576 sage: phi = FiniteDimensionalEJAOperator(J,K,jordan_isom_matrix)
3577 sage: all( phi(x*y) == phi(x)*phi(y)
3578 ....: for x in J.gens()
3579 ....: for y in J.gens() )
3581 sage: x,y = J.random_elements(2)
3582 sage: phi(x*y) == phi(x)*phi(y)
3587 Random elements should satisfy the same conditions that the basis
3590 sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ,
3591 ....: orthonormalize=False)
3592 sage: x,y = K.random_elements(2)
3594 sage: x = x.to_matrix()
3595 sage: y = y.to_matrix()
3596 sage: z = z.to_matrix()
3597 sage: all( e.is_skew_symmetric() for e in (x,y,z) )
3599 sage: J = -K.one().to_matrix()
3600 sage: all( e*J == J*e.conjugate() for e in (x,y,z) )
3603 The power law in Faraut & Koranyi's II.7.a is satisfied.
3604 We're in a subalgebra of theirs, but powers are still
3607 sage: K = ComplexSkewSymmetricEJA.random_instance(field=QQ,
3608 ....: orthonormalize=False)
3609 sage: x = K.random_element()
3610 sage: k = ZZ.random_element(5)
3612 sage: J = -K.one().to_matrix()
3613 sage: expected = K(-J*(J*x.to_matrix())^k)
3614 sage: actual == expected
3619 def _max_random_instance_size(max_dimension
):
3620 # Obtained by solving d = 2n^2 - n, which comes from noticing
3621 # that, in 2x2 block form, any element of this algebra has a
3622 # free skew-symmetric top-left block, a Hermitian top-right
3623 # block, and two bottom blocks that are determined by the top.
3624 # The ZZ-int-ZZ thing is just "floor."
3625 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/4 + 1/4))
3628 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
3630 Return a random instance of this type of algebra.
3632 class_max_d
= cls
._max
_random
_instance
_dimension
()
3633 if (max_dimension
is None or max_dimension
> class_max_d
):
3634 max_dimension
= class_max_d
3635 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
3636 n
= ZZ
.random_element(max_size
+ 1)
3637 return cls(n
, **kwargs
)
3640 def _denormalized_basis(A
):
3644 sage: from mjo.hurwitz import ComplexMatrixAlgebra
3645 sage: from mjo.eja.eja_algebra import ComplexSkewSymmetricEJA
3649 The basis elements are all skew-Hermitian::
3651 sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension()
3652 sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max)
3653 sage: n = ZZ.random_element(n_max + 1)
3654 sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ)
3655 sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A)
3656 sage: all( M.is_skew_symmetric() for M in B)
3659 The basis elements ``b`` all satisfy ``b*J == J*b.conjugate()``,
3660 as in the definition of the algebra::
3662 sage: d_max = ComplexSkewSymmetricEJA._max_random_instance_dimension()
3663 sage: n_max = ComplexSkewSymmetricEJA._max_random_instance_size(d_max)
3664 sage: n = ZZ.random_element(n_max + 1)
3665 sage: A = ComplexMatrixAlgebra(2*n, scalars=QQ)
3666 sage: I_n = matrix.identity(ZZ, n)
3667 sage: J = matrix.block(ZZ, 2, 2, (0, I_n, -I_n, 0), subdivide=False)
3668 sage: J = A.from_list(J.rows())
3669 sage: B = ComplexSkewSymmetricEJA._denormalized_basis(A)
3670 sage: all( b*J == J*b.conjugate() for b in B )
3674 es
= A
.entry_algebra_gens()
3675 gen
= lambda A
,m
: A
.monomial(m
)
3679 # The size of the blocks. We're going to treat these thing as
3680 # 2x2 block matrices,
3683 # [ -x2-conj x1-conj ]
3685 # where x1 is skew-symmetric and x2 is Hermitian.
3689 # We only loop through the top half of the matrix, because the
3690 # bottom can be constructed from the top.
3692 # First do the top-left block, which is skew-symmetric.
3693 # We can compute the bottom-right block in the process.
3694 for j
in range(i
+1):
3696 # Skew-symmetry implies zeros for (i == j).
3698 # Top-left block's entry.
3699 E_ij
= gen(A
, (i
,j
,e
))
3700 E_ij
-= gen(A
, (j
,i
,e
))
3702 # Bottom-right block's entry.
3703 F_ij
= gen(A
, (i
+m
,j
+m
,e
)).conjugate()
3704 F_ij
-= gen(A
, (j
+m
,i
+m
,e
)).conjugate()
3706 basis
.append(E_ij
+ F_ij
)
3708 # Now do the top-right block, which is Hermitian, and compute
3709 # the bottom-left block along the way.
3710 for j
in range(m
,i
+m
+1):
3712 # Hermitian matrices have real diagonal entries.
3713 # Top-right block's entry.
3714 E_ii
= gen(A
, (i
,j
,es
[0]))
3716 # Bottom-left block's entry. Don't conjugate
3718 E_ii
-= gen(A
, (i
+m
,j
-m
,es
[0]))
3722 # Top-right block's entry. BEWARE! We're not
3723 # reflecting across the main diagonal as in
3724 # (i,j)~(j,i). We're only reflecting across
3725 # the diagonal for the top-right block.
3726 E_ij
= gen(A
, (i
,j
,e
))
3728 # Shift it back to non-offset coords, transpose,
3729 # conjugate, and put it back:
3731 # (i,j) -> (i,j-m) -> (j-m, i) -> (j-m, i+m)
3732 E_ij
+= gen(A
, (j
-m
,i
+m
,e
)).conjugate()
3734 # Bottom-left's block's below-diagonal entry.
3735 # Just shift the top-right coords down m and
3737 F_ij
= -gen(A
, (i
+m
,j
-m
,e
)).conjugate()
3738 F_ij
+= -gen(A
, (j
,i
,e
)) # double-conjugate cancels
3740 basis
.append(E_ij
+ F_ij
)
3742 return tuple( basis
)
3746 def _J_matrix(matrix_space
):
3747 n
= matrix_space
.nrows() // 2
3748 F
= matrix_space
.base_ring()
3749 I_n
= matrix
.identity(F
, n
)
3750 J
= matrix
.block(F
, 2, 2, (0, I_n
, -I_n
, 0), subdivide
=False)
3751 return matrix_space
.from_list(J
.rows())
3754 return ComplexSkewSymmetricEJA
._J
_matrix
(self
.matrix_space())
3756 def __init__(self
, n
, field
=AA
, **kwargs
):
3757 # New code; always check the axioms.
3758 #if "check_axioms" not in kwargs: kwargs["check_axioms"] = False
3760 from mjo
.hurwitz
import ComplexMatrixAlgebra
3761 A
= ComplexMatrixAlgebra(2*n
, scalars
=field
)
3762 J
= ComplexSkewSymmetricEJA
._J
_matrix
(A
)
3764 def jordan_product(X
,Y
):
3765 return (X
*J
*Y
+ Y
*J
*X
)/2
3767 def inner_product(X
,Y
):
3768 return (X
.conjugate_transpose()*Y
).trace().real()
3770 super().__init
__(self
._denormalized
_basis
(A
),
3777 # This algebra is conjectured (by me) to be isomorphic to
3778 # the quaternion Hermitian EJA of size n, and the rank
3779 # would follow from that.
3780 #self.rank.set_cache(n)
3781 self
.one
.set_cache( self(-J
) )