]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_element.py
eja: drop "FiniteDimensional" prefix everywhere.
[sage.d.git] / mjo / eja / eja_element.py
index 0983ea0607a438464b3e9ab6245a344ba59a5aef..8af3b77698457486af3d5b64eb52bf6f6e7b5ea8 100644 (file)
@@ -3,11 +3,11 @@ from sage.misc.cachefunc import cached_method
 from sage.modules.free_module import VectorSpace
 from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
 
-from mjo.eja.eja_operator import FiniteDimensionalEJAOperator
+from mjo.eja.eja_operator import EJAOperator
 from mjo.eja.eja_utils import _scale
 
 
-class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
+class EJAElement(IndexedFreeModuleElement):
     """
     An element of a Euclidean Jordan algebra.
     """
@@ -795,7 +795,23 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
 
         ALGORITHM:
 
-        .........
+        First we handle the special cases where the algebra is
+        trivial, this element is zero, or the dimension of the algebra
+        is one and this element is not zero. With those out of the
+        way, we may assume that ``self`` is nonzero, the algebra is
+        nontrivial, and that the dimension of the algebra is at least
+        two.
+
+        Beginning with the algebra's unit element (power zero), we add
+        successive (basis representations of) powers of this element
+        to a matrix, row-reducing at each step. After row-reducing, we
+        check the rank of the matrix. If adding a row and row-reducing
+        does not increase the rank of the matrix at any point, the row
+        we've just added lives in the span of the previous ones; thus
+        the corresponding power of ``self`` lives in the span of its
+        lesser powers. When that happens, the degree of the minimal
+        polynomial is the rank of the matrix; if it never happens, the
+        degree must be the dimension of the entire space.
 
         SETUP::
 
@@ -838,7 +854,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             sage: x = random_eja().random_element()
             sage: x.degree() == x.minimal_polynomial().degree()
             True
-
         """
         n = self.parent().dimension()
 
@@ -1147,7 +1162,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         P = self.parent()
         left_mult_by_self = lambda y: self*y
         L = P.module_morphism(function=left_mult_by_self, codomain=P)
-        return FiniteDimensionalEJAOperator(P, P, L.matrix() )
+        return EJAOperator(P, P, L.matrix() )
 
 
     def quadratic_representation(self, other=None):
@@ -1440,7 +1455,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         if self.is_nilpotent():
             raise ValueError("this only works with non-nilpotent elements!")
 
-        J = self.subalgebra_generated_by()
+        # The subalgebra is transient (we return an element of the
+        # superalgebra, i.e. this algebra) so why bother
+        # orthonormalizing?
+        J = self.subalgebra_generated_by(orthonormalize=False)
         u = J(self)
 
         # The image of the matrix of left-u^m-multiplication
@@ -1461,14 +1479,12 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         # subspace... or do we? Can't we just solve, knowing that
         # A(c) = u^(s+1) should have a solution in the big space,
         # too?
-        #
-        # Beware, solve_right() means that we're using COLUMN vectors.
-        # Our FiniteDimensionalAlgebraElement superclass uses rows.
         u_next = u**(s+1)
         A = u_next.operator().matrix()
         c = J.from_vector(A.solve_right(u_next.to_vector()))
 
-        # Now c is the idempotent we want, but it still lives in the subalgebra.
+        # Now c is the idempotent we want, but it still lives in
+        # the subalgebra.
         return c.superalgebra_element()
 
 
@@ -1542,6 +1558,102 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         # we want the negative of THAT for the trace.
         return -p(*self.to_vector())
 
+    def operator_inner_product(self, other):
+        r"""
+        Return the operator inner product of myself and ``other``.
+
+        The "operator inner product," whose name is not standard, is
+        defined be the usual linear-algebraic trace of the
+        ``(x*y).operator()``.
+
+        Proposition III.1.5 in Faraut and Korányi shows that on any
+        Euclidean Jordan algebra, this is another associative inner
+        product under which the cone of squares is symmetric.
+
+        This *probably* works even if the basis hasn't been
+        orthonormalized because the eigenvalues of the corresponding
+        matrix don't change when the basis does (they're preserved by
+        any similarity transformation).
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (JordanSpinEJA,
+            ....:                                  RealSymmetricEJA,
+            ....:                                  ComplexHermitianEJA,
+            ....:                                  random_eja)
+
+        EXAMPLES:
+
+        Proposition III.4.2 of Faraut and Korányi shows that on a
+        simple algebra of rank `r` and dimension `n`, this inner
+        product is `n/r` times the canonical
+        :meth:`trace_inner_product`::
+
+            sage: J = JordanSpinEJA(4, field=QQ)
+            sage: x,y = J.random_elements(2)
+            sage: n = J.dimension()
+            sage: r = J.rank()
+            sage: actual = x.operator_inner_product(y)
+            sage: expected = (n/r)*x.trace_inner_product(y)
+            sage: actual == expected
+            True
+
+        ::
+
+            sage: J = RealSymmetricEJA(3)
+            sage: x,y = J.random_elements(2)
+            sage: n = J.dimension()
+            sage: r = J.rank()
+            sage: actual = x.operator_inner_product(y)
+            sage: expected = (n/r)*x.trace_inner_product(y)
+            sage: actual == expected
+            True
+
+        ::
+
+            sage: J = ComplexHermitianEJA(3, field=QQ, orthonormalize=False)
+            sage: x,y = J.random_elements(2)
+            sage: n = J.dimension()
+            sage: r = J.rank()
+            sage: actual = x.operator_inner_product(y)
+            sage: expected = (n/r)*x.trace_inner_product(y)
+            sage: actual == expected
+            True
+
+        TESTS:
+
+        The operator inner product is commutative, bilinear, and
+        associative::
+
+            sage: J = random_eja()
+            sage: x,y,z = J.random_elements(3)
+            sage: # commutative
+            sage: x.operator_inner_product(y) == y.operator_inner_product(x)
+            True
+            sage: # bilinear
+            sage: a = J.base_ring().random_element()
+            sage: actual = (a*(x+z)).operator_inner_product(y)
+            sage: expected = ( a*x.operator_inner_product(y) +
+            ....:              a*z.operator_inner_product(y) )
+            sage: actual == expected
+            True
+            sage: actual = x.operator_inner_product(a*(y+z))
+            sage: expected = ( a*x.operator_inner_product(y) +
+            ....:              a*x.operator_inner_product(z) )
+            sage: actual == expected
+            True
+            sage: # associative
+            sage: actual = (x*y).operator_inner_product(z)
+            sage: expected = y.operator_inner_product(x*z)
+            sage: actual == expected
+            True
+
+        """
+        if not other in self.parent():
+            raise TypeError("'other' must live in the same algebra")
+
+        return (self*other).operator().matrix().trace()
+
 
     def trace_inner_product(self, other):
         """
@@ -1561,7 +1673,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             sage: x.trace_inner_product(y) == y.trace_inner_product(x)
             True
             sage: # bilinear
-            sage: a = J.base_ring().random_element();
+            sage: a = J.base_ring().random_element()
             sage: actual = (a*(x+z)).trace_inner_product(y)
             sage: expected = ( a*x.trace_inner_product(y) +
             ....:              a*z.trace_inner_product(y) )
@@ -1610,19 +1722,18 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         return self.trace_inner_product(self).sqrt()
 
 
-class CartesianProductEJAElement(FiniteDimensionalEJAElement):
-    def det(self):
-        r"""
-        Compute the determinant of this product-element using the
-        determianants of its factors.
-
-        This result Follows from the spectral decomposition of (say)
-        the pair `(x,y)` in terms of the Jordan frame `\left\{ (c_1,
-        0),(c_2, 0),...,(0,d_1),(0,d_2),... \right\}.
-        """
-        from sage.misc.misc_c import prod
-        return prod( f.det() for f in self.cartesian_factors() )
+class CartesianProductParentEJAElement(EJAElement):
+    r"""
+    An intermediate class for elements that have a Cartesian
+    product as their parent algebra.
 
+    This is needed because the ``to_matrix`` method (which gives you a
+    representation from the superalgebra) needs to do special stuff
+    for Cartesian products. Specifically, an EJA subalgebra of a
+    Cartesian product EJA will not itself be a Cartesian product (it
+    has its own basis) -- but we want ``to_matrix()`` to be able to
+    give us a Cartesian product representation.
+    """
     def to_matrix(self):
         # An override is necessary to call our custom _scale().
         B = self.parent().matrix_basis()
@@ -1630,7 +1741,20 @@ class CartesianProductEJAElement(FiniteDimensionalEJAElement):
 
         # Aaaaand linear combinations don't work in Cartesian
         # product spaces, even though they provide a method with
-        # that name. This is hidden behind an "if" because the
+        # that name. This is hidden in a subclass because the
         # _scale() function is slow.
         pairs = zip(B, self.to_vector())
         return W.sum( _scale(b, alpha) for (b,alpha) in pairs )
+
+class CartesianProductEJAElement(CartesianProductParentEJAElement):
+    def det(self):
+        r"""
+        Compute the determinant of this product-element using the
+        determianants of its factors.
+
+        This result Follows from the spectral decomposition of (say)
+        the pair `(x,y)` in terms of the Jordan frame `\left\{ (c_1,
+        0),(c_2, 0),...,(0,d_1),(0,d_2),... \right\}.
+        """
+        from sage.misc.misc_c import prod
+        return prod( f.det() for f in self.cartesian_factors() )