2 Representations and constructions for Euclidean Jordan algebras.
4 A Euclidean Jordan algebra is a Jordan algebra that has some
7 1. It is finite-dimensional.
8 2. Its scalar field is the real numbers.
9 3a. An inner product is defined on it, and...
10 3b. That inner product is compatible with the Jordan product
11 in the sense that `<x*y,z> = <y,x*z>` for all elements
12 `x,y,z` in the algebra.
14 Every Euclidean Jordan algebra is formally-real: for any two elements
15 `x` and `y` in the algebra, `x^{2} + y^{2} = 0` implies that `x = y =
16 0`. Conversely, every finite-dimensional formally-real Jordan algebra
17 can be made into a Euclidean Jordan algebra with an appropriate choice
20 Formally-real Jordan algebras were originally studied as a framework
21 for quantum mechanics. Today, Euclidean Jordan algebras are crucial in
22 symmetric cone optimization, since every symmetric cone arises as the
23 cone of squares in some Euclidean Jordan algebra.
25 It is known that every Euclidean Jordan algebra decomposes into an
26 orthogonal direct sum (essentially, a Cartesian product) of simple
27 algebras, and that moreover, up to Jordan-algebra isomorphism, there
28 are only five families of simple algebras. We provide constructions
29 for these simple algebras:
31 * :class:`BilinearFormEJA`
32 * :class:`RealSymmetricEJA`
33 * :class:`ComplexHermitianEJA`
34 * :class:`QuaternionHermitianEJA`
35 * :class:`OctonionHermitianEJA`
37 In addition to these, we provide two other example constructions,
39 * :class:`JordanSpinEJA`
40 * :class:`HadamardEJA`
44 The Jordan spin algebra is a bilinear form algebra where the bilinear
45 form is the identity. The Hadamard EJA is simply a Cartesian product
46 of one-dimensional spin algebras. The Albert EJA is simply a special
47 case of the :class:`OctonionHermitianEJA` where the matrices are
48 three-by-three and the resulting space has dimension 27. And
49 last/least, the trivial EJA is exactly what you think it is; it could
50 also be obtained by constructing a dimension-zero instance of any of
51 the other algebras. Cartesian products of these are also supported
52 using the usual ``cartesian_product()`` function; as a result, we
53 support (up to isomorphism) all Euclidean Jordan algebras.
55 At a minimum, the following are required to construct a Euclidean
58 * A basis of matrices, column vectors, or MatrixAlgebra elements
59 * A Jordan product defined on the basis
60 * Its inner product defined on the basis
62 The real numbers form a Euclidean Jordan algebra when both the Jordan
63 and inner products are the usual multiplication. We use this as our
64 example, and demonstrate a few ways to construct an EJA.
66 First, we can use one-by-one SageMath matrices with algebraic real
67 entries to represent real numbers. We define the Jordan and inner
68 products to be essentially real-number multiplication, with the only
69 difference being that the Jordan product again returns a one-by-one
70 matrix, whereas the inner product must return a scalar. Our basis for
71 the one-by-one matrices is of course the set consisting of a single
72 matrix with its sole entry non-zero::
74 sage: from mjo.eja.eja_algebra import FiniteDimensionalEJA
75 sage: jp = lambda X,Y: X*Y
76 sage: ip = lambda X,Y: X[0,0]*Y[0,0]
77 sage: b1 = matrix(AA, [[1]])
78 sage: J1 = FiniteDimensionalEJA((b1,), jp, ip)
80 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
82 In fact, any positive scalar multiple of that inner-product would work::
84 sage: ip2 = lambda X,Y: 16*ip(X,Y)
85 sage: J2 = FiniteDimensionalEJA((b1,), jp, ip2)
87 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
89 But beware that your basis will be orthonormalized _with respect to the
90 given inner-product_ unless you pass ``orthonormalize=False`` to the
91 constructor. For example::
93 sage: J3 = FiniteDimensionalEJA((b1,), jp, ip2, orthonormalize=False)
95 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
97 To see the difference, you can take the first and only basis element
98 of the resulting algebra, and ask for it to be converted back into
101 sage: J1.basis()[0].to_matrix()
103 sage: J2.basis()[0].to_matrix()
105 sage: J3.basis()[0].to_matrix()
108 Since square roots are used in that process, the default scalar field
109 that we use is the field of algebraic real numbers, ``AA``. You can
110 also Use rational numbers, but only if you either pass
111 ``orthonormalize=False`` or know that orthonormalizing your basis
112 won't stray beyond the rational numbers. The example above would
113 have worked only because ``sqrt(16) == 4`` is rational.
115 Another option for your basis is to use elemebts of a
116 :class:`MatrixAlgebra`::
118 sage: from mjo.matrix_algebra import MatrixAlgebra
119 sage: A = MatrixAlgebra(1,AA,AA)
120 sage: J4 = FiniteDimensionalEJA(A.gens(), jp, ip)
122 Euclidean Jordan algebra of dimension 1 over Algebraic Real Field
123 sage: J4.basis()[0].to_matrix()
128 An easier way to view the entire EJA basis in its original (but
129 perhaps orthonormalized) matrix form is to use the ``matrix_basis``
132 sage: J4.matrix_basis()
137 In particular, a :class:`MatrixAlgebra` is needed to work around the
138 fact that matrices in SageMath must have entries in the same
139 (commutative and associative) ring as its scalars. There are many
140 Euclidean Jordan algebras whose elements are matrices that violate
141 those assumptions. The complex, quaternion, and octonion Hermitian
142 matrices all have entries in a ring (the complex numbers, quaternions,
143 or octonions...) that differs from the algebra's scalar ring (the real
144 numbers). Quaternions are also non-commutative; the octonions are
145 neither commutative nor associative.
149 sage: from mjo.eja.eja_algebra import random_eja
154 Euclidean Jordan algebra of dimension...
157 from sage
.algebras
.quatalg
.quaternion_algebra
import QuaternionAlgebra
158 from sage
.categories
.magmatic_algebras
import MagmaticAlgebras
159 from sage
.categories
.sets_cat
import cartesian_product
160 from sage
.combinat
.free_module
import CombinatorialFreeModule
161 from sage
.matrix
.constructor
import matrix
162 from sage
.matrix
.matrix_space
import MatrixSpace
163 from sage
.misc
.cachefunc
import cached_method
164 from sage
.misc
.table
import table
165 from sage
.modules
.free_module
import FreeModule
, VectorSpace
166 from sage
.rings
.all
import (ZZ
, QQ
, AA
, QQbar
, RR
, RLF
, CLF
,
169 from mjo
.eja
.eja_element
import FiniteDimensionalEJAElement
170 from mjo
.eja
.eja_operator
import FiniteDimensionalEJAOperator
171 from mjo
.eja
.eja_utils
import _all2list
173 def EuclideanJordanAlgebras(field
):
175 The category of Euclidean Jordan algebras over ``field``, which
176 must be a subfield of the real numbers. For now this is just a
177 convenient wrapper around all of the other category axioms that
180 category
= MagmaticAlgebras(field
).FiniteDimensional()
181 category
= category
.WithBasis().Unital().Commutative()
184 class FiniteDimensionalEJA(CombinatorialFreeModule
):
186 A finite-dimensional Euclidean Jordan algebra.
190 - ``basis`` -- a tuple; a tuple of basis elements in "matrix
191 form," which must be the same form as the arguments to
192 ``jordan_product`` and ``inner_product``. In reality, "matrix
193 form" can be either vectors, matrices, or a Cartesian product
194 (ordered tuple) of vectors or matrices. All of these would
195 ideally be vector spaces in sage with no special-casing
196 needed; but in reality we turn vectors into column-matrices
197 and Cartesian products `(a,b)` into column matrices
198 `(a,b)^{T}` after converting `a` and `b` themselves.
200 - ``jordan_product`` -- a function; afunction of two ``basis``
201 elements (in matrix form) that returns their jordan product,
202 also in matrix form; this will be applied to ``basis`` to
203 compute a multiplication table for the algebra.
205 - ``inner_product`` -- a function; a function of two ``basis``
206 elements (in matrix form) that returns their inner
207 product. This will be applied to ``basis`` to compute an
208 inner-product table (basically a matrix) for this algebra.
210 - ``matrix_space`` -- the space that your matrix basis lives in,
211 or ``None`` (the default). So long as your basis does not have
212 length zero you can omit this. But in trivial algebras, it is
215 - ``field`` -- a subfield of the reals (default: ``AA``); the scalar
216 field for the algebra.
218 - ``orthonormalize`` -- boolean (default: ``True``); whether or
219 not to orthonormalize the basis. Doing so is expensive and
220 generally rules out using the rationals as your ``field``, but
221 is required for spectral decompositions.
225 sage: from mjo.eja.eja_algebra import random_eja
229 We should compute that an element subalgebra is associative even
230 if we circumvent the element method::
232 sage: set_random_seed()
233 sage: J = random_eja(field=QQ,orthonormalize=False)
234 sage: x = J.random_element()
235 sage: A = x.subalgebra_generated_by(orthonormalize=False)
236 sage: basis = tuple(b.superalgebra_element() for b in A.basis())
237 sage: J.subalgebra(basis, orthonormalize=False).is_associative()
240 Element
= FiniteDimensionalEJAElement
243 def _check_input_field(field
):
244 if not field
.is_subring(RR
):
245 # Note: this does return true for the real algebraic
246 # field, the rationals, and any quadratic field where
247 # we've specified a real embedding.
248 raise ValueError("scalar field is not real")
251 def _check_input_axioms(basis
, jordan_product
, inner_product
):
252 if not all( jordan_product(bi
,bj
) == jordan_product(bj
,bi
)
255 raise ValueError("Jordan product is not commutative")
257 if not all( inner_product(bi
,bj
) == inner_product(bj
,bi
)
260 raise ValueError("inner-product is not commutative")
277 self
._check
_input
_field
(field
)
280 # Check commutativity of the Jordan and inner-products.
281 # This has to be done before we build the multiplication
282 # and inner-product tables/matrices, because we take
283 # advantage of symmetry in the process.
284 self
._check
_input
_axioms
(basis
, jordan_product
, inner_product
)
287 # All zero- and one-dimensional algebras are just the real
288 # numbers with (some positive multiples of) the usual
289 # multiplication as its Jordan and inner-product.
291 if associative
is None:
292 # We should figure it out. As with check_axioms, we have to do
293 # this without the help of the _jordan_product_is_associative()
294 # method because we need to know the category before we
295 # initialize the algebra.
296 associative
= all( jordan_product(jordan_product(bi
,bj
),bk
)
298 jordan_product(bi
,jordan_product(bj
,bk
))
303 category
= EuclideanJordanAlgebras(field
)
306 # Element subalgebras can take advantage of this.
307 category
= category
.Associative()
309 # Call the superclass constructor so that we can use its from_vector()
310 # method to build our multiplication table.
311 CombinatorialFreeModule
.__init
__(self
,
318 # Now comes all of the hard work. We'll be constructing an
319 # ambient vector space V that our (vectorized) basis lives in,
320 # as well as a subspace W of V spanned by those (vectorized)
321 # basis elements. The W-coordinates are the coefficients that
322 # we see in things like x = 1*b1 + 2*b2.
326 degree
= len(_all2list(basis
[0]))
328 # Build an ambient space that fits our matrix basis when
329 # written out as "long vectors."
330 V
= VectorSpace(field
, degree
)
332 # The matrix that will hold the orthonormal -> unorthonormal
333 # coordinate transformation. Default to an identity matrix of
334 # the appropriate size to avoid special cases for None
336 self
._deortho
_matrix
= matrix
.identity(field
,n
)
339 # Save a copy of the un-orthonormalized basis for later.
340 # Convert it to ambient V (vector) coordinates while we're
341 # at it, because we'd have to do it later anyway.
342 deortho_vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
344 from mjo
.eja
.eja_utils
import gram_schmidt
345 basis
= tuple(gram_schmidt(basis
, inner_product
))
347 # Save the (possibly orthonormalized) matrix basis for
348 # later, as well as the space that its elements live in.
349 # In most cases we can deduce the matrix space, but when
350 # n == 0 (that is, there are no basis elements) we cannot.
351 self
._matrix
_basis
= basis
352 if matrix_space
is None:
353 self
._matrix
_space
= self
._matrix
_basis
[0].parent()
355 self
._matrix
_space
= matrix_space
357 # Now create the vector space for the algebra, which will have
358 # its own set of non-ambient coordinates (in terms of the
360 vector_basis
= tuple( V(_all2list(b
)) for b
in basis
)
362 # Save the span of our matrix basis (when written out as long
363 # vectors) because otherwise we'll have to reconstruct it
364 # every time we want to coerce a matrix into the algebra.
365 self
._matrix
_span
= V
.span_of_basis( vector_basis
, check
=check_axioms
)
368 # Now "self._matrix_span" is the vector space of our
369 # algebra coordinates. The variables "X1", "X2",... refer
370 # to the entries of vectors in self._matrix_span. Thus to
371 # convert back and forth between the orthonormal
372 # coordinates and the given ones, we need to stick the
373 # original basis in self._matrix_span.
374 U
= V
.span_of_basis( deortho_vector_basis
, check
=check_axioms
)
375 self
._deortho
_matrix
= matrix
.column( U
.coordinate_vector(q
)
376 for q
in vector_basis
)
379 # Now we actually compute the multiplication and inner-product
380 # tables/matrices using the possibly-orthonormalized basis.
381 self
._inner
_product
_matrix
= matrix
.identity(field
, n
)
383 self
._multiplication
_table
= [ [zed
for j
in range(i
+1)]
386 # Note: the Jordan and inner-products are defined in terms
387 # of the ambient basis. It's important that their arguments
388 # are in ambient coordinates as well.
391 # ortho basis w.r.t. ambient coords
395 # The jordan product returns a matrixy answer, so we
396 # have to convert it to the algebra coordinates.
397 elt
= jordan_product(q_i
, q_j
)
398 elt
= self
._matrix
_span
.coordinate_vector(V(_all2list(elt
)))
399 self
._multiplication
_table
[i
][j
] = self
.from_vector(elt
)
401 if not orthonormalize
:
402 # If we're orthonormalizing the basis with respect
403 # to an inner-product, then the inner-product
404 # matrix with respect to the resulting basis is
405 # just going to be the identity.
406 ip
= inner_product(q_i
, q_j
)
407 self
._inner
_product
_matrix
[i
,j
] = ip
408 self
._inner
_product
_matrix
[j
,i
] = ip
410 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
411 self
._inner
_product
_matrix
.set_immutable()
414 if not self
._is
_jordanian
():
415 raise ValueError("Jordan identity does not hold")
416 if not self
._inner
_product
_is
_associative
():
417 raise ValueError("inner product is not associative")
420 def _coerce_map_from_base_ring(self
):
422 Disable the map from the base ring into the algebra.
424 Performing a nonsense conversion like this automatically
425 is counterpedagogical. The fallback is to try the usual
426 element constructor, which should also fail.
430 sage: from mjo.eja.eja_algebra import random_eja
434 sage: set_random_seed()
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: set_random_seed()
460 sage: J = random_eja()
461 sage: n = J.dimension()
464 sage: bi_bj = J.zero()*J.zero()
466 ....: i = ZZ.random_element(n)
467 ....: j = ZZ.random_element(n)
468 ....: bi = J.monomial(i)
469 ....: bj = J.monomial(j)
470 ....: bi_bj = J.product_on_basis(i,j)
475 # We only stored the lower-triangular portion of the
476 # multiplication table.
478 return self
._multiplication
_table
[i
][j
]
480 return self
._multiplication
_table
[j
][i
]
482 def inner_product(self
, x
, y
):
484 The inner product associated with this Euclidean Jordan algebra.
486 Defaults to the trace inner product, but can be overridden by
487 subclasses if they are sure that the necessary properties are
492 sage: from mjo.eja.eja_algebra import (random_eja,
494 ....: BilinearFormEJA)
498 Our inner product is "associative," which means the following for
499 a symmetric bilinear form::
501 sage: set_random_seed()
502 sage: J = random_eja()
503 sage: x,y,z = J.random_elements(3)
504 sage: (x*y).inner_product(z) == y.inner_product(x*z)
509 Ensure that this is the usual inner product for the algebras
512 sage: set_random_seed()
513 sage: J = HadamardEJA.random_instance()
514 sage: x,y = J.random_elements(2)
515 sage: actual = x.inner_product(y)
516 sage: expected = x.to_vector().inner_product(y.to_vector())
517 sage: actual == expected
520 Ensure that this is one-half of the trace inner-product in a
521 BilinearFormEJA that isn't just the reals (when ``n`` isn't
522 one). This is in Faraut and Koranyi, and also my "On the
525 sage: set_random_seed()
526 sage: J = BilinearFormEJA.random_instance()
527 sage: n = J.dimension()
528 sage: x = J.random_element()
529 sage: y = J.random_element()
530 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
534 B
= self
._inner
_product
_matrix
535 return (B
*x
.to_vector()).inner_product(y
.to_vector())
538 def is_associative(self
):
540 Return whether or not this algebra's Jordan product is associative.
544 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
548 sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False)
549 sage: J.is_associative()
551 sage: x = sum(J.gens())
552 sage: A = x.subalgebra_generated_by(orthonormalize=False)
553 sage: A.is_associative()
557 return "Associative" in self
.category().axioms()
559 def _is_commutative(self
):
561 Whether or not this algebra's multiplication table is commutative.
563 This method should of course always return ``True``, unless
564 this algebra was constructed with ``check_axioms=False`` and
565 passed an invalid multiplication table.
567 return all( x
*y
== y
*x
for x
in self
.gens() for y
in self
.gens() )
569 def _is_jordanian(self
):
571 Whether or not this algebra's multiplication table respects the
572 Jordan identity `(x^{2})(xy) = x(x^{2}y)`.
574 We only check one arrangement of `x` and `y`, so for a
575 ``True`` result to be truly true, you should also check
576 :meth:`_is_commutative`. This method should of course always
577 return ``True``, unless this algebra was constructed with
578 ``check_axioms=False`` and passed an invalid multiplication table.
580 return all( (self
.monomial(i
)**2)*(self
.monomial(i
)*self
.monomial(j
))
582 (self
.monomial(i
))*((self
.monomial(i
)**2)*self
.monomial(j
))
583 for i
in range(self
.dimension())
584 for j
in range(self
.dimension()) )
586 def _jordan_product_is_associative(self
):
588 Return whether or not this algebra's Jordan product is
589 associative; that is, whether or not `x*(y*z) = (x*y)*z`
592 This method should agree with :meth:`is_associative` unless
593 you lied about the value of the ``associative`` parameter
594 when you constructed the algebra.
598 sage: from mjo.eja.eja_algebra import (random_eja,
599 ....: RealSymmetricEJA,
600 ....: ComplexHermitianEJA,
601 ....: QuaternionHermitianEJA)
605 sage: J = RealSymmetricEJA(4, orthonormalize=False)
606 sage: J._jordan_product_is_associative()
608 sage: x = sum(J.gens())
609 sage: A = x.subalgebra_generated_by()
610 sage: A._jordan_product_is_associative()
615 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
616 sage: J._jordan_product_is_associative()
618 sage: x = sum(J.gens())
619 sage: A = x.subalgebra_generated_by(orthonormalize=False)
620 sage: A._jordan_product_is_associative()
625 sage: J = QuaternionHermitianEJA(2)
626 sage: J._jordan_product_is_associative()
628 sage: x = sum(J.gens())
629 sage: A = x.subalgebra_generated_by()
630 sage: A._jordan_product_is_associative()
635 The values we've presupplied to the constructors agree with
638 sage: set_random_seed()
639 sage: J = random_eja()
640 sage: J.is_associative() == J._jordan_product_is_associative()
646 # Used to check whether or not something is zero.
649 # I don't know of any examples that make this magnitude
650 # necessary because I don't know how to make an
651 # associative algebra when the element subalgebra
652 # construction is unreliable (as it is over RDF; we can't
653 # find the degree of an element because we can't compute
654 # the rank of a matrix). But even multiplication of floats
655 # is non-associative, so *some* epsilon is needed... let's
656 # just take the one from _inner_product_is_associative?
659 for i
in range(self
.dimension()):
660 for j
in range(self
.dimension()):
661 for k
in range(self
.dimension()):
665 diff
= (x
*y
)*z
- x
*(y
*z
)
667 if diff
.norm() > epsilon
:
672 def _inner_product_is_associative(self
):
674 Return whether or not this algebra's inner product `B` is
675 associative; that is, whether or not `B(xy,z) = B(x,yz)`.
677 This method should of course always return ``True``, unless
678 this algebra was constructed with ``check_axioms=False`` and
679 passed an invalid Jordan or inner-product.
683 # Used to check whether or not something is zero.
686 # This choice is sufficient to allow the construction of
687 # QuaternionHermitianEJA(2, field=RDF) with check_axioms=True.
690 for i
in range(self
.dimension()):
691 for j
in range(self
.dimension()):
692 for k
in range(self
.dimension()):
696 diff
= (x
*y
).inner_product(z
) - x
.inner_product(y
*z
)
698 if diff
.abs() > epsilon
:
703 def _element_constructor_(self
, elt
):
705 Construct an element of this algebra or a subalgebra from its
706 EJA element, vector, or matrix representation.
708 This gets called only after the parent element _call_ method
709 fails to find a coercion for the argument.
713 sage: from mjo.eja.eja_algebra import (random_eja,
716 ....: RealSymmetricEJA)
720 The identity in `S^n` is converted to the identity in the EJA::
722 sage: J = RealSymmetricEJA(3)
723 sage: I = matrix.identity(QQ,3)
724 sage: J(I) == J.one()
727 This skew-symmetric matrix can't be represented in the EJA::
729 sage: J = RealSymmetricEJA(3)
730 sage: A = matrix(QQ,3, lambda i,j: i-j)
732 Traceback (most recent call last):
734 ValueError: not an element of this algebra
736 Tuples work as well, provided that the matrix basis for the
737 algebra consists of them::
739 sage: J1 = HadamardEJA(3)
740 sage: J2 = RealSymmetricEJA(2)
741 sage: J = cartesian_product([J1,J2])
742 sage: J( (J1.matrix_basis()[1], J2.matrix_basis()[2]) )
745 Subalgebra elements are embedded into the superalgebra::
747 sage: J = JordanSpinEJA(3)
750 sage: x = sum(J.gens())
751 sage: A = x.subalgebra_generated_by()
757 Ensure that we can convert any element back and forth
758 faithfully between its matrix and algebra representations::
760 sage: set_random_seed()
761 sage: J = random_eja()
762 sage: x = J.random_element()
763 sage: J(x.to_matrix()) == x
766 We cannot coerce elements between algebras just because their
767 matrix representations are compatible::
769 sage: J1 = HadamardEJA(3)
770 sage: J2 = JordanSpinEJA(3)
772 Traceback (most recent call last):
774 ValueError: not an element of this algebra
776 Traceback (most recent call last):
778 ValueError: not an element of this algebra
781 msg
= "not an element of this algebra"
782 if elt
in self
.base_ring():
783 # Ensure that no base ring -> algebra coercion is performed
784 # by this method. There's some stupidity in sage that would
785 # otherwise propagate to this method; for example, sage thinks
786 # that the integer 3 belongs to the space of 2-by-2 matrices.
787 raise ValueError(msg
)
789 if hasattr(elt
, 'superalgebra_element'):
790 # Handle subalgebra elements
791 if elt
.parent().superalgebra() == self
:
792 return elt
.superalgebra_element()
794 if hasattr(elt
, 'sparse_vector'):
795 # Convert a vector into a column-matrix. We check for
796 # "sparse_vector" and not "column" because matrices also
797 # have a "column" method.
800 if elt
not in self
.matrix_space():
801 raise ValueError(msg
)
803 # Thanks for nothing! Matrix spaces aren't vector spaces in
804 # Sage, so we have to figure out its matrix-basis coordinates
805 # ourselves. We use the basis space's ring instead of the
806 # element's ring because the basis space might be an algebraic
807 # closure whereas the base ring of the 3-by-3 identity matrix
808 # could be QQ instead of QQbar.
810 # And, we also have to handle Cartesian product bases (when
811 # the matrix basis consists of tuples) here. The "good news"
812 # is that we're already converting everything to long vectors,
813 # and that strategy works for tuples as well.
815 elt
= self
._matrix
_span
.ambient_vector_space()(_all2list(elt
))
818 coords
= self
._matrix
_span
.coordinate_vector(elt
)
819 except ArithmeticError: # vector is not in free module
820 raise ValueError(msg
)
822 return self
.from_vector(coords
)
826 Return a string representation of ``self``.
830 sage: from mjo.eja.eja_algebra import JordanSpinEJA
834 Ensure that it says what we think it says::
836 sage: JordanSpinEJA(2, field=AA)
837 Euclidean Jordan algebra of dimension 2 over Algebraic Real Field
838 sage: JordanSpinEJA(3, field=RDF)
839 Euclidean Jordan algebra of dimension 3 over Real Double Field
842 fmt
= "Euclidean Jordan algebra of dimension {} over {}"
843 return fmt
.format(self
.dimension(), self
.base_ring())
847 def characteristic_polynomial_of(self
):
849 Return the algebra's "characteristic polynomial of" function,
850 which is itself a multivariate polynomial that, when evaluated
851 at the coordinates of some algebra element, returns that
852 element's characteristic polynomial.
854 The resulting polynomial has `n+1` variables, where `n` is the
855 dimension of this algebra. The first `n` variables correspond to
856 the coordinates of an algebra element: when evaluated at the
857 coordinates of an algebra element with respect to a certain
858 basis, the result is a univariate polynomial (in the one
859 remaining variable ``t``), namely the characteristic polynomial
864 sage: from mjo.eja.eja_algebra import JordanSpinEJA, TrivialEJA
868 The characteristic polynomial in the spin algebra is given in
869 Alizadeh, Example 11.11::
871 sage: J = JordanSpinEJA(3)
872 sage: p = J.characteristic_polynomial_of(); p
873 X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2
874 sage: xvec = J.one().to_vector()
878 By definition, the characteristic polynomial is a monic
879 degree-zero polynomial in a rank-zero algebra. Note that
880 Cayley-Hamilton is indeed satisfied since the polynomial
881 ``1`` evaluates to the identity element of the algebra on
884 sage: J = TrivialEJA()
885 sage: J.characteristic_polynomial_of()
892 # The list of coefficient polynomials a_0, a_1, a_2, ..., a_(r-1).
893 a
= self
._charpoly
_coefficients
()
895 # We go to a bit of trouble here to reorder the
896 # indeterminates, so that it's easier to evaluate the
897 # characteristic polynomial at x's coordinates and get back
898 # something in terms of t, which is what we want.
899 S
= PolynomialRing(self
.base_ring(),'t')
903 S
= PolynomialRing(S
, R
.variable_names())
906 return (t
**r
+ sum( a
[k
]*(t
**k
) for k
in range(r
) ))
908 def coordinate_polynomial_ring(self
):
910 The multivariate polynomial ring in which this algebra's
911 :meth:`characteristic_polynomial_of` lives.
915 sage: from mjo.eja.eja_algebra import (HadamardEJA,
916 ....: RealSymmetricEJA)
920 sage: J = HadamardEJA(2)
921 sage: J.coordinate_polynomial_ring()
922 Multivariate Polynomial Ring in X1, X2...
923 sage: J = RealSymmetricEJA(3,field=QQ,orthonormalize=False)
924 sage: J.coordinate_polynomial_ring()
925 Multivariate Polynomial Ring in X1, X2, X3, X4, X5, X6...
928 var_names
= tuple( "X%d" % z
for z
in range(1, self
.dimension()+1) )
929 return PolynomialRing(self
.base_ring(), var_names
)
931 def inner_product(self
, x
, y
):
933 The inner product associated with this Euclidean Jordan algebra.
935 Defaults to the trace inner product, but can be overridden by
936 subclasses if they are sure that the necessary properties are
941 sage: from mjo.eja.eja_algebra import (random_eja,
943 ....: BilinearFormEJA)
947 Our inner product is "associative," which means the following for
948 a symmetric bilinear form::
950 sage: set_random_seed()
951 sage: J = random_eja()
952 sage: x,y,z = J.random_elements(3)
953 sage: (x*y).inner_product(z) == y.inner_product(x*z)
958 Ensure that this is the usual inner product for the algebras
961 sage: set_random_seed()
962 sage: J = HadamardEJA.random_instance()
963 sage: x,y = J.random_elements(2)
964 sage: actual = x.inner_product(y)
965 sage: expected = x.to_vector().inner_product(y.to_vector())
966 sage: actual == expected
969 Ensure that this is one-half of the trace inner-product in a
970 BilinearFormEJA that isn't just the reals (when ``n`` isn't
971 one). This is in Faraut and Koranyi, and also my "On the
974 sage: set_random_seed()
975 sage: J = BilinearFormEJA.random_instance()
976 sage: n = J.dimension()
977 sage: x = J.random_element()
978 sage: y = J.random_element()
979 sage: (n == 1) or (x.inner_product(y) == (x*y).trace()/2)
982 B
= self
._inner
_product
_matrix
983 return (B
*x
.to_vector()).inner_product(y
.to_vector())
986 def is_trivial(self
):
988 Return whether or not this algebra is trivial.
990 A trivial algebra contains only the zero element.
994 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
999 sage: J = ComplexHermitianEJA(3)
1000 sage: J.is_trivial()
1005 sage: J = TrivialEJA()
1006 sage: J.is_trivial()
1010 return self
.dimension() == 0
1013 def multiplication_table(self
):
1015 Return a visual representation of this algebra's multiplication
1016 table (on basis elements).
1020 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1024 sage: J = JordanSpinEJA(4)
1025 sage: J.multiplication_table()
1026 +----++----+----+----+----+
1027 | * || b0 | b1 | b2 | b3 |
1028 +====++====+====+====+====+
1029 | b0 || b0 | b1 | b2 | b3 |
1030 +----++----+----+----+----+
1031 | b1 || b1 | b0 | 0 | 0 |
1032 +----++----+----+----+----+
1033 | b2 || b2 | 0 | b0 | 0 |
1034 +----++----+----+----+----+
1035 | b3 || b3 | 0 | 0 | b0 |
1036 +----++----+----+----+----+
1039 n
= self
.dimension()
1040 # Prepend the header row.
1041 M
= [["*"] + list(self
.gens())]
1043 # And to each subsequent row, prepend an entry that belongs to
1044 # the left-side "header column."
1045 M
+= [ [self
.monomial(i
)] + [ self
.monomial(i
)*self
.monomial(j
)
1049 return table(M
, header_row
=True, header_column
=True, frame
=True)
1052 def matrix_basis(self
):
1054 Return an (often more natural) representation of this algebras
1055 basis as an ordered tuple of matrices.
1057 Every finite-dimensional Euclidean Jordan Algebra is a, up to
1058 Jordan isomorphism, a direct sum of five simple
1059 algebras---four of which comprise Hermitian matrices. And the
1060 last type of algebra can of course be thought of as `n`-by-`1`
1061 column matrices (ambiguusly called column vectors) to avoid
1062 special cases. As a result, matrices (and column vectors) are
1063 a natural representation format for Euclidean Jordan algebra
1066 But, when we construct an algebra from a basis of matrices,
1067 those matrix representations are lost in favor of coordinate
1068 vectors *with respect to* that basis. We could eventually
1069 convert back if we tried hard enough, but having the original
1070 representations handy is valuable enough that we simply store
1071 them and return them from this method.
1073 Why implement this for non-matrix algebras? Avoiding special
1074 cases for the :class:`BilinearFormEJA` pays with simplicity in
1075 its own right. But mainly, we would like to be able to assume
1076 that elements of a :class:`CartesianProductEJA` can be displayed
1077 nicely, without having to have special classes for direct sums
1078 one of whose components was a matrix algebra.
1082 sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
1083 ....: RealSymmetricEJA)
1087 sage: J = RealSymmetricEJA(2)
1089 Finite family {0: b0, 1: b1, 2: b2}
1090 sage: J.matrix_basis()
1092 [1 0] [ 0 0.7071067811865475?] [0 0]
1093 [0 0], [0.7071067811865475? 0], [0 1]
1098 sage: J = JordanSpinEJA(2)
1100 Finite family {0: b0, 1: b1}
1101 sage: J.matrix_basis()
1107 return self
._matrix
_basis
1110 def matrix_space(self
):
1112 Return the matrix space in which this algebra's elements live, if
1113 we think of them as matrices (including column vectors of the
1116 "By default" this will be an `n`-by-`1` column-matrix space,
1117 except when the algebra is trivial. There it's `n`-by-`n`
1118 (where `n` is zero), to ensure that two elements of the matrix
1119 space (empty matrices) can be multiplied. For algebras of
1120 matrices, this returns the space in which their
1121 real embeddings live.
1125 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
1126 ....: JordanSpinEJA,
1127 ....: QuaternionHermitianEJA,
1132 By default, the matrix representation is just a column-matrix
1133 equivalent to the vector representation::
1135 sage: J = JordanSpinEJA(3)
1136 sage: J.matrix_space()
1137 Full MatrixSpace of 3 by 1 dense matrices over Algebraic
1140 The matrix representation in the trivial algebra is
1141 zero-by-zero instead of the usual `n`-by-one::
1143 sage: J = TrivialEJA()
1144 sage: J.matrix_space()
1145 Full MatrixSpace of 0 by 0 dense matrices over Algebraic
1148 The matrix space for complex/quaternion Hermitian matrix EJA
1149 is the space in which their real-embeddings live, not the
1150 original complex/quaternion matrix space::
1152 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
1153 sage: J.matrix_space()
1154 Module of 2 by 2 matrices with entries in Algebraic Field over
1155 the scalar ring Rational Field
1156 sage: J = QuaternionHermitianEJA(1,field=QQ,orthonormalize=False)
1157 sage: J.matrix_space()
1158 Module of 1 by 1 matrices with entries in Quaternion
1159 Algebra (-1, -1) with base ring Rational Field over
1160 the scalar ring Rational Field
1163 return self
._matrix
_space
1169 Return the unit element of this algebra.
1173 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1178 We can compute unit element in the Hadamard EJA::
1180 sage: J = HadamardEJA(5)
1182 b0 + b1 + b2 + b3 + b4
1184 The unit element in the Hadamard EJA is inherited in the
1185 subalgebras generated by its elements::
1187 sage: J = HadamardEJA(5)
1189 b0 + b1 + b2 + b3 + b4
1190 sage: x = sum(J.gens())
1191 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1194 sage: A.one().superalgebra_element()
1195 b0 + b1 + b2 + b3 + b4
1199 The identity element acts like the identity, regardless of
1200 whether or not we orthonormalize::
1202 sage: set_random_seed()
1203 sage: J = random_eja()
1204 sage: x = J.random_element()
1205 sage: J.one()*x == x and x*J.one() == x
1207 sage: A = x.subalgebra_generated_by()
1208 sage: y = A.random_element()
1209 sage: A.one()*y == y and y*A.one() == y
1214 sage: set_random_seed()
1215 sage: J = random_eja(field=QQ, orthonormalize=False)
1216 sage: x = J.random_element()
1217 sage: J.one()*x == x and x*J.one() == x
1219 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1220 sage: y = A.random_element()
1221 sage: A.one()*y == y and y*A.one() == y
1224 The matrix of the unit element's operator is the identity,
1225 regardless of the base field and whether or not we
1228 sage: set_random_seed()
1229 sage: J = random_eja()
1230 sage: actual = J.one().operator().matrix()
1231 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1232 sage: actual == expected
1234 sage: x = J.random_element()
1235 sage: A = x.subalgebra_generated_by()
1236 sage: actual = A.one().operator().matrix()
1237 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1238 sage: actual == expected
1243 sage: set_random_seed()
1244 sage: J = random_eja(field=QQ, orthonormalize=False)
1245 sage: actual = J.one().operator().matrix()
1246 sage: expected = matrix.identity(J.base_ring(), J.dimension())
1247 sage: actual == expected
1249 sage: x = J.random_element()
1250 sage: A = x.subalgebra_generated_by(orthonormalize=False)
1251 sage: actual = A.one().operator().matrix()
1252 sage: expected = matrix.identity(A.base_ring(), A.dimension())
1253 sage: actual == expected
1256 Ensure that the cached unit element (often precomputed by
1257 hand) agrees with the computed one::
1259 sage: set_random_seed()
1260 sage: J = random_eja()
1261 sage: cached = J.one()
1262 sage: J.one.clear_cache()
1263 sage: J.one() == cached
1268 sage: set_random_seed()
1269 sage: J = random_eja(field=QQ, orthonormalize=False)
1270 sage: cached = J.one()
1271 sage: J.one.clear_cache()
1272 sage: J.one() == cached
1276 # We can brute-force compute the matrices of the operators
1277 # that correspond to the basis elements of this algebra.
1278 # If some linear combination of those basis elements is the
1279 # algebra identity, then the same linear combination of
1280 # their matrices has to be the identity matrix.
1282 # Of course, matrices aren't vectors in sage, so we have to
1283 # appeal to the "long vectors" isometry.
1285 V
= VectorSpace(self
.base_ring(), self
.dimension()**2)
1286 oper_vecs
= [ V(g
.operator().matrix().list()) for g
in self
.gens() ]
1288 # Now we use basic linear algebra to find the coefficients,
1289 # of the matrices-as-vectors-linear-combination, which should
1290 # work for the original algebra basis too.
1291 A
= matrix(self
.base_ring(), oper_vecs
)
1293 # We used the isometry on the left-hand side already, but we
1294 # still need to do it for the right-hand side. Recall that we
1295 # wanted something that summed to the identity matrix.
1296 b
= V( matrix
.identity(self
.base_ring(), self
.dimension()).list() )
1298 # Now if there's an identity element in the algebra, this
1299 # should work. We solve on the left to avoid having to
1300 # transpose the matrix "A".
1301 return self
.from_vector(A
.solve_left(b
))
1304 def peirce_decomposition(self
, c
):
1306 The Peirce decomposition of this algebra relative to the
1309 In the future, this can be extended to a complete system of
1310 orthogonal idempotents.
1314 - ``c`` -- an idempotent of this algebra.
1318 A triple (J0, J5, J1) containing two subalgebras and one subspace
1321 - ``J0`` -- the algebra on the eigenspace of ``c.operator()``
1322 corresponding to the eigenvalue zero.
1324 - ``J5`` -- the eigenspace (NOT a subalgebra) of ``c.operator()``
1325 corresponding to the eigenvalue one-half.
1327 - ``J1`` -- the algebra on the eigenspace of ``c.operator()``
1328 corresponding to the eigenvalue one.
1330 These are the only possible eigenspaces for that operator, and this
1331 algebra is a direct sum of them. The spaces ``J0`` and ``J1`` are
1332 orthogonal, and are subalgebras of this algebra with the appropriate
1337 sage: from mjo.eja.eja_algebra import random_eja, RealSymmetricEJA
1341 The canonical example comes from the symmetric matrices, which
1342 decompose into diagonal and off-diagonal parts::
1344 sage: J = RealSymmetricEJA(3)
1345 sage: C = matrix(QQ, [ [1,0,0],
1349 sage: J0,J5,J1 = J.peirce_decomposition(c)
1351 Euclidean Jordan algebra of dimension 1...
1353 Vector space of degree 6 and dimension 2...
1355 Euclidean Jordan algebra of dimension 3...
1356 sage: J0.one().to_matrix()
1360 sage: orig_df = AA.options.display_format
1361 sage: AA.options.display_format = 'radical'
1362 sage: J.from_vector(J5.basis()[0]).to_matrix()
1366 sage: J.from_vector(J5.basis()[1]).to_matrix()
1370 sage: AA.options.display_format = orig_df
1371 sage: J1.one().to_matrix()
1378 Every algebra decomposes trivially with respect to its identity
1381 sage: set_random_seed()
1382 sage: J = random_eja()
1383 sage: J0,J5,J1 = J.peirce_decomposition(J.one())
1384 sage: J0.dimension() == 0 and J5.dimension() == 0
1386 sage: J1.superalgebra() == J and J1.dimension() == J.dimension()
1389 The decomposition is into eigenspaces, and its components are
1390 therefore necessarily orthogonal. Moreover, the identity
1391 elements in the two subalgebras are the projections onto their
1392 respective subspaces of the superalgebra's identity element::
1394 sage: set_random_seed()
1395 sage: J = random_eja()
1396 sage: x = J.random_element()
1397 sage: if not J.is_trivial():
1398 ....: while x.is_nilpotent():
1399 ....: x = J.random_element()
1400 sage: c = x.subalgebra_idempotent()
1401 sage: J0,J5,J1 = J.peirce_decomposition(c)
1403 sage: for (w,y,z) in zip(J0.basis(), J5.basis(), J1.basis()):
1404 ....: w = w.superalgebra_element()
1405 ....: y = J.from_vector(y)
1406 ....: z = z.superalgebra_element()
1407 ....: ipsum += w.inner_product(y).abs()
1408 ....: ipsum += w.inner_product(z).abs()
1409 ....: ipsum += y.inner_product(z).abs()
1412 sage: J1(c) == J1.one()
1414 sage: J0(J.one() - c) == J0.one()
1418 if not c
.is_idempotent():
1419 raise ValueError("element is not idempotent: %s" % c
)
1421 # Default these to what they should be if they turn out to be
1422 # trivial, because eigenspaces_left() won't return eigenvalues
1423 # corresponding to trivial spaces (e.g. it returns only the
1424 # eigenspace corresponding to lambda=1 if you take the
1425 # decomposition relative to the identity element).
1426 trivial
= self
.subalgebra((), check_axioms
=False)
1427 J0
= trivial
# eigenvalue zero
1428 J5
= VectorSpace(self
.base_ring(), 0) # eigenvalue one-half
1429 J1
= trivial
# eigenvalue one
1431 for (eigval
, eigspace
) in c
.operator().matrix().right_eigenspaces():
1432 if eigval
== ~
(self
.base_ring()(2)):
1435 gens
= tuple( self
.from_vector(b
) for b
in eigspace
.basis() )
1436 subalg
= self
.subalgebra(gens
, check_axioms
=False)
1442 raise ValueError("unexpected eigenvalue: %s" % eigval
)
1447 def random_element(self
, thorough
=False):
1449 Return a random element of this algebra.
1451 Our algebra superclass method only returns a linear
1452 combination of at most two basis elements. We instead
1453 want the vector space "random element" method that
1454 returns a more diverse selection.
1458 - ``thorough`` -- (boolean; default False) whether or not we
1459 should generate irrational coefficients for the random
1460 element when our base ring is irrational; this slows the
1461 algebra operations to a crawl, but any truly random method
1465 # For a general base ring... maybe we can trust this to do the
1466 # right thing? Unlikely, but.
1467 V
= self
.vector_space()
1468 v
= V
.random_element()
1470 if self
.base_ring() is AA
:
1471 # The "random element" method of the algebraic reals is
1472 # stupid at the moment, and only returns integers between
1473 # -2 and 2, inclusive:
1475 # https://trac.sagemath.org/ticket/30875
1477 # Instead, we implement our own "random vector" method,
1478 # and then coerce that into the algebra. We use the vector
1479 # space degree here instead of the dimension because a
1480 # subalgebra could (for example) be spanned by only two
1481 # vectors, each with five coordinates. We need to
1482 # generate all five coordinates.
1484 v
*= QQbar
.random_element().real()
1486 v
*= QQ
.random_element()
1488 return self
.from_vector(V
.coordinate_vector(v
))
1490 def random_elements(self
, count
, thorough
=False):
1492 Return ``count`` random elements as a tuple.
1496 - ``thorough`` -- (boolean; default False) whether or not we
1497 should generate irrational coefficients for the random
1498 elements when our base ring is irrational; this slows the
1499 algebra operations to a crawl, but any truly random method
1504 sage: from mjo.eja.eja_algebra import JordanSpinEJA
1508 sage: J = JordanSpinEJA(3)
1509 sage: x,y,z = J.random_elements(3)
1510 sage: all( [ x in J, y in J, z in J ])
1512 sage: len( J.random_elements(10) ) == 10
1516 return tuple( self
.random_element(thorough
)
1517 for idx
in range(count
) )
1521 def _charpoly_coefficients(self
):
1523 The `r` polynomial coefficients of the "characteristic polynomial
1528 sage: from mjo.eja.eja_algebra import random_eja
1532 The theory shows that these are all homogeneous polynomials of
1535 sage: set_random_seed()
1536 sage: J = random_eja()
1537 sage: all(p.is_homogeneous() for p in J._charpoly_coefficients())
1541 n
= self
.dimension()
1542 R
= self
.coordinate_polynomial_ring()
1544 F
= R
.fraction_field()
1547 # From a result in my book, these are the entries of the
1548 # basis representation of L_x.
1549 return sum( vars[k
]*self
.monomial(k
).operator().matrix()[i
,j
]
1552 L_x
= matrix(F
, n
, n
, L_x_i_j
)
1555 if self
.rank
.is_in_cache():
1557 # There's no need to pad the system with redundant
1558 # columns if we *know* they'll be redundant.
1561 # Compute an extra power in case the rank is equal to
1562 # the dimension (otherwise, we would stop at x^(r-1)).
1563 x_powers
= [ (L_x
**k
)*self
.one().to_vector()
1564 for k
in range(n
+1) ]
1565 A
= matrix
.column(F
, x_powers
[:n
])
1566 AE
= A
.extended_echelon_form()
1573 # The theory says that only the first "r" coefficients are
1574 # nonzero, and they actually live in the original polynomial
1575 # ring and not the fraction field. We negate them because in
1576 # the actual characteristic polynomial, they get moved to the
1577 # other side where x^r lives. We don't bother to trim A_rref
1578 # down to a square matrix and solve the resulting system,
1579 # because the upper-left r-by-r portion of A_rref is
1580 # guaranteed to be the identity matrix, so e.g.
1582 # A_rref.solve_right(Y)
1584 # would just be returning Y.
1585 return (-E
*b
)[:r
].change_ring(R
)
1590 Return the rank of this EJA.
1592 This is a cached method because we know the rank a priori for
1593 all of the algebras we can construct. Thus we can avoid the
1594 expensive ``_charpoly_coefficients()`` call unless we truly
1595 need to compute the whole characteristic polynomial.
1599 sage: from mjo.eja.eja_algebra import (HadamardEJA,
1600 ....: JordanSpinEJA,
1601 ....: RealSymmetricEJA,
1602 ....: ComplexHermitianEJA,
1603 ....: QuaternionHermitianEJA,
1608 The rank of the Jordan spin algebra is always two::
1610 sage: JordanSpinEJA(2).rank()
1612 sage: JordanSpinEJA(3).rank()
1614 sage: JordanSpinEJA(4).rank()
1617 The rank of the `n`-by-`n` Hermitian real, complex, or
1618 quaternion matrices is `n`::
1620 sage: RealSymmetricEJA(4).rank()
1622 sage: ComplexHermitianEJA(3).rank()
1624 sage: QuaternionHermitianEJA(2).rank()
1629 Ensure that every EJA that we know how to construct has a
1630 positive integer rank, unless the algebra is trivial in
1631 which case its rank will be zero::
1633 sage: set_random_seed()
1634 sage: J = random_eja()
1638 sage: r > 0 or (r == 0 and J.is_trivial())
1641 Ensure that computing the rank actually works, since the ranks
1642 of all simple algebras are known and will be cached by default::
1644 sage: set_random_seed() # long time
1645 sage: J = random_eja() # long time
1646 sage: cached = J.rank() # long time
1647 sage: J.rank.clear_cache() # long time
1648 sage: J.rank() == cached # long time
1652 return len(self
._charpoly
_coefficients
())
1655 def subalgebra(self
, basis
, **kwargs
):
1657 Create a subalgebra of this algebra from the given basis.
1659 from mjo
.eja
.eja_subalgebra
import FiniteDimensionalEJASubalgebra
1660 return FiniteDimensionalEJASubalgebra(self
, basis
, **kwargs
)
1663 def vector_space(self
):
1665 Return the vector space that underlies this algebra.
1669 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
1673 sage: J = RealSymmetricEJA(2)
1674 sage: J.vector_space()
1675 Vector space of dimension 3 over...
1678 return self
.zero().to_vector().parent().ambient_vector_space()
1682 class RationalBasisEJA(FiniteDimensionalEJA
):
1684 Algebras whose supplied basis elements have all rational entries.
1688 sage: from mjo.eja.eja_algebra import BilinearFormEJA
1692 The supplied basis is orthonormalized by default::
1694 sage: B = matrix(QQ, [[1, 0, 0], [0, 25, -32], [0, -32, 41]])
1695 sage: J = BilinearFormEJA(B)
1696 sage: J.matrix_basis()
1713 # Abuse the check_field parameter to check that the entries of
1714 # out basis (in ambient coordinates) are in the field QQ.
1715 # Use _all2list to get the vector coordinates of octonion
1716 # entries and not the octonions themselves (which are not
1718 if not all( all(b_i
in QQ
for b_i
in _all2list(b
))
1720 raise TypeError("basis not rational")
1722 super().__init
__(basis
,
1726 check_field
=check_field
,
1729 self
._rational
_algebra
= None
1731 # There's no point in constructing the extra algebra if this
1732 # one is already rational.
1734 # Note: the same Jordan and inner-products work here,
1735 # because they are necessarily defined with respect to
1736 # ambient coordinates and not any particular basis.
1737 self
._rational
_algebra
= FiniteDimensionalEJA(
1742 matrix_space
=self
.matrix_space(),
1743 associative
=self
.is_associative(),
1744 orthonormalize
=False,
1748 def rational_algebra(self
):
1749 # Using None as a flag here (rather than just assigning "self"
1750 # to self._rational_algebra by default) feels a little bit
1751 # more sane to me in a garbage-collected environment.
1752 if self
._rational
_algebra
is None:
1755 return self
._rational
_algebra
1758 def _charpoly_coefficients(self
):
1762 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
1763 ....: JordanSpinEJA)
1767 The base ring of the resulting polynomial coefficients is what
1768 it should be, and not the rationals (unless the algebra was
1769 already over the rationals)::
1771 sage: J = JordanSpinEJA(3)
1772 sage: J._charpoly_coefficients()
1773 (X1^2 - X2^2 - X3^2, -2*X1)
1774 sage: a0 = J._charpoly_coefficients()[0]
1776 Algebraic Real Field
1777 sage: a0.base_ring()
1778 Algebraic Real Field
1781 if self
.rational_algebra() is self
:
1782 # Bypass the hijinks if they won't benefit us.
1783 return super()._charpoly
_coefficients
()
1785 # Do the computation over the rationals. The answer will be
1786 # the same, because all we've done is a change of basis.
1787 # Then, change back from QQ to our real base ring
1788 a
= ( a_i
.change_ring(self
.base_ring())
1789 for a_i
in self
.rational_algebra()._charpoly
_coefficients
() )
1791 # Otherwise, convert the coordinate variables back to the
1792 # deorthonormalized ones.
1793 R
= self
.coordinate_polynomial_ring()
1794 from sage
.modules
.free_module_element
import vector
1795 X
= vector(R
, R
.gens())
1796 BX
= self
._deortho
_matrix
*X
1798 subs_dict
= { X[i]: BX[i] for i in range(len(X)) }
1799 return tuple( a_i
.subs(subs_dict
) for a_i
in a
)
1801 class ConcreteEJA(FiniteDimensionalEJA
):
1803 A class for the Euclidean Jordan algebras that we know by name.
1805 These are the Jordan algebras whose basis, multiplication table,
1806 rank, and so on are known a priori. More to the point, they are
1807 the Euclidean Jordan algebras for which we are able to conjure up
1808 a "random instance."
1812 sage: from mjo.eja.eja_algebra import ConcreteEJA
1816 Our basis is normalized with respect to the algebra's inner
1817 product, unless we specify otherwise::
1819 sage: set_random_seed()
1820 sage: J = ConcreteEJA.random_instance()
1821 sage: all( b.norm() == 1 for b in J.gens() )
1824 Since our basis is orthonormal with respect to the algebra's inner
1825 product, and since we know that this algebra is an EJA, any
1826 left-multiplication operator's matrix will be symmetric because
1827 natural->EJA basis representation is an isometry and within the
1828 EJA the operator is self-adjoint by the Jordan axiom::
1830 sage: set_random_seed()
1831 sage: J = ConcreteEJA.random_instance()
1832 sage: x = J.random_element()
1833 sage: x.operator().is_self_adjoint()
1838 def _max_random_instance_dimension():
1840 The maximum dimension of any random instance. Ten dimensions seems
1841 to be about the point where everything takes a turn for the
1842 worse. And dimension ten (but not nine) allows the 4-by-4 real
1843 Hermitian matrices, the 2-by-2 quaternion Hermitian matrices,
1844 and the 2-by-2 octonion Hermitian matrices.
1849 def _max_random_instance_size(max_dimension
):
1851 Return an integer "size" that is an upper bound on the size of
1852 this algebra when it is used in a random test case. This size
1853 (which can be passed to the algebra's constructor) is itself
1854 based on the ``max_dimension`` parameter.
1856 This method must be implemented in each subclass.
1858 raise NotImplementedError
1861 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
1863 Return a random instance of this type of algebra whose dimension
1864 is less than or equal to the lesser of ``max_dimension`` and
1865 the value returned by ``_max_random_instance_dimension()``. If
1866 the dimension bound is omitted, then only the
1867 ``_max_random_instance_dimension()`` is used as a bound.
1869 This method should be implemented in each subclass.
1873 sage: from mjo.eja.eja_algebra import ConcreteEJA
1877 Both the class bound and the ``max_dimension`` argument are upper
1878 bounds on the dimension of the algebra returned::
1880 sage: from sage.misc.prandom import choice
1881 sage: eja_class = choice(ConcreteEJA.__subclasses__())
1882 sage: class_max_d = eja_class._max_random_instance_dimension()
1883 sage: J = eja_class.random_instance(max_dimension=20,
1885 ....: orthonormalize=False)
1886 sage: J.dimension() <= class_max_d
1888 sage: J = eja_class.random_instance(max_dimension=2,
1890 ....: orthonormalize=False)
1891 sage: J.dimension() <= 2
1895 from sage
.misc
.prandom
import choice
1896 eja_class
= choice(cls
.__subclasses
__())
1898 # These all bubble up to the RationalBasisEJA superclass
1899 # constructor, so any (kw)args valid there are also valid
1901 return eja_class
.random_instance(max_dimension
, *args
, **kwargs
)
1904 class MatrixEJA(FiniteDimensionalEJA
):
1906 def _denormalized_basis(A
):
1908 Returns a basis for the space of complex Hermitian n-by-n matrices.
1910 Why do we embed these? Basically, because all of numerical linear
1911 algebra assumes that you're working with vectors consisting of `n`
1912 entries from a field and scalars from the same field. There's no way
1913 to tell SageMath that (for example) the vectors contain complex
1914 numbers, while the scalar field is real.
1918 sage: from mjo.hurwitz import (ComplexMatrixAlgebra,
1919 ....: QuaternionMatrixAlgebra,
1920 ....: OctonionMatrixAlgebra)
1921 sage: from mjo.eja.eja_algebra import MatrixEJA
1925 sage: set_random_seed()
1926 sage: n = ZZ.random_element(1,5)
1927 sage: A = MatrixSpace(QQ, n)
1928 sage: B = MatrixEJA._denormalized_basis(A)
1929 sage: all( M.is_hermitian() for M in B)
1934 sage: set_random_seed()
1935 sage: n = ZZ.random_element(1,5)
1936 sage: A = ComplexMatrixAlgebra(n, scalars=QQ)
1937 sage: B = MatrixEJA._denormalized_basis(A)
1938 sage: all( M.is_hermitian() for M in B)
1943 sage: set_random_seed()
1944 sage: n = ZZ.random_element(1,5)
1945 sage: A = QuaternionMatrixAlgebra(n, scalars=QQ)
1946 sage: B = MatrixEJA._denormalized_basis(A)
1947 sage: all( M.is_hermitian() for M in B )
1952 sage: set_random_seed()
1953 sage: n = ZZ.random_element(1,5)
1954 sage: A = OctonionMatrixAlgebra(n, scalars=QQ)
1955 sage: B = MatrixEJA._denormalized_basis(A)
1956 sage: all( M.is_hermitian() for M in B )
1960 # These work for real MatrixSpace, whose monomials only have
1961 # two coordinates (because the last one would always be "1").
1962 es
= A
.base_ring().gens()
1963 gen
= lambda A
,m
: A
.monomial(m
[:2])
1965 if hasattr(A
, 'entry_algebra_gens'):
1966 # We've got a MatrixAlgebra, and its monomials will have
1967 # three coordinates.
1968 es
= A
.entry_algebra_gens()
1969 gen
= lambda A
,m
: A
.monomial(m
)
1972 for i
in range(A
.nrows()):
1973 for j
in range(i
+1):
1975 E_ii
= gen(A
, (i
,j
,es
[0]))
1979 E_ij
= gen(A
, (i
,j
,e
))
1980 E_ij
+= E_ij
.conjugate_transpose()
1983 return tuple( basis
)
1986 def jordan_product(X
,Y
):
1987 return (X
*Y
+ Y
*X
)/2
1990 def trace_inner_product(X
,Y
):
1992 A trace inner-product for matrices that aren't embedded in the
1993 reals. It takes MATRICES as arguments, not EJA elements.
1997 sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
1998 ....: ComplexHermitianEJA,
1999 ....: QuaternionHermitianEJA,
2000 ....: OctonionHermitianEJA)
2004 sage: J = RealSymmetricEJA(2,field=QQ,orthonormalize=False)
2005 sage: I = J.one().to_matrix()
2006 sage: J.trace_inner_product(I, -I)
2011 sage: J = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
2012 sage: I = J.one().to_matrix()
2013 sage: J.trace_inner_product(I, -I)
2018 sage: J = QuaternionHermitianEJA(2,field=QQ,orthonormalize=False)
2019 sage: I = J.one().to_matrix()
2020 sage: J.trace_inner_product(I, -I)
2025 sage: J = OctonionHermitianEJA(2,field=QQ,orthonormalize=False)
2026 sage: I = J.one().to_matrix()
2027 sage: J.trace_inner_product(I, -I)
2032 if hasattr(tr
, 'coefficient'):
2033 # Works for octonions, and has to come first because they
2034 # also have a "real()" method that doesn't return an
2035 # element of the scalar ring.
2036 return tr
.coefficient(0)
2037 elif hasattr(tr
, 'coefficient_tuple'):
2038 # Works for quaternions.
2039 return tr
.coefficient_tuple()[0]
2041 # Works for real and complex numbers.
2045 def __init__(self
, matrix_space
, **kwargs
):
2046 # We know this is a valid EJA, but will double-check
2047 # if the user passes check_axioms=True.
2048 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2050 super().__init
__(self
._denormalized
_basis
(matrix_space
),
2051 self
.jordan_product
,
2052 self
.trace_inner_product
,
2053 field
=matrix_space
.base_ring(),
2054 matrix_space
=matrix_space
,
2057 self
.rank
.set_cache(matrix_space
.nrows())
2058 self
.one
.set_cache( self(matrix_space
.one()) )
2060 class RealSymmetricEJA(MatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2062 The rank-n simple EJA consisting of real symmetric n-by-n
2063 matrices, the usual symmetric Jordan product, and the trace inner
2064 product. It has dimension `(n^2 + n)/2` over the reals.
2068 sage: from mjo.eja.eja_algebra import RealSymmetricEJA
2072 sage: J = RealSymmetricEJA(2)
2073 sage: b0, b1, b2 = J.gens()
2081 In theory, our "field" can be any subfield of the reals::
2083 sage: RealSymmetricEJA(2, field=RDF, check_axioms=True)
2084 Euclidean Jordan algebra of dimension 3 over Real Double Field
2085 sage: RealSymmetricEJA(2, field=RR, check_axioms=True)
2086 Euclidean Jordan algebra of dimension 3 over Real Field with
2087 53 bits of precision
2091 The dimension of this algebra is `(n^2 + n) / 2`::
2093 sage: set_random_seed()
2094 sage: d = RealSymmetricEJA._max_random_instance_dimension()
2095 sage: n = RealSymmetricEJA._max_random_instance_size(d)
2096 sage: J = RealSymmetricEJA(n)
2097 sage: J.dimension() == (n^2 + n)/2
2100 The Jordan multiplication is what we think it is::
2102 sage: set_random_seed()
2103 sage: J = RealSymmetricEJA.random_instance()
2104 sage: x,y = J.random_elements(2)
2105 sage: actual = (x*y).to_matrix()
2106 sage: X = x.to_matrix()
2107 sage: Y = y.to_matrix()
2108 sage: expected = (X*Y + Y*X)/2
2109 sage: actual == expected
2111 sage: J(expected) == x*y
2114 We can change the generator prefix::
2116 sage: RealSymmetricEJA(3, prefix='q').gens()
2117 (q0, q1, q2, q3, q4, q5)
2119 We can construct the (trivial) algebra of rank zero::
2121 sage: RealSymmetricEJA(0)
2122 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2126 def _max_random_instance_size(max_dimension
):
2127 # Obtained by solving d = (n^2 + n)/2.
2128 # The ZZ-int-ZZ thing is just "floor."
2129 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/2 - 1/2))
2132 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2134 Return a random instance of this type of algebra.
2136 class_max_d
= cls
._max
_random
_instance
_dimension
()
2137 if (max_dimension
is None or max_dimension
> class_max_d
):
2138 max_dimension
= class_max_d
2139 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2140 n
= ZZ
.random_element(max_size
+ 1)
2141 return cls(n
, **kwargs
)
2143 def __init__(self
, n
, field
=AA
, **kwargs
):
2144 A
= MatrixSpace(field
, n
)
2145 super().__init
__(A
, **kwargs
)
2147 from mjo
.eja
.eja_cache
import real_symmetric_eja_coeffs
2148 a
= real_symmetric_eja_coeffs(self
)
2150 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2154 class ComplexHermitianEJA(MatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2156 The rank-n simple EJA consisting of complex Hermitian n-by-n
2157 matrices over the real numbers, the usual symmetric Jordan product,
2158 and the real-part-of-trace inner product. It has dimension `n^2` over
2163 sage: from mjo.eja.eja_algebra import ComplexHermitianEJA
2167 In theory, our "field" can be any subfield of the reals, but we
2168 can't use inexact real fields at the moment because SageMath
2169 doesn't know how to convert their elements into complex numbers,
2170 or even into algebraic reals::
2173 Traceback (most recent call last):
2175 TypeError: Illegal initializer for algebraic number
2177 Traceback (most recent call last):
2179 TypeError: Illegal initializer for algebraic number
2181 This causes the following error when we try to scale a matrix of
2182 complex numbers by an inexact real number::
2184 sage: ComplexHermitianEJA(2,field=RR)
2185 Traceback (most recent call last):
2187 TypeError: Unable to coerce entries (=(1.00000000000000,
2188 -0.000000000000000)) to coefficients in Algebraic Real Field
2192 The dimension of this algebra is `n^2`::
2194 sage: set_random_seed()
2195 sage: d = ComplexHermitianEJA._max_random_instance_dimension()
2196 sage: n = ComplexHermitianEJA._max_random_instance_size(d)
2197 sage: J = ComplexHermitianEJA(n)
2198 sage: J.dimension() == n^2
2201 The Jordan multiplication is what we think it is::
2203 sage: set_random_seed()
2204 sage: J = ComplexHermitianEJA.random_instance()
2205 sage: x,y = J.random_elements(2)
2206 sage: actual = (x*y).to_matrix()
2207 sage: X = x.to_matrix()
2208 sage: Y = y.to_matrix()
2209 sage: expected = (X*Y + Y*X)/2
2210 sage: actual == expected
2212 sage: J(expected) == x*y
2215 We can change the generator prefix::
2217 sage: ComplexHermitianEJA(2, prefix='z').gens()
2220 We can construct the (trivial) algebra of rank zero::
2222 sage: ComplexHermitianEJA(0)
2223 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2226 def __init__(self
, n
, field
=AA
, **kwargs
):
2227 from mjo
.hurwitz
import ComplexMatrixAlgebra
2228 A
= ComplexMatrixAlgebra(n
, scalars
=field
)
2229 super().__init
__(A
, **kwargs
)
2231 from mjo
.eja
.eja_cache
import complex_hermitian_eja_coeffs
2232 a
= complex_hermitian_eja_coeffs(self
)
2234 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2237 def _max_random_instance_size(max_dimension
):
2238 # Obtained by solving d = n^2.
2239 # The ZZ-int-ZZ thing is just "floor."
2240 return ZZ(int(ZZ(max_dimension
).sqrt()))
2243 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2245 Return a random instance of this type of algebra.
2247 class_max_d
= cls
._max
_random
_instance
_dimension
()
2248 if (max_dimension
is None or max_dimension
> class_max_d
):
2249 max_dimension
= class_max_d
2250 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2251 n
= ZZ
.random_element(max_size
+ 1)
2252 return cls(n
, **kwargs
)
2255 class QuaternionHermitianEJA(MatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2257 The rank-n simple EJA consisting of self-adjoint n-by-n quaternion
2258 matrices, the usual symmetric Jordan product, and the
2259 real-part-of-trace inner product. It has dimension `2n^2 - n` over
2264 sage: from mjo.eja.eja_algebra import QuaternionHermitianEJA
2268 In theory, our "field" can be any subfield of the reals::
2270 sage: QuaternionHermitianEJA(2, field=RDF, check_axioms=True)
2271 Euclidean Jordan algebra of dimension 6 over Real Double Field
2272 sage: QuaternionHermitianEJA(2, field=RR, check_axioms=True)
2273 Euclidean Jordan algebra of dimension 6 over Real Field with
2274 53 bits of precision
2278 The dimension of this algebra is `2*n^2 - n`::
2280 sage: set_random_seed()
2281 sage: d = QuaternionHermitianEJA._max_random_instance_dimension()
2282 sage: n = QuaternionHermitianEJA._max_random_instance_size(d)
2283 sage: J = QuaternionHermitianEJA(n)
2284 sage: J.dimension() == 2*(n^2) - n
2287 The Jordan multiplication is what we think it is::
2289 sage: set_random_seed()
2290 sage: J = QuaternionHermitianEJA.random_instance()
2291 sage: x,y = J.random_elements(2)
2292 sage: actual = (x*y).to_matrix()
2293 sage: X = x.to_matrix()
2294 sage: Y = y.to_matrix()
2295 sage: expected = (X*Y + Y*X)/2
2296 sage: actual == expected
2298 sage: J(expected) == x*y
2301 We can change the generator prefix::
2303 sage: QuaternionHermitianEJA(2, prefix='a').gens()
2304 (a0, a1, a2, a3, a4, a5)
2306 We can construct the (trivial) algebra of rank zero::
2308 sage: QuaternionHermitianEJA(0)
2309 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2312 def __init__(self
, n
, field
=AA
, **kwargs
):
2313 from mjo
.hurwitz
import QuaternionMatrixAlgebra
2314 A
= QuaternionMatrixAlgebra(n
, scalars
=field
)
2315 super().__init
__(A
, **kwargs
)
2317 from mjo
.eja
.eja_cache
import quaternion_hermitian_eja_coeffs
2318 a
= quaternion_hermitian_eja_coeffs(self
)
2320 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2325 def _max_random_instance_size(max_dimension
):
2327 The maximum rank of a random QuaternionHermitianEJA.
2329 # Obtained by solving d = 2n^2 - n.
2330 # The ZZ-int-ZZ thing is just "floor."
2331 return ZZ(int(ZZ(8*max_dimension
+ 1).sqrt()/4 + 1/4))
2334 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2336 Return a random instance of this type of algebra.
2338 class_max_d
= cls
._max
_random
_instance
_dimension
()
2339 if (max_dimension
is None or max_dimension
> class_max_d
):
2340 max_dimension
= class_max_d
2341 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2342 n
= ZZ
.random_element(max_size
+ 1)
2343 return cls(n
, **kwargs
)
2345 class OctonionHermitianEJA(MatrixEJA
, RationalBasisEJA
, ConcreteEJA
):
2349 sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA,
2350 ....: OctonionHermitianEJA)
2351 sage: from mjo.hurwitz import Octonions, OctonionMatrixAlgebra
2355 The 3-by-3 algebra satisfies the axioms of an EJA::
2357 sage: OctonionHermitianEJA(3, # long time
2358 ....: field=QQ, # long time
2359 ....: orthonormalize=False, # long time
2360 ....: check_axioms=True) # long time
2361 Euclidean Jordan algebra of dimension 27 over Rational Field
2363 After a change-of-basis, the 2-by-2 algebra has the same
2364 multiplication table as the ten-dimensional Jordan spin algebra::
2366 sage: A = OctonionMatrixAlgebra(2,Octonions(QQ),QQ)
2367 sage: b = OctonionHermitianEJA._denormalized_basis(A)
2368 sage: basis = (b[0] + b[9],) + b[1:9] + (b[0] - b[9],)
2369 sage: jp = OctonionHermitianEJA.jordan_product
2370 sage: ip = OctonionHermitianEJA.trace_inner_product
2371 sage: J = FiniteDimensionalEJA(basis,
2375 ....: orthonormalize=False)
2376 sage: J.multiplication_table()
2377 +----++----+----+----+----+----+----+----+----+----+----+
2378 | * || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 |
2379 +====++====+====+====+====+====+====+====+====+====+====+
2380 | b0 || b0 | b1 | b2 | b3 | b4 | b5 | b6 | b7 | b8 | b9 |
2381 +----++----+----+----+----+----+----+----+----+----+----+
2382 | b1 || b1 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
2383 +----++----+----+----+----+----+----+----+----+----+----+
2384 | b2 || b2 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
2385 +----++----+----+----+----+----+----+----+----+----+----+
2386 | b3 || b3 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 | 0 |
2387 +----++----+----+----+----+----+----+----+----+----+----+
2388 | b4 || b4 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 | 0 |
2389 +----++----+----+----+----+----+----+----+----+----+----+
2390 | b5 || b5 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 | 0 |
2391 +----++----+----+----+----+----+----+----+----+----+----+
2392 | b6 || b6 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 | 0 |
2393 +----++----+----+----+----+----+----+----+----+----+----+
2394 | b7 || b7 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 | 0 |
2395 +----++----+----+----+----+----+----+----+----+----+----+
2396 | b8 || b8 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 | 0 |
2397 +----++----+----+----+----+----+----+----+----+----+----+
2398 | b9 || b9 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | b0 |
2399 +----++----+----+----+----+----+----+----+----+----+----+
2403 We can actually construct the 27-dimensional Albert algebra,
2404 and we get the right unit element if we recompute it::
2406 sage: J = OctonionHermitianEJA(3, # long time
2407 ....: field=QQ, # long time
2408 ....: orthonormalize=False) # long time
2409 sage: J.one.clear_cache() # long time
2410 sage: J.one() # long time
2412 sage: J.one().to_matrix() # long time
2421 The 2-by-2 algebra is isomorphic to the ten-dimensional Jordan
2422 spin algebra, but just to be sure, we recompute its rank::
2424 sage: J = OctonionHermitianEJA(2, # long time
2425 ....: field=QQ, # long time
2426 ....: orthonormalize=False) # long time
2427 sage: J.rank.clear_cache() # long time
2428 sage: J.rank() # long time
2433 def _max_random_instance_size(max_dimension
):
2435 The maximum rank of a random QuaternionHermitianEJA.
2437 # There's certainly a formula for this, but with only four
2438 # cases to worry about, I'm not that motivated to derive it.
2439 if max_dimension
>= 27:
2441 elif max_dimension
>= 10:
2443 elif max_dimension
>= 1:
2449 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2451 Return a random instance of this type of algebra.
2453 class_max_d
= cls
._max
_random
_instance
_dimension
()
2454 if (max_dimension
is None or max_dimension
> class_max_d
):
2455 max_dimension
= class_max_d
2456 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2457 n
= ZZ
.random_element(max_size
+ 1)
2458 return cls(n
, **kwargs
)
2460 def __init__(self
, n
, field
=AA
, **kwargs
):
2462 # Otherwise we don't get an EJA.
2463 raise ValueError("n cannot exceed 3")
2465 # We know this is a valid EJA, but will double-check
2466 # if the user passes check_axioms=True.
2467 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2469 from mjo
.hurwitz
import OctonionMatrixAlgebra
2470 A
= OctonionMatrixAlgebra(n
, scalars
=field
)
2471 super().__init
__(A
, **kwargs
)
2473 from mjo
.eja
.eja_cache
import octonion_hermitian_eja_coeffs
2474 a
= octonion_hermitian_eja_coeffs(self
)
2476 self
.rational_algebra()._charpoly
_coefficients
.set_cache(a
)
2479 class AlbertEJA(OctonionHermitianEJA
):
2481 The Albert algebra is the algebra of three-by-three Hermitian
2482 matrices whose entries are octonions.
2486 sage: from mjo.eja.eja_algebra import AlbertEJA
2490 sage: AlbertEJA(field=QQ, orthonormalize=False)
2491 Euclidean Jordan algebra of dimension 27 over Rational Field
2492 sage: AlbertEJA() # long time
2493 Euclidean Jordan algebra of dimension 27 over Algebraic Real Field
2496 def __init__(self
, *args
, **kwargs
):
2497 super().__init
__(3, *args
, **kwargs
)
2500 class HadamardEJA(RationalBasisEJA
, ConcreteEJA
):
2502 Return the Euclidean Jordan algebra on `R^n` with the Hadamard
2503 (pointwise real-number multiplication) Jordan product and the
2504 usual inner-product.
2506 This is nothing more than the Cartesian product of ``n`` copies of
2507 the one-dimensional Jordan spin algebra, and is the most common
2508 example of a non-simple Euclidean Jordan algebra.
2512 sage: from mjo.eja.eja_algebra import HadamardEJA
2516 This multiplication table can be verified by hand::
2518 sage: J = HadamardEJA(3)
2519 sage: b0,b1,b2 = J.gens()
2535 We can change the generator prefix::
2537 sage: HadamardEJA(3, prefix='r').gens()
2540 def __init__(self
, n
, field
=AA
, **kwargs
):
2541 MS
= MatrixSpace(field
, n
, 1)
2544 jordan_product
= lambda x
,y
: x
2545 inner_product
= lambda x
,y
: x
2547 def jordan_product(x
,y
):
2548 return MS( xi
*yi
for (xi
,yi
) in zip(x
,y
) )
2550 def inner_product(x
,y
):
2553 # New defaults for keyword arguments. Don't orthonormalize
2554 # because our basis is already orthonormal with respect to our
2555 # inner-product. Don't check the axioms, because we know this
2556 # is a valid EJA... but do double-check if the user passes
2557 # check_axioms=True. Note: we DON'T override the "check_field"
2558 # default here, because the user can pass in a field!
2559 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2560 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2562 column_basis
= tuple( MS(b
) for b
in FreeModule(field
, n
).basis() )
2563 super().__init
__(column_basis
,
2570 self
.rank
.set_cache(n
)
2572 self
.one
.set_cache( self
.sum(self
.gens()) )
2575 def _max_random_instance_dimension():
2577 There's no reason to go higher than five here. That's
2578 enough to get the point across.
2583 def _max_random_instance_size(max_dimension
):
2585 The maximum size (=dimension) of a random HadamardEJA.
2587 return max_dimension
2590 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2592 Return a random instance of this type of algebra.
2594 class_max_d
= cls
._max
_random
_instance
_dimension
()
2595 if (max_dimension
is None or max_dimension
> class_max_d
):
2596 max_dimension
= class_max_d
2597 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2598 n
= ZZ
.random_element(max_size
+ 1)
2599 return cls(n
, **kwargs
)
2602 class BilinearFormEJA(RationalBasisEJA
, ConcreteEJA
):
2604 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2605 with the half-trace inner product and jordan product ``x*y =
2606 (<Bx,y>,y_bar>, x0*y_bar + y0*x_bar)`` where `B = 1 \times B22` is
2607 a symmetric positive-definite "bilinear form" matrix. Its
2608 dimension is the size of `B`, and it has rank two in dimensions
2609 larger than two. It reduces to the ``JordanSpinEJA`` when `B` is
2610 the identity matrix of order ``n``.
2612 We insist that the one-by-one upper-left identity block of `B` be
2613 passed in as well so that we can be passed a matrix of size zero
2614 to construct a trivial algebra.
2618 sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
2619 ....: JordanSpinEJA)
2623 When no bilinear form is specified, the identity matrix is used,
2624 and the resulting algebra is the Jordan spin algebra::
2626 sage: B = matrix.identity(AA,3)
2627 sage: J0 = BilinearFormEJA(B)
2628 sage: J1 = JordanSpinEJA(3)
2629 sage: J0.multiplication_table() == J0.multiplication_table()
2632 An error is raised if the matrix `B` does not correspond to a
2633 positive-definite bilinear form::
2635 sage: B = matrix.random(QQ,2,3)
2636 sage: J = BilinearFormEJA(B)
2637 Traceback (most recent call last):
2639 ValueError: bilinear form is not positive-definite
2640 sage: B = matrix.zero(QQ,3)
2641 sage: J = BilinearFormEJA(B)
2642 Traceback (most recent call last):
2644 ValueError: bilinear form is not positive-definite
2648 We can create a zero-dimensional algebra::
2650 sage: B = matrix.identity(AA,0)
2651 sage: J = BilinearFormEJA(B)
2655 We can check the multiplication condition given in the Jordan, von
2656 Neumann, and Wigner paper (and also discussed on my "On the
2657 symmetry..." paper). Note that this relies heavily on the standard
2658 choice of basis, as does anything utilizing the bilinear form
2659 matrix. We opt not to orthonormalize the basis, because if we
2660 did, we would have to normalize the `s_{i}` in a similar manner::
2662 sage: set_random_seed()
2663 sage: n = ZZ.random_element(5)
2664 sage: M = matrix.random(QQ, max(0,n-1), algorithm='unimodular')
2665 sage: B11 = matrix.identity(QQ,1)
2666 sage: B22 = M.transpose()*M
2667 sage: B = block_matrix(2,2,[ [B11,0 ],
2669 sage: J = BilinearFormEJA(B, orthonormalize=False)
2670 sage: eis = VectorSpace(M.base_ring(), M.ncols()).basis()
2671 sage: V = J.vector_space()
2672 sage: sis = [ J( V([0] + (M.inverse()*ei).list()).column() )
2673 ....: for ei in eis ]
2674 sage: actual = [ sis[i]*sis[j]
2675 ....: for i in range(n-1)
2676 ....: for j in range(n-1) ]
2677 sage: expected = [ J.one() if i == j else J.zero()
2678 ....: for i in range(n-1)
2679 ....: for j in range(n-1) ]
2680 sage: actual == expected
2684 def __init__(self
, B
, field
=AA
, **kwargs
):
2685 # The matrix "B" is supplied by the user in most cases,
2686 # so it makes sense to check whether or not its positive-
2687 # definite unless we are specifically asked not to...
2688 if ("check_axioms" not in kwargs
) or kwargs
["check_axioms"]:
2689 if not B
.is_positive_definite():
2690 raise ValueError("bilinear form is not positive-definite")
2692 # However, all of the other data for this EJA is computed
2693 # by us in manner that guarantees the axioms are
2694 # satisfied. So, again, unless we are specifically asked to
2695 # verify things, we'll skip the rest of the checks.
2696 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2699 MS
= MatrixSpace(field
, n
, 1)
2701 def inner_product(x
,y
):
2702 return (y
.T
*B
*x
)[0,0]
2704 def jordan_product(x
,y
):
2709 z0
= inner_product(y
,x
)
2710 zbar
= y0
*xbar
+ x0
*ybar
2711 return MS([z0
] + zbar
.list())
2713 column_basis
= tuple( MS(b
) for b
in FreeModule(field
, n
).basis() )
2715 # TODO: I haven't actually checked this, but it seems legit.
2720 super().__init
__(column_basis
,
2725 associative
=associative
,
2728 # The rank of this algebra is two, unless we're in a
2729 # one-dimensional ambient space (because the rank is bounded
2730 # by the ambient dimension).
2731 self
.rank
.set_cache(min(n
,2))
2733 self
.one
.set_cache( self
.zero() )
2735 self
.one
.set_cache( self
.monomial(0) )
2738 def _max_random_instance_dimension():
2740 There's no reason to go higher than five here. That's
2741 enough to get the point across.
2746 def _max_random_instance_size(max_dimension
):
2748 The maximum size (=dimension) of a random BilinearFormEJA.
2750 return max_dimension
2753 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2755 Return a random instance of this algebra.
2757 class_max_d
= cls
._max
_random
_instance
_dimension
()
2758 if (max_dimension
is None or max_dimension
> class_max_d
):
2759 max_dimension
= class_max_d
2760 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2761 n
= ZZ
.random_element(max_size
+ 1)
2764 B
= matrix
.identity(ZZ
, n
)
2765 return cls(B
, **kwargs
)
2767 B11
= matrix
.identity(ZZ
, 1)
2768 M
= matrix
.random(ZZ
, n
-1)
2769 I
= matrix
.identity(ZZ
, n
-1)
2771 while alpha
.is_zero():
2772 alpha
= ZZ
.random_element().abs()
2774 B22
= M
.transpose()*M
+ alpha
*I
2776 from sage
.matrix
.special
import block_matrix
2777 B
= block_matrix(2,2, [ [B11
, ZZ(0) ],
2780 return cls(B
, **kwargs
)
2783 class JordanSpinEJA(BilinearFormEJA
):
2785 The rank-2 simple EJA consisting of real vectors ``x=(x0, x_bar)``
2786 with the usual inner product and jordan product ``x*y =
2787 (<x,y>, x0*y_bar + y0*x_bar)``. It has dimension `n` over
2792 sage: from mjo.eja.eja_algebra import JordanSpinEJA
2796 This multiplication table can be verified by hand::
2798 sage: J = JordanSpinEJA(4)
2799 sage: b0,b1,b2,b3 = J.gens()
2815 We can change the generator prefix::
2817 sage: JordanSpinEJA(2, prefix='B').gens()
2822 Ensure that we have the usual inner product on `R^n`::
2824 sage: set_random_seed()
2825 sage: J = JordanSpinEJA.random_instance()
2826 sage: x,y = J.random_elements(2)
2827 sage: actual = x.inner_product(y)
2828 sage: expected = x.to_vector().inner_product(y.to_vector())
2829 sage: actual == expected
2833 def __init__(self
, n
, *args
, **kwargs
):
2834 # This is a special case of the BilinearFormEJA with the
2835 # identity matrix as its bilinear form.
2836 B
= matrix
.identity(ZZ
, n
)
2838 # Don't orthonormalize because our basis is already
2839 # orthonormal with respect to our inner-product.
2840 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2842 # But also don't pass check_field=False here, because the user
2843 # can pass in a field!
2844 super().__init
__(B
, *args
, **kwargs
)
2847 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2849 Return a random instance of this type of algebra.
2851 Needed here to override the implementation for ``BilinearFormEJA``.
2853 class_max_d
= cls
._max
_random
_instance
_dimension
()
2854 if (max_dimension
is None or max_dimension
> class_max_d
):
2855 max_dimension
= class_max_d
2856 max_size
= cls
._max
_random
_instance
_size
(max_dimension
)
2857 n
= ZZ
.random_element(max_size
+ 1)
2858 return cls(n
, **kwargs
)
2861 class TrivialEJA(RationalBasisEJA
, ConcreteEJA
):
2863 The trivial Euclidean Jordan algebra consisting of only a zero element.
2867 sage: from mjo.eja.eja_algebra import TrivialEJA
2871 sage: J = TrivialEJA()
2878 sage: 7*J.one()*12*J.one()
2880 sage: J.one().inner_product(J.one())
2882 sage: J.one().norm()
2884 sage: J.one().subalgebra_generated_by()
2885 Euclidean Jordan algebra of dimension 0 over Algebraic Real Field
2890 def __init__(self
, field
=AA
, **kwargs
):
2891 jordan_product
= lambda x
,y
: x
2892 inner_product
= lambda x
,y
: field
.zero()
2894 MS
= MatrixSpace(field
,0)
2896 # New defaults for keyword arguments
2897 if "orthonormalize" not in kwargs
: kwargs
["orthonormalize"] = False
2898 if "check_axioms" not in kwargs
: kwargs
["check_axioms"] = False
2900 super().__init
__(basis
,
2908 # The rank is zero using my definition, namely the dimension of the
2909 # largest subalgebra generated by any element.
2910 self
.rank
.set_cache(0)
2911 self
.one
.set_cache( self
.zero() )
2914 def random_instance(cls
, max_dimension
=None, *args
, **kwargs
):
2915 # We don't take a "size" argument so the superclass method is
2916 # inappropriate for us. The ``max_dimension`` argument is
2917 # included so that if this method is called generically with a
2918 # ``max_dimension=<whatever>`` argument, we don't try to pass
2919 # it on to the algebra constructor.
2920 return cls(**kwargs
)
2923 class CartesianProductEJA(FiniteDimensionalEJA
):
2925 The external (orthogonal) direct sum of two or more Euclidean
2926 Jordan algebras. Every Euclidean Jordan algebra decomposes into an
2927 orthogonal direct sum of simple Euclidean Jordan algebras which is
2928 then isometric to a Cartesian product, so no generality is lost by
2929 providing only this construction.
2933 sage: from mjo.eja.eja_algebra import (random_eja,
2934 ....: CartesianProductEJA,
2935 ....: ComplexHermitianEJA,
2937 ....: JordanSpinEJA,
2938 ....: RealSymmetricEJA)
2942 The Jordan product is inherited from our factors and implemented by
2943 our CombinatorialFreeModule Cartesian product superclass::
2945 sage: set_random_seed()
2946 sage: J1 = HadamardEJA(2)
2947 sage: J2 = RealSymmetricEJA(2)
2948 sage: J = cartesian_product([J1,J2])
2949 sage: x,y = J.random_elements(2)
2953 The ability to retrieve the original factors is implemented by our
2954 CombinatorialFreeModule Cartesian product superclass::
2956 sage: J1 = HadamardEJA(2, field=QQ)
2957 sage: J2 = JordanSpinEJA(3, field=QQ)
2958 sage: J = cartesian_product([J1,J2])
2959 sage: J.cartesian_factors()
2960 (Euclidean Jordan algebra of dimension 2 over Rational Field,
2961 Euclidean Jordan algebra of dimension 3 over Rational Field)
2963 You can provide more than two factors::
2965 sage: J1 = HadamardEJA(2)
2966 sage: J2 = JordanSpinEJA(3)
2967 sage: J3 = RealSymmetricEJA(3)
2968 sage: cartesian_product([J1,J2,J3])
2969 Euclidean Jordan algebra of dimension 2 over Algebraic Real
2970 Field (+) Euclidean Jordan algebra of dimension 3 over Algebraic
2971 Real Field (+) Euclidean Jordan algebra of dimension 6 over
2972 Algebraic Real Field
2974 Rank is additive on a Cartesian product::
2976 sage: J1 = HadamardEJA(1)
2977 sage: J2 = RealSymmetricEJA(2)
2978 sage: J = cartesian_product([J1,J2])
2979 sage: J1.rank.clear_cache()
2980 sage: J2.rank.clear_cache()
2981 sage: J.rank.clear_cache()
2984 sage: J.rank() == J1.rank() + J2.rank()
2987 The same rank computation works over the rationals, with whatever
2990 sage: J1 = HadamardEJA(1, field=QQ, orthonormalize=False)
2991 sage: J2 = RealSymmetricEJA(2, field=QQ, orthonormalize=False)
2992 sage: J = cartesian_product([J1,J2])
2993 sage: J1.rank.clear_cache()
2994 sage: J2.rank.clear_cache()
2995 sage: J.rank.clear_cache()
2998 sage: J.rank() == J1.rank() + J2.rank()
3001 The product algebra will be associative if and only if all of its
3002 components are associative::
3004 sage: J1 = HadamardEJA(2)
3005 sage: J1.is_associative()
3007 sage: J2 = HadamardEJA(3)
3008 sage: J2.is_associative()
3010 sage: J3 = RealSymmetricEJA(3)
3011 sage: J3.is_associative()
3013 sage: CP1 = cartesian_product([J1,J2])
3014 sage: CP1.is_associative()
3016 sage: CP2 = cartesian_product([J1,J3])
3017 sage: CP2.is_associative()
3020 Cartesian products of Cartesian products work::
3022 sage: J1 = JordanSpinEJA(1)
3023 sage: J2 = JordanSpinEJA(1)
3024 sage: J3 = JordanSpinEJA(1)
3025 sage: J = cartesian_product([J1,cartesian_product([J2,J3])])
3026 sage: J.multiplication_table()
3027 +----++----+----+----+
3028 | * || b0 | b1 | b2 |
3029 +====++====+====+====+
3030 | b0 || b0 | 0 | 0 |
3031 +----++----+----+----+
3032 | b1 || 0 | b1 | 0 |
3033 +----++----+----+----+
3034 | b2 || 0 | 0 | b2 |
3035 +----++----+----+----+
3036 sage: HadamardEJA(3).multiplication_table()
3037 +----++----+----+----+
3038 | * || b0 | b1 | b2 |
3039 +====++====+====+====+
3040 | b0 || b0 | 0 | 0 |
3041 +----++----+----+----+
3042 | b1 || 0 | b1 | 0 |
3043 +----++----+----+----+
3044 | b2 || 0 | 0 | b2 |
3045 +----++----+----+----+
3047 The "matrix space" of a Cartesian product always consists of
3048 ordered pairs (or triples, or...) whose components are the
3049 matrix spaces of its factors::
3051 sage: J1 = HadamardEJA(2)
3052 sage: J2 = ComplexHermitianEJA(2)
3053 sage: J = cartesian_product([J1,J2])
3054 sage: J.matrix_space()
3055 The Cartesian product of (Full MatrixSpace of 2 by 1 dense
3056 matrices over Algebraic Real Field, Module of 2 by 2 matrices
3057 with entries in Algebraic Field over the scalar ring Algebraic
3059 sage: J.one().to_matrix()[0]
3062 sage: J.one().to_matrix()[1]
3071 All factors must share the same base field::
3073 sage: J1 = HadamardEJA(2, field=QQ)
3074 sage: J2 = RealSymmetricEJA(2)
3075 sage: CartesianProductEJA((J1,J2))
3076 Traceback (most recent call last):
3078 ValueError: all factors must share the same base field
3080 The cached unit element is the same one that would be computed::
3082 sage: set_random_seed() # long time
3083 sage: J1 = random_eja() # long time
3084 sage: J2 = random_eja() # long time
3085 sage: J = cartesian_product([J1,J2]) # long time
3086 sage: actual = J.one() # long time
3087 sage: J.one.clear_cache() # long time
3088 sage: expected = J.one() # long time
3089 sage: actual == expected # long time
3092 def __init__(self
, factors
, **kwargs
):
3097 self
._sets
= factors
3099 field
= factors
[0].base_ring()
3100 if not all( J
.base_ring() == field
for J
in factors
):
3101 raise ValueError("all factors must share the same base field")
3103 # Figure out the category to use.
3104 associative
= all( f
.is_associative() for f
in factors
)
3105 category
= EuclideanJordanAlgebras(field
)
3106 if associative
: category
= category
.Associative()
3107 category
= category
.join([category
, category
.CartesianProducts()])
3109 # Compute my matrix space. We don't simply use the
3110 # ``cartesian_product()`` functor here because it acts
3111 # differently on SageMath MatrixSpaces and our custom
3112 # MatrixAlgebras, which are CombinatorialFreeModules. We
3113 # always want the result to be represented (and indexed) as an
3114 # ordered tuple. This category isn't perfect, but is good
3115 # enough for what we need to do.
3116 MS_cat
= MagmaticAlgebras(field
).FiniteDimensional().WithBasis()
3117 MS_cat
= MS_cat
.Unital().CartesianProducts()
3118 MS_factors
= tuple( J
.matrix_space() for J
in factors
)
3119 from sage
.sets
.cartesian_product
import CartesianProduct
3120 self
._matrix
_space
= CartesianProduct(MS_factors
, MS_cat
)
3122 self
._matrix
_basis
= []
3123 zero
= self
._matrix
_space
.zero()
3125 for b
in factors
[i
].matrix_basis():
3128 self
._matrix
_basis
.append(z
)
3130 self
._matrix
_basis
= tuple( self
._matrix
_space
(b
)
3131 for b
in self
._matrix
_basis
)
3132 n
= len(self
._matrix
_basis
)
3134 # We already have what we need for the super-superclass constructor.
3135 CombinatorialFreeModule
.__init
__(self
,
3142 # Now create the vector space for the algebra, which will have
3143 # its own set of non-ambient coordinates (in terms of the
3145 degree
= sum( f
._matrix
_span
.ambient_vector_space().degree()
3147 V
= VectorSpace(field
, degree
)
3148 vector_basis
= tuple( V(_all2list(b
)) for b
in self
._matrix
_basis
)
3150 # Save the span of our matrix basis (when written out as long
3151 # vectors) because otherwise we'll have to reconstruct it
3152 # every time we want to coerce a matrix into the algebra.
3153 self
._matrix
_span
= V
.span_of_basis( vector_basis
, check
=False)
3155 # Since we don't (re)orthonormalize the basis, the FDEJA
3156 # constructor is going to set self._deortho_matrix to the
3157 # identity matrix. Here we set it to the correct value using
3158 # the deortho matrices from our factors.
3159 self
._deortho
_matrix
= matrix
.block_diagonal(
3160 [J
._deortho
_matrix
for J
in factors
]
3163 self
._inner
_product
_matrix
= matrix
.block_diagonal(
3164 [J
._inner
_product
_matrix
for J
in factors
]
3166 self
._inner
_product
_matrix
._cache
= {'hermitian': True}
3167 self
._inner
_product
_matrix
.set_immutable()
3169 # Building the multiplication table is a bit more tricky
3170 # because we have to embed the entries of the factors'
3171 # multiplication tables into the product EJA.
3173 self
._multiplication
_table
= [ [zed
for j
in range(i
+1)]
3176 # Keep track of an offset that tallies the dimensions of all
3177 # previous factors. If the second factor is dim=2 and if the
3178 # first one is dim=3, then we want to skip the first 3x3 block
3179 # when copying the multiplication table for the second factor.
3182 phi_f
= self
.cartesian_embedding(f
)
3183 factor_dim
= factors
[f
].dimension()
3184 for i
in range(factor_dim
):
3185 for j
in range(i
+1):
3186 f_ij
= factors
[f
]._multiplication
_table
[i
][j
]
3188 self
._multiplication
_table
[offset
+i
][offset
+j
] = e
3189 offset
+= factor_dim
3191 self
.rank
.set_cache(sum(J
.rank() for J
in factors
))
3192 ones
= tuple(J
.one().to_matrix() for J
in factors
)
3193 self
.one
.set_cache(self(ones
))
3195 def _sets_keys(self
):
3200 sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
3201 ....: RealSymmetricEJA)
3205 The superclass uses ``_sets_keys()`` to implement its
3206 ``cartesian_factors()`` method::
3208 sage: J1 = RealSymmetricEJA(2,
3210 ....: orthonormalize=False,
3212 sage: J2 = ComplexHermitianEJA(2,field=QQ,orthonormalize=False)
3213 sage: J = cartesian_product([J1,J2])
3214 sage: x = sum(i*J.gens()[i] for i in range(len(J.gens())))
3215 sage: x.cartesian_factors()
3216 (a1 + 2*a2, 3*b0 + 4*b1 + 5*b2 + 6*b3)
3219 # Copy/pasted from CombinatorialFreeModule_CartesianProduct,
3220 # but returning a tuple instead of a list.
3221 return tuple(range(len(self
.cartesian_factors())))
3223 def cartesian_factors(self
):
3224 # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
3227 def cartesian_factor(self
, i
):
3229 Return the ``i``th factor of this algebra.
3231 return self
._sets
[i
]
3234 # Copy/pasted from CombinatorialFreeModule_CartesianProduct.
3235 from sage
.categories
.cartesian_product
import cartesian_product
3236 return cartesian_product
.symbol
.join("%s" % factor
3237 for factor
in self
._sets
)
3241 def cartesian_projection(self
, i
):
3245 sage: from mjo.eja.eja_algebra import (random_eja,
3246 ....: JordanSpinEJA,
3248 ....: RealSymmetricEJA,
3249 ....: ComplexHermitianEJA)
3253 The projection morphisms are Euclidean Jordan algebra
3256 sage: J1 = HadamardEJA(2)
3257 sage: J2 = RealSymmetricEJA(2)
3258 sage: J = cartesian_product([J1,J2])
3259 sage: J.cartesian_projection(0)
3260 Linear operator between finite-dimensional Euclidean Jordan
3261 algebras represented by the matrix:
3264 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3265 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3266 Algebraic Real Field
3267 Codomain: Euclidean Jordan algebra of dimension 2 over Algebraic
3269 sage: J.cartesian_projection(1)
3270 Linear operator between finite-dimensional Euclidean Jordan
3271 algebras represented by the matrix:
3275 Domain: Euclidean Jordan algebra of dimension 2 over Algebraic
3276 Real Field (+) Euclidean Jordan algebra of dimension 3 over
3277 Algebraic Real Field
3278 Codomain: Euclidean Jordan algebra of dimension 3 over Algebraic
3281 The projections work the way you'd expect on the vector
3282 representation of an element::
3284 sage: J1 = JordanSpinEJA(2)
3285 sage: J2 = ComplexHermitianEJA(2)
3286 sage: J = cartesian_product([J1,J2])
3287 sage: pi_left = J.cartesian_projection(0)
3288 sage: pi_right = J.cartesian_projection(1)
3289 sage: pi_left(J.one()).to_vector()
3291 sage: pi_right(J.one()).to_vector()
3293 sage: J.one().to_vector()
3298 The answer never changes::
3300 sage: set_random_seed()
3301 sage: J1 = random_eja()
3302 sage: J2 = random_eja()
3303 sage: J = cartesian_product([J1,J2])
3304 sage: P0 = J.cartesian_projection(0)
3305 sage: P1 = J.cartesian_projection(0)
3310 offset
= sum( self
.cartesian_factor(k
).dimension()
3312 Ji
= self
.cartesian_factor(i
)
3313 Pi
= self
._module
_morphism
(lambda j
: Ji
.monomial(j
- offset
),
3316 return FiniteDimensionalEJAOperator(self
,Ji
,Pi
.matrix())
3319 def cartesian_embedding(self
, i
):
3323 sage: from mjo.eja.eja_algebra import (random_eja,
3324 ....: JordanSpinEJA,
3326 ....: RealSymmetricEJA)
3330 The embedding morphisms are Euclidean Jordan algebra
3333 sage: J1 = HadamardEJA(2)
3334 sage: J2 = RealSymmetricEJA(2)
3335 sage: J = cartesian_product([J1,J2])
3336 sage: J.cartesian_embedding(0)
3337 Linear operator between finite-dimensional Euclidean Jordan
3338 algebras represented by the matrix:
3344 Domain: Euclidean Jordan algebra of dimension 2 over
3345 Algebraic Real Field
3346 Codomain: Euclidean Jordan algebra of dimension 2 over
3347 Algebraic Real Field (+) Euclidean Jordan algebra of
3348 dimension 3 over Algebraic Real Field
3349 sage: J.cartesian_embedding(1)
3350 Linear operator between finite-dimensional Euclidean Jordan
3351 algebras represented by the matrix:
3357 Domain: Euclidean Jordan algebra of dimension 3 over
3358 Algebraic Real Field
3359 Codomain: Euclidean Jordan algebra of dimension 2 over
3360 Algebraic Real Field (+) Euclidean Jordan algebra of
3361 dimension 3 over Algebraic Real Field
3363 The embeddings work the way you'd expect on the vector
3364 representation of an element::
3366 sage: J1 = JordanSpinEJA(3)
3367 sage: J2 = RealSymmetricEJA(2)
3368 sage: J = cartesian_product([J1,J2])
3369 sage: iota_left = J.cartesian_embedding(0)
3370 sage: iota_right = J.cartesian_embedding(1)
3371 sage: iota_left(J1.zero()) == J.zero()
3373 sage: iota_right(J2.zero()) == J.zero()
3375 sage: J1.one().to_vector()
3377 sage: iota_left(J1.one()).to_vector()
3379 sage: J2.one().to_vector()
3381 sage: iota_right(J2.one()).to_vector()
3383 sage: J.one().to_vector()
3388 The answer never changes::
3390 sage: set_random_seed()
3391 sage: J1 = random_eja()
3392 sage: J2 = random_eja()
3393 sage: J = cartesian_product([J1,J2])
3394 sage: E0 = J.cartesian_embedding(0)
3395 sage: E1 = J.cartesian_embedding(0)
3399 Composing a projection with the corresponding inclusion should
3400 produce the identity map, and mismatching them should produce
3403 sage: set_random_seed()
3404 sage: J1 = random_eja()
3405 sage: J2 = random_eja()
3406 sage: J = cartesian_product([J1,J2])
3407 sage: iota_left = J.cartesian_embedding(0)
3408 sage: iota_right = J.cartesian_embedding(1)
3409 sage: pi_left = J.cartesian_projection(0)
3410 sage: pi_right = J.cartesian_projection(1)
3411 sage: pi_left*iota_left == J1.one().operator()
3413 sage: pi_right*iota_right == J2.one().operator()
3415 sage: (pi_left*iota_right).is_zero()
3417 sage: (pi_right*iota_left).is_zero()
3421 offset
= sum( self
.cartesian_factor(k
).dimension()
3423 Ji
= self
.cartesian_factor(i
)
3424 Ei
= Ji
._module
_morphism
(lambda j
: self
.monomial(j
+ offset
),
3426 return FiniteDimensionalEJAOperator(Ji
,self
,Ei
.matrix())
3430 FiniteDimensionalEJA
.CartesianProduct
= CartesianProductEJA
3432 class RationalBasisCartesianProductEJA(CartesianProductEJA
,
3435 A separate class for products of algebras for which we know a
3440 sage: from mjo.eja.eja_algebra import (FiniteDimensionalEJA,
3442 ....: JordanSpinEJA,
3443 ....: RealSymmetricEJA)
3447 This gives us fast characteristic polynomial computations in
3448 product algebras, too::
3451 sage: J1 = JordanSpinEJA(2)
3452 sage: J2 = RealSymmetricEJA(3)
3453 sage: J = cartesian_product([J1,J2])
3454 sage: J.characteristic_polynomial_of().degree()
3461 The ``cartesian_product()`` function only uses the first factor to
3462 decide where the result will live; thus we have to be careful to
3463 check that all factors do indeed have a ``rational_algebra()`` method
3464 before we construct an algebra that claims to have a rational basis::
3466 sage: J1 = HadamardEJA(2)
3467 sage: jp = lambda X,Y: X*Y
3468 sage: ip = lambda X,Y: X[0,0]*Y[0,0]
3469 sage: b1 = matrix(QQ, [[1]])
3470 sage: J2 = FiniteDimensionalEJA((b1,), jp, ip)
3471 sage: cartesian_product([J2,J1]) # factor one not RationalBasisEJA
3472 Euclidean Jordan algebra of dimension 1 over Algebraic Real
3473 Field (+) Euclidean Jordan algebra of dimension 2 over Algebraic
3475 sage: cartesian_product([J1,J2]) # factor one is RationalBasisEJA
3476 Traceback (most recent call last):
3478 ValueError: factor not a RationalBasisEJA
3481 def __init__(self
, algebras
, **kwargs
):
3482 if not all( hasattr(r
, "rational_algebra") for r
in algebras
):
3483 raise ValueError("factor not a RationalBasisEJA")
3485 CartesianProductEJA
.__init
__(self
, algebras
, **kwargs
)
3488 def rational_algebra(self
):
3489 if self
.base_ring() is QQ
:
3492 return cartesian_product([
3493 r
.rational_algebra() for r
in self
.cartesian_factors()
3497 RationalBasisEJA
.CartesianProduct
= RationalBasisCartesianProductEJA
3499 def random_eja(max_dimension
=None, *args
, **kwargs
):
3504 sage: from mjo.eja.eja_algebra import random_eja
3508 sage: set_random_seed()
3509 sage: n = ZZ.random_element(1,5)
3510 sage: J = random_eja(max_dimension=n, field=QQ, orthonormalize=False)
3511 sage: J.dimension() <= n
3515 # Use the ConcreteEJA default as the total upper bound (regardless
3516 # of any whether or not any individual factors set a lower limit).
3517 if max_dimension
is None:
3518 max_dimension
= ConcreteEJA
._max
_random
_instance
_dimension
()
3519 J1
= ConcreteEJA
.random_instance(max_dimension
, *args
, **kwargs
)
3522 # Roll the dice to see if we attempt a Cartesian product.
3523 dice_roll
= ZZ
.random_element(len(ConcreteEJA
.__subclasses
__()) + 1)
3524 new_max_dimension
= max_dimension
- J1
.dimension()
3525 if new_max_dimension
== 0 or dice_roll
!= 0:
3526 # If it's already as big as we're willing to tolerate, just
3527 # return it and don't worry about Cartesian products.
3530 # Use random_eja() again so we can get more than two factors
3531 # if the sub-call also Decides on a cartesian product.
3532 J2
= random_eja(new_max_dimension
, *args
, **kwargs
)
3533 return cartesian_product([J1
,J2
])