-# -*- coding: utf-8 -*-
-
from sage.matrix.constructor import matrix
from sage.modules.free_module import VectorSpace
from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
-# TODO: make this unnecessary somehow.
-from sage.misc.lazy_import import lazy_import
-lazy_import('mjo.eja.eja_algebra', 'FiniteDimensionalEuclideanJordanAlgebra')
-lazy_import('mjo.eja.eja_element_subalgebra',
- 'FiniteDimensionalEuclideanJordanElementSubalgebra')
from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
from mjo.eja.eja_utils import _mat2vec
We should always get back an element of the algebra::
sage: set_random_seed()
- sage: p = PolynomialRing(QQ, 't').random_element()
+ sage: p = PolynomialRing(AA, 't').random_element()
sage: J = random_eja()
sage: x = J.random_element()
sage: x.apply_univariate_polynomial(p) in J
True
"""
- p = self.parent().characteristic_polynomial()
+ p = self.parent().characteristic_polynomial_of()
return p(*self.to_vector())
SETUP::
sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+ ....: TrivialEJA,
....: random_eja)
EXAMPLES::
sage: x.det()
-1
+ The determinant of the sole element in the rank-zero trivial
+ algebra is ``1``, by three paths of reasoning. First, its
+ characteristic polynomial is a constant ``1``, so the constant
+ term in that polynomial is ``1``. Second, the characteristic
+ polynomial evaluated at zero is again ``1``. And finally, the
+ (empty) product of its eigenvalues is likewise just unity::
+
+ sage: J = TrivialEJA()
+ sage: J.zero().det()
+ 1
+
TESTS:
An element is invertible if and only if its determinant is
sage: x,y = J.random_elements(2)
sage: (x*y).det() == x.det()*y.det()
True
-
"""
P = self.parent()
r = P.rank()
- p = P._charpoly_coeff(0)
- # The _charpoly_coeff function already adds the factor of
- # -1 to ensure that _charpoly_coeff(0) is really what
- # appears in front of t^{0} in the charpoly. However,
- # we want (-1)^r times THAT for the determinant.
+
+ if r == 0:
+ # Special case, since we don't get the a0=1
+ # coefficient when the rank of the algebra
+ # is zero.
+ return P.base_ring().one()
+
+ p = P._charpoly_coefficients()[0]
+ # The _charpoly_coeff function already adds the factor of -1
+ # to ensure that _charpoly_coefficients()[0] is really what
+ # appears in front of t^{0} in the charpoly. However, we want
+ # (-1)^r times THAT for the determinant.
return ((-1)**r)*p(*self.to_vector())
sage: while not x.is_invertible():
....: x = J.random_element()
sage: x_vec = x.to_vector()
- sage: x0 = x_vec[0]
+ sage: x0 = x_vec[:1]
sage: x_bar = x_vec[1:]
- sage: coeff = ~(x0^2 - x_bar.inner_product(x_bar))
- sage: inv_vec = x_vec.parent()([x0] + (-x_bar).list())
- sage: x_inverse = coeff*inv_vec
+ sage: coeff = x0.inner_product(x0) - x_bar.inner_product(x_bar)
+ sage: x_inverse = x_vec.parent()(x0.list() + (-x_bar).list())
+ sage: if not coeff.is_zero(): x_inverse = x_inverse/coeff
sage: x.inverse() == J.from_vector(x_inverse)
True
if not self.is_invertible():
raise ValueError("element is not invertible")
+ if self.parent()._charpoly_coefficients.is_in_cache():
+ # We can invert using our charpoly if it will be fast to
+ # compute. If the coefficients are cached, our rank had
+ # better be too!
+ r = self.parent().rank()
+ a = self.characteristic_polynomial().coefficients(sparse=False)
+ return (-1)**(r+1)*sum(a[i+1]*self**i for i in range(r))/self.det()
+
return (~self.quadratic_representation())(self)
whether or not the paren't algebra's zero element is a root
of this element's minimal polynomial.
+ That is... unless the coefficients of our algebra's
+ "characteristic polynomial of" function are already cached!
+ In that case, we just use the determinant (which will be fast
+ as a result).
+
Beware that we can't use the superclass method, because it
relies on the algebra being associative.
else:
return False
+ if self.parent()._charpoly_coefficients.is_in_cache():
+ # The determinant will be quicker than computing the minimal
+ # polynomial from scratch, most likely.
+ return (not self.det().is_zero())
+
# In fact, we only need to know if the constant term is non-zero,
# so we can pass in the field's zero element instead.
zero = self.base_ring().zero()
The spectral decomposition of a non-regular element should always
contain at least one non-minimal idempotent::
- sage: J = RealSymmetricEJA(3, AA)
+ sage: J = RealSymmetricEJA(3)
sage: x = sum(J.gens())
sage: x.is_regular()
False
On the other hand, the spectral decomposition of a regular
element should always be in terms of minimal idempotents::
- sage: J = JordanSpinEJA(4, AA)
+ sage: J = JordanSpinEJA(4)
sage: x = sum( i*J.gens()[i] for i in range(len(J.gens())) )
sage: x.is_regular()
True
sage: set_random_seed()
sage: J = JordanSpinEJA.random_instance()
+ sage: n = J.dimension()
sage: x = J.random_element()
- sage: x == x.coefficient(0)*J.one() or x.degree() == 2
+ sage: x.degree() == min(n,2) or (x == x.coefficient(0)*J.one())
True
TESTS:
TESTS:
The minimal polynomial of the identity and zero elements are
- always the same::
+ always the same, except in trivial algebras where the minimal
+ polynomial of the unit/zero element is ``1``::
sage: set_random_seed()
- sage: J = random_eja(nontrivial=True)
- sage: J.one().minimal_polynomial()
+ sage: J = random_eja()
+ sage: mu = J.one().minimal_polynomial()
+ sage: t = mu.parent().gen()
+ sage: mu + int(J.is_trivial())*(t-2)
t - 1
- sage: J.zero().minimal_polynomial()
+ sage: mu = J.zero().minimal_polynomial()
+ sage: t = mu.parent().gen()
+ sage: mu + int(J.is_trivial())*(t-1)
t
The degree of an element is (by one definition) the degree
two here so that said elements actually exist::
sage: set_random_seed()
- sage: n_max = max(2, JordanSpinEJA._max_test_case_size())
+ sage: n_max = max(2, JordanSpinEJA._max_random_instance_size())
sage: n = ZZ.random_element(2, n_max)
sage: J = JordanSpinEJA(n)
sage: y = J.random_element()
and in particular, a re-scaling of the basis::
sage: set_random_seed()
- sage: n_max = RealSymmetricEJA._max_test_case_size()
+ sage: n_max = RealSymmetricEJA._max_random_instance_size()
sage: n = ZZ.random_element(1, n_max)
- sage: J1 = RealSymmetricEJA(n,QQ)
- sage: J2 = RealSymmetricEJA(n,QQ,normalize_basis=False)
- sage: X = random_matrix(QQ,n)
+ sage: J1 = RealSymmetricEJA(n)
+ sage: J2 = RealSymmetricEJA(n,normalize_basis=False)
+ sage: X = random_matrix(AA,n)
sage: X = X*X.transpose()
sage: x1 = J1(X)
sage: x2 = J2(X)
# in the "normal" case without us having to think about it.
return self.operator().minimal_polynomial()
- A = self.subalgebra_generated_by()
+ A = self.subalgebra_generated_by(orthonormalize_basis=False)
return A(self).operator().minimal_polynomial()
- def natural_representation(self):
+ def to_matrix(self):
"""
- Return a more-natural representation of this element.
+ Return an (often more natural) representation of this element as a
+ matrix.
- Every finite-dimensional Euclidean Jordan Algebra is a
- direct sum of five simple algebras, four of which comprise
- Hermitian matrices. This method returns the original
- "natural" representation of this element as a Hermitian
- matrix, if it has one. If not, you get the usual representation.
+ Every finite-dimensional Euclidean Jordan Algebra is a direct
+ sum of five simple algebras, four of which comprise Hermitian
+ matrices. This method returns a "natural" matrix
+ representation of this element as either a Hermitian matrix or
+ column vector.
SETUP::
sage: J = ComplexHermitianEJA(3)
sage: J.one()
e0 + e3 + e8
- sage: J.one().natural_representation()
+ sage: J.one().to_matrix()
[1 0 0 0 0 0]
[0 1 0 0 0 0]
[0 0 1 0 0 0]
sage: J = QuaternionHermitianEJA(3)
sage: J.one()
e0 + e5 + e14
- sage: J.one().natural_representation()
+ sage: J.one().to_matrix()
[1 0 0 0 0 0 0 0 0 0 0 0]
[0 1 0 0 0 0 0 0 0 0 0 0]
[0 0 1 0 0 0 0 0 0 0 0 0]
[0 0 0 0 0 0 0 0 0 1 0 0]
[0 0 0 0 0 0 0 0 0 0 1 0]
[0 0 0 0 0 0 0 0 0 0 0 1]
-
"""
- B = self.parent().natural_basis()
- W = self.parent().natural_basis_space()
- return W.linear_combination(zip(B,self.to_vector()))
+ B = self.parent().matrix_basis()
+ W = self.parent().matrix_space()
+
+ # This is just a manual "from_vector()", but of course
+ # matrix spaces aren't vector spaces in sage, so they
+ # don't have a from_vector() method.
+ return W.linear_combination( zip(B, self.to_vector()) )
def norm(self):
sage: J = HadamardEJA(2)
sage: x = sum(J.gens())
sage: x.norm()
- sqrt(2)
+ 1.414213562373095?
::
sage: set_random_seed()
sage: x = JordanSpinEJA.random_instance().random_element()
sage: x_vec = x.to_vector()
+ sage: Q = matrix.identity(x.base_ring(), 0)
sage: n = x_vec.degree()
- sage: x0 = x_vec[0]
- sage: x_bar = x_vec[1:]
- sage: A = matrix(QQ, 1, [x_vec.inner_product(x_vec)])
- sage: B = 2*x0*x_bar.row()
- sage: C = 2*x0*x_bar.column()
- sage: D = matrix.identity(QQ, n-1)
- sage: D = (x0^2 - x_bar.inner_product(x_bar))*D
- sage: D = D + 2*x_bar.tensor_product(x_bar)
- sage: Q = matrix.block(2,2,[A,B,C,D])
+ sage: if n > 0:
+ ....: x0 = x_vec[0]
+ ....: x_bar = x_vec[1:]
+ ....: A = matrix(x.base_ring(), 1, [x_vec.inner_product(x_vec)])
+ ....: B = 2*x0*x_bar.row()
+ ....: C = 2*x0*x_bar.column()
+ ....: D = matrix.identity(x.base_ring(), n-1)
+ ....: D = (x0^2 - x_bar.inner_product(x_bar))*D
+ ....: D = D + 2*x_bar.tensor_product(x_bar)
+ ....: Q = matrix.block(2,2,[A,B,C,D])
sage: Q == x.quadratic_representation().matrix()
True
The spectral decomposition of the identity is ``1`` times itself,
and the spectral decomposition of zero is ``0`` times the identity::
- sage: J = RealSymmetricEJA(3,AA)
+ sage: J = RealSymmetricEJA(3)
sage: J.one()
e0 + e2 + e5
sage: J.one().spectral_decomposition()
TESTS::
- sage: J = RealSymmetricEJA(4,AA)
+ sage: J = RealSymmetricEJA(4)
sage: x = sum(J.gens())
sage: sd = x.spectral_decomposition()
sage: l0 = sd[0][0]
sage: l0*c0 + l1*c1 == x
True
+ The spectral decomposition should work in subalgebras, too::
+
+ sage: J = RealSymmetricEJA(4)
+ sage: (e0, e1, e2, e3, e4, e5, e6, e7, e8, e9) = J.gens()
+ sage: A = 2*e5 - 2*e8
+ sage: (lambda1, c1) = A.spectral_decomposition()[1]
+ sage: (J0, J5, J1) = J.peirce_decomposition(c1)
+ sage: (f0, f1, f2) = J1.gens()
+ sage: f0.spectral_decomposition()
+ [(0, f2), (1, f0)]
+
"""
- P = self.parent()
A = self.subalgebra_generated_by(orthonormalize_basis=True)
result = []
for (evalue, proj) in A(self).operator().spectral_decomposition():
True
"""
+ from mjo.eja.eja_element_subalgebra import FiniteDimensionalEuclideanJordanElementSubalgebra
return FiniteDimensionalEuclideanJordanElementSubalgebra(self, orthonormalize_basis)
TESTS:
Ensure that we can find an idempotent in a non-trivial algebra
- where there are non-nilpotent elements::
+ where there are non-nilpotent elements, or that we get the dumb
+ solution in the trivial algebra::
sage: set_random_seed()
- sage: J = random_eja(nontrivial=True)
+ sage: J = random_eja()
sage: x = J.random_element()
- sage: while x.is_nilpotent():
+ sage: while x.is_nilpotent() and not J.is_trivial():
....: x = J.random_element()
sage: c = x.subalgebra_idempotent()
sage: c^2 == c
# the trace is an empty sum.
return P.base_ring().zero()
- p = P._charpoly_coeff(r-1)
+ p = P._charpoly_coefficients()[r-1]
# The _charpoly_coeff function already adds the factor of
# -1 to ensure that _charpoly_coeff(r-1) is really what
# appears in front of t^{r-1} in the charpoly. However,
sage: J = HadamardEJA(2)
sage: x = sum(J.gens())
sage: x.trace_norm()
- sqrt(2)
+ 1.414213562373095?
::
sage: J = JordanSpinEJA(4)
sage: x = sum(J.gens())
sage: x.trace_norm()
- 2*sqrt(2)
+ 2.828427124746190?
"""
return self.trace_inner_product(self).sqrt()
projects / sage.d.git / blobdiff