]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_element.py
eja: move eja_subalgebra to eja_element_subalgebra.
[sage.d.git] / mjo / eja / eja_element.py
index b8db12e9827de13720c8ac1dae97ec424edbd59c..5b8bc1e98a58ef03234e8a966a75eefd859c461c 100644 (file)
@@ -1,3 +1,7 @@
+# -*- coding: utf-8 -*-
+
+from itertools import izip
+
 from sage.matrix.constructor import matrix
 from sage.modules.free_module import VectorSpace
 from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
@@ -5,7 +9,7 @@ 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_subalgebra',
+lazy_import('mjo.eja.eja_element_subalgebra',
             'FiniteDimensionalEuclideanJordanElementSubalgebra')
 from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
 from mjo.eja.eja_utils import _mat2vec
@@ -32,7 +36,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         Return ``self`` raised to the power ``n``.
 
         Jordan algebras are always power-associative; see for
-        example Faraut and Koranyi, Proposition II.1.2 (ii).
+        example Faraut and Korányi, Proposition II.1.2 (ii).
 
         We have to override this because our superclass uses row
         vectors instead of column vectors! We, on the other hand,
@@ -373,7 +377,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             True
 
         Ensure that the determinant is multiplicative on an associative
-        subalgebra as in Faraut and Koranyi's Proposition II.2.2::
+        subalgebra as in Faraut and Korányi's Proposition II.2.2::
 
             sage: set_random_seed()
             sage: J = random_eja().random_element().subalgebra_generated_by()
@@ -403,7 +407,8 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
         SETUP::
 
-            sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+            sage: from mjo.eja.eja_algebra import (ComplexHermitianEJA,
+            ....:                                  JordanSpinEJA,
             ....:                                  random_eja)
 
         EXAMPLES:
@@ -425,6 +430,13 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: x.inverse() == J.from_vector(x_inverse)
             True
 
+        Trying to invert a non-invertible element throws an error:
+
+            sage: JordanSpinEJA(3).zero().inverse()
+            Traceback (most recent call last):
+            ...
+            ValueError: element is not invertible
+
         TESTS:
 
         The identity element is its own inverse::
@@ -450,13 +462,32 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: (not x.is_invertible()) or (x.inverse().inverse() == x)
             True
 
-        The zero element is never invertible::
+        Proposition II.2.3 in Faraut and Korányi says that the inverse
+        of an element is the inverse of its left-multiplication operator
+        applied to the algebra's identity, when that inverse exists::
 
             sage: set_random_seed()
-            sage: J = random_eja().zero().inverse()
-            Traceback (most recent call last):
-            ...
-            ValueError: element is not invertible
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: (not x.operator().is_invertible()) or (
+            ....:    x.operator().inverse()(J.one()) == x.inverse() )
+            True
+
+        Proposition II.2.4 in Faraut and Korányi gives a formula for
+        the inverse based on the characteristic polynomial and the
+        Cayley-Hamilton theorem for Euclidean Jordan algebras::
+
+            sage: set_random_seed()
+            sage: J = ComplexHermitianEJA(3)
+            sage: x = J.random_element()
+            sage: while not x.is_invertible():
+            ....:     x = J.random_element()
+            sage: r = J.rank()
+            sage: a = x.characteristic_polynomial().coefficients(sparse=False)
+            sage: expected  = (-1)^(r+1)/x.det()
+            sage: expected *= sum( a[i+1]*x^i for i in range(r) )
+            sage: x.inverse() == expected
+            True
 
         """
         if not self.is_invertible():
@@ -515,6 +546,96 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         return not (p(zero) == zero)
 
 
+    def is_minimal_idempotent(self):
+        """
+        Return whether or not this element is a minimal idempotent.
+
+
+        An element of a Euclidean Jordan algebra is a minimal idempotent
+        if it :meth:`is_idempotent` and if its Peirce subalgebra
+        corresponding to the eigenvalue ``1`` has dimension ``1`` (Baes,
+        Proposition 2.7.17).
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+            ....:                                  RealSymmetricEJA,
+            ....:                                  random_eja)
+
+        WARNING::
+
+        This method is sloooooow.
+
+        EXAMPLES:
+
+        The spectral decomposition of a non-regular element should always
+        contain at least one non-minimal idempotent::
+
+            sage: J = RealSymmetricEJA(3, AA)
+            sage: x = sum(J.gens())
+            sage: x.is_regular()
+            False
+            sage: [ c.is_minimal_idempotent()
+            ....:   for (l,c) in x.spectral_decomposition() ]
+            [False, True]
+
+        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: x = sum( i*J.gens()[i] for i in range(len(J.gens())) )
+            sage: x.is_regular()
+            True
+            sage: [ c.is_minimal_idempotent()
+            ....:   for (l,c) in x.spectral_decomposition() ]
+            [True, True]
+
+        TESTS:
+
+        The identity element is minimal only in an EJA of rank one::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: J.rank() == 1 or not J.one().is_minimal_idempotent()
+            True
+
+        A non-idempotent cannot be a minimal idempotent::
+
+            sage: set_random_seed()
+            sage: J = JordanSpinEJA(4)
+            sage: x = J.random_element()
+            sage: (not x.is_idempotent()) and x.is_minimal_idempotent()
+            False
+
+        Proposition 2.7.19 in Baes says that an element is a minimal
+        idempotent if and only if it's idempotent with trace equal to
+        unity::
+
+            sage: set_random_seed()
+            sage: J = JordanSpinEJA(4)
+            sage: x = J.random_element()
+            sage: expected = (x.is_idempotent() and x.trace() == 1)
+            sage: actual = x.is_minimal_idempotent()
+            sage: actual == expected
+            True
+
+        """
+        # TODO: when the Peirce decomposition is implemented for real,
+        # we can use that instead of finding this eigenspace manually.
+        #
+        # Trivial eigenspaces don't appear in the list, so we default to the
+        # trivial one and override it if there's a nontrivial space in the
+        # list.
+        if not self.is_idempotent():
+            return False
+
+        J1 = VectorSpace(self.parent().base_ring(), 0)
+        for (eigval, eigspace) in self.operator().matrix().left_eigenspaces():
+            if eigval == 1:
+                J1 = eigspace
+        return (J1.dimension() == 1)
+
+
     def is_nilpotent(self):
         """
         Return whether or not some power of this element is zero.
@@ -542,10 +663,11 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
         TESTS:
 
-        The identity element is never nilpotent::
+        The identity element is never nilpotent, except in a trivial EJA::
 
             sage: set_random_seed()
-            sage: random_eja().one().is_nilpotent()
+            sage: J = random_eja()
+            sage: J.one().is_nilpotent() and not J.is_trivial()
             False
 
         The additive identity is always nilpotent::
@@ -589,11 +711,11 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         TESTS:
 
         The zero element should never be regular, unless the parent
-        algebra has dimension one::
+        algebra has dimension less than or equal to one::
 
             sage: set_random_seed()
             sage: J = random_eja()
-            sage: J.dimension() == 1 or not J.zero().is_regular()
+            sage: J.dimension() <= 1 or not J.zero().is_regular()
             True
 
         The unit element isn't regular unless the algebra happens to
@@ -601,7 +723,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
             sage: set_random_seed()
             sage: J = random_eja()
-            sage: J.dimension() == 1 or not J.one().is_regular()
+            sage: J.dimension() <= 1 or not J.one().is_regular()
             True
 
         """
@@ -645,14 +767,17 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
         TESTS:
 
-        The zero and unit elements are both of degree one::
+        The zero and unit elements are both of degree one in nontrivial
+        algebras::
 
             sage: set_random_seed()
             sage: J = random_eja()
-            sage: J.zero().degree()
-            1
-            sage: J.one().degree()
-            1
+            sage: d = J.zero().degree()
+            sage: (J.is_trivial() and d == 0) or d == 1
+            True
+            sage: d = J.one().degree()
+            sage: (J.is_trivial() and d == 0) or d == 1
+            True
 
         Our implementation agrees with the definition::
 
@@ -696,15 +821,29 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
             ....:                                  RealSymmetricEJA,
+            ....:                                  TrivialEJA,
             ....:                                  random_eja)
 
+        EXAMPLES:
+
+        Keeping in mind that the polynomial ``1`` evaluates the identity
+        element (also the zero element) of the trivial algebra, it is clear
+        that the polynomial ``1`` is the minimal polynomial of the only
+        element in a trivial algebra::
+
+            sage: J = TrivialEJA()
+            sage: J.one().minimal_polynomial()
+            1
+            sage: J.zero().minimal_polynomial()
+            1
+
         TESTS:
 
         The minimal polynomial of the identity and zero elements are
         always the same::
 
             sage: set_random_seed()
-            sage: J = random_eja()
+            sage: J = random_eja(nontrivial=True)
             sage: J.one().minimal_polynomial()
             t - 1
             sage: J.zero().minimal_polynomial()
@@ -754,7 +893,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: n_max = RealSymmetricEJA._max_test_case_size()
             sage: n = ZZ.random_element(1, n_max)
             sage: J1 = RealSymmetricEJA(n,QQ)
-            sage: J2 = RealSymmetricEJA(n,QQ,False)
+            sage: J2 = RealSymmetricEJA(n,QQ,normalize_basis=False)
             sage: X = random_matrix(QQ,n)
             sage: X = X*X.transpose()
             sage: x1 = J1(X)
@@ -830,7 +969,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         """
         B = self.parent().natural_basis()
         W = self.parent().natural_basis_space()
-        return W.linear_combination(zip(B,self.to_vector()))
+        return W.linear_combination(izip(B,self.to_vector()))
 
 
     def norm(self):
@@ -1008,12 +1147,82 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
 
 
+    def spectral_decomposition(self):
+        """
+        Return the unique spectral decomposition of this element.
+
+        ALGORITHM:
+
+        Following Faraut and Korányi's Theorem III.1.1, we restrict this
+        element's left-multiplication-by operator to the subalgebra it
+        generates. We then compute the spectral decomposition of that
+        operator, and the spectral projectors we get back must be the
+        left-multiplication-by operators for the idempotents we
+        seek. Thus applying them to the identity element gives us those
+        idempotents.
+
+        Since the eigenvalues are required to be distinct, we take
+        the spectral decomposition of the zero element to be zero
+        times the identity element of the algebra (which is idempotent,
+        obviously).
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import RealSymmetricEJA
+
+        EXAMPLES:
+
+        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.one()
+            e0 + e2 + e5
+            sage: J.one().spectral_decomposition()
+            [(1, e0 + e2 + e5)]
+            sage: J.zero().spectral_decomposition()
+            [(0, e0 + e2 + e5)]
+
+        TESTS::
+
+            sage: J = RealSymmetricEJA(4,AA)
+            sage: x = sum(J.gens())
+            sage: sd = x.spectral_decomposition()
+            sage: l0 = sd[0][0]
+            sage: l1 = sd[1][0]
+            sage: c0 = sd[0][1]
+            sage: c1 = sd[1][1]
+            sage: c0.inner_product(c1) == 0
+            True
+            sage: c0.is_idempotent()
+            True
+            sage: c1.is_idempotent()
+            True
+            sage: c0 + c1 == J.one()
+            True
+            sage: l0*c0 + l1*c1 == x
+            True
 
-    def subalgebra_generated_by(self):
+        """
+        P = self.parent()
+        A = self.subalgebra_generated_by(orthonormalize_basis=True)
+        result = []
+        for (evalue, proj) in A(self).operator().spectral_decomposition():
+            result.append( (evalue, proj(A.one()).superalgebra_element()) )
+        return result
+
+    def subalgebra_generated_by(self, orthonormalize_basis=False):
         """
         Return the associative subalgebra of the parent EJA generated
         by this element.
 
+        Since our parent algebra is unital, we want "subalgebra" to mean
+        "unital subalgebra" as well; thus the subalgebra that an element
+        generates will itself be a Euclidean Jordan algebra after
+        restricting the algebra operations appropriately. This is the
+        subalgebra that Faraut and Korányi work with in section II.2, for
+        example.
+
         SETUP::
 
             sage: from mjo.eja.eja_algebra import random_eja
@@ -1038,17 +1247,18 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             sage: A(x^2) == A(x)*A(x)
             True
 
-        The subalgebra generated by the zero element is trivial::
+        By definition, the subalgebra generated by the zero element is
+        the one-dimensional algebra generated by the identity
+        element... unless the original algebra was trivial, in which
+        case the subalgebra is trivial too::
 
             sage: set_random_seed()
             sage: A = random_eja().zero().subalgebra_generated_by()
-            sage: A
-            Euclidean Jordan algebra of dimension 0 over...
-            sage: A.one()
-            0
+            sage: (A.is_trivial() and A.dimension() == 0) or A.dimension() == 1
+            True
 
         """
-        return FiniteDimensionalEuclideanJordanElementSubalgebra(self)
+        return FiniteDimensionalEuclideanJordanElementSubalgebra(self, orthonormalize_basis)
 
 
     def subalgebra_idempotent(self):
@@ -1111,14 +1321,23 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         """
         Return my trace, the sum of my eigenvalues.
 
+        In a trivial algebra, however you want to look at it, the trace is
+        an empty sum for which we declare the result to be zero.
+
         SETUP::
 
             sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
             ....:                                  RealCartesianProductEJA,
+            ....:                                  TrivialEJA,
             ....:                                  random_eja)
 
         EXAMPLES::
 
+            sage: J = TrivialEJA()
+            sage: J.zero().trace()
+            0
+
+        ::
             sage: J = JordanSpinEJA(3)
             sage: x = sum(J.gens())
             sage: x.trace()
@@ -1142,6 +1361,12 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
         """
         P = self.parent()
         r = P.rank()
+
+        if r == 0:
+            # Special case for the trivial algebra where
+            # the trace is an empty sum.
+            return P.base_ring().zero()
+
         p = P._charpoly_coeff(r-1)
         # The _charpoly_coeff function already adds the factor of
         # -1 to ensure that _charpoly_coeff(r-1) is really what
@@ -1160,8 +1385,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
 
         TESTS:
 
-        The trace inner product is commutative, bilinear, and satisfies
-        the Jordan axiom:
+        The trace inner product is commutative, bilinear, and associative::
 
             sage: set_random_seed()
             sage: J = random_eja()
@@ -1181,7 +1405,7 @@ class FiniteDimensionalEuclideanJordanAlgebraElement(IndexedFreeModuleElement):
             ....:              a*x.trace_inner_product(z) )
             sage: actual == expected
             True
-            sage: # jordan axiom
+            sage: # associative
             sage: (x*y).trace_inner_product(z) == y.trace_inner_product(x*z)
             True