]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_operator.py
eja: add is_self_adjoint() for operators.
[sage.d.git] / mjo / eja / eja_operator.py
index 0e898b59fb9d8a99090da6f4baec3e49df51d9cd..0b52f555d51f58341fe56f5b63984f88cdf99da0 100644 (file)
@@ -3,12 +3,48 @@ from sage.categories.all import FreeModules
 from sage.categories.map import Map
 
 class FiniteDimensionalEuclideanJordanAlgebraOperator(Map):
+    r"""
+    An operator between two finite-dimensional Euclidean Jordan algebras.
+
+    SETUP::
+
+        sage: from mjo.eja.eja_algebra import HadamardEJA
+        sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
+
+    EXAMPLES:
+
+    The domain and codomain must be finite-dimensional Euclidean
+    Jordan algebras; if either is not, then an error is raised::
+
+        sage: J = HadamardEJA(3)
+        sage: V = VectorSpace(J.base_ring(), 3)
+        sage: M = matrix.identity(J.base_ring(), 3)
+        sage: FiniteDimensionalEuclideanJordanAlgebraOperator(V,J,M)
+        Traceback (most recent call last):
+        ...
+        TypeError: domain must be a finite-dimensional Euclidean
+        Jordan algebra
+        sage: FiniteDimensionalEuclideanJordanAlgebraOperator(J,V,M)
+        Traceback (most recent call last):
+        ...
+        TypeError: codomain must be a finite-dimensional Euclidean
+        Jordan algebra
+
+    """
+
     def __init__(self, domain_eja, codomain_eja, mat):
-        # if not (
-        #   isinstance(domain_eja, FiniteDimensionalEuclideanJordanAlgebra) and
-        #   isinstance(codomain_eja, FiniteDimensionalEuclideanJordanAlgebra) ):
-        #     raise ValueError('(co)domains must be finite-dimensional Euclidean '
-        #                      'Jordan algebras')
+        from mjo.eja.eja_algebra import FiniteDimensionalEuclideanJordanAlgebra
+
+        # I guess we should check this, because otherwise you could
+        # pass in pretty much anything algebraish.
+        if not isinstance(domain_eja,
+                          FiniteDimensionalEuclideanJordanAlgebra):
+            raise TypeError('domain must be a finite-dimensional '
+                            'Euclidean Jordan algebra')
+        if not isinstance(codomain_eja,
+                          FiniteDimensionalEuclideanJordanAlgebra):
+            raise TypeError('codomain must be a finite-dimensional '
+                            'Euclidean Jordan algebra')
 
         F = domain_eja.base_ring()
         if not (F == codomain_eja.base_ring()):
@@ -117,17 +153,17 @@ class FiniteDimensionalEuclideanJordanAlgebraOperator(Map):
             sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
             sage: from mjo.eja.eja_algebra import (
             ....:   JordanSpinEJA,
-            ....:   RealCartesianProductEJA,
+            ....:   HadamardEJA,
             ....:   RealSymmetricEJA)
 
         EXAMPLES::
 
             sage: J1 = JordanSpinEJA(3)
-            sage: J2 = RealCartesianProductEJA(2)
+            sage: J2 = HadamardEJA(2)
             sage: J3 = RealSymmetricEJA(1)
-            sage: mat1 = matrix(QQ, [[1,2,3],
+            sage: mat1 = matrix(AA, [[1,2,3],
             ....:                    [4,5,6]])
-            sage: mat2 = matrix(QQ, [[7,8]])
+            sage: mat2 = matrix(AA, [[7,8]])
             sage: g = FiniteDimensionalEuclideanJordanAlgebraOperator(J1,
             ....:                                                     J2,
             ....:                                                     mat1)
@@ -139,9 +175,9 @@ class FiniteDimensionalEuclideanJordanAlgebraOperator(Map):
             algebras represented by the matrix:
             [39 54 69]
             Domain: Euclidean Jordan algebra of dimension 3 over
-            Rational Field
+            Algebraic Real Field
             Codomain: Euclidean Jordan algebra of dimension 1 over
-            Rational Field
+            Algebraic Real Field
 
         """
         return FiniteDimensionalEuclideanJordanAlgebraOperator(
@@ -341,9 +377,9 @@ class FiniteDimensionalEuclideanJordanAlgebraOperator(Map):
             [1 0]
             [0 1]
             Domain: Euclidean Jordan algebra of dimension 2 over
-            Rational Field
+            Algebraic Real Field
             Codomain: Euclidean Jordan algebra of dimension 2 over
-            Rational Field
+            Algebraic Real Field
 
         """
         msg = ("Linear operator between finite-dimensional Euclidean Jordan "
@@ -383,6 +419,171 @@ class FiniteDimensionalEuclideanJordanAlgebraOperator(Map):
         return (self + (-other))
 
 
+    def is_self_adjoint(self):
+        r"""
+        Return whether or not this operator is self-adjoint.
+
+        At least in Sage, the fact that the base field is real means
+        that the algebra elements have to be real as well (this is why
+        we embed the complex numbers and quaternions). As a result, the
+        matrix of this operator will contain only real entries, and it
+        suffices to check only symmetry, not conjugate symmetry.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (JordanSpinEJA)
+
+        EXAMPLES::
+
+            sage: J = JordanSpinEJA(4)
+            sage: J.one().operator().is_self_adjoint()
+            True
+
+        """
+        return self.matrix().is_symmetric()
+
+
+    def is_zero(self):
+        r"""
+        Return whether or not this map is the zero operator.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator
+            sage: from mjo.eja.eja_algebra import (random_eja,
+            ....:                                  JordanSpinEJA,
+            ....:                                  RealSymmetricEJA)
+
+        EXAMPLES::
+
+            sage: J1 = JordanSpinEJA(2)
+            sage: J2 = RealSymmetricEJA(2)
+            sage: R = J1.base_ring()
+            sage: M = matrix(R, [ [0, 0],
+            ....:                 [0, 0],
+            ....:                 [0, 0] ])
+            sage: L = FiniteDimensionalEuclideanJordanAlgebraOperator(J1,J2,M)
+            sage: L.is_zero()
+            True
+            sage: M = matrix(R, [ [0, 0],
+            ....:                 [0, 1],
+            ....:                 [0, 0] ])
+            sage: L = FiniteDimensionalEuclideanJordanAlgebraOperator(J1,J2,M)
+            sage: L.is_zero()
+            False
+
+        TESTS:
+
+        The left-multiplication-by-zero operation on a given algebra
+        is its zero map::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: J.zero().operator().is_zero()
+            True
+
+        """
+        return self.matrix().is_zero()
+
+
+    def inverse(self):
+        """
+        Return the inverse of this operator, if it exists.
+
+        The reason this method is not simply an alias for the built-in
+        :meth:`__invert__` is that the built-in inversion is a bit magic
+        since it's intended to be a unary operator. If we alias ``inverse``
+        to ``__invert__``, then we wind up having to call e.g. ``A.inverse``
+        without parentheses.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import RealSymmetricEJA, random_eja
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: x = sum(J.gens())
+            sage: x.operator().inverse().matrix()
+            [3/2  -1 1/2]
+            [ -1   2  -1]
+            [1/2  -1 3/2]
+            sage: x.operator().matrix().inverse()
+            [3/2  -1 1/2]
+            [ -1   2  -1]
+            [1/2  -1 3/2]
+
+        TESTS:
+
+        The identity operator is its own inverse::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: idJ = J.one().operator()
+            sage: idJ.inverse() == idJ
+            True
+
+        The inverse of the inverse is the operator we started with::
+
+            sage: set_random_seed()
+            sage: x = random_eja().random_element()
+            sage: L = x.operator()
+            sage: not L.is_invertible() or (L.inverse().inverse() == L)
+            True
+
+        """
+        return ~self
+
+
+    def is_invertible(self):
+        """
+        Return whether or not this operator is invertible.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (RealSymmetricEJA,
+            ....:                                  TrivialEJA,
+            ....:                                  random_eja)
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: x = sum(J.gens())
+            sage: x.operator().matrix()
+            [  1 1/2   0]
+            [1/2   1 1/2]
+            [  0 1/2   1]
+            sage: x.operator().matrix().is_invertible()
+            True
+            sage: x.operator().is_invertible()
+            True
+
+        The zero operator is invertible in a trivial algebra::
+
+            sage: J = TrivialEJA()
+            sage: J.zero().operator().is_invertible()
+            True
+
+        TESTS:
+
+        The identity operator is always invertible::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: J.one().operator().is_invertible()
+            True
+
+        The zero operator is never invertible in a nontrivial algebra::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: not J.is_trivial() and J.zero().operator().is_invertible()
+            False
+
+        """
+        return self.matrix().is_invertible()
+
+
     def matrix(self):
         """
         Return the matrix representation of this operator with respect
@@ -426,3 +627,62 @@ class FiniteDimensionalEuclideanJordanAlgebraOperator(Map):
         """
         # The matrix method returns a polynomial in 'x' but want one in 't'.
         return self.matrix().minimal_polynomial().change_variable_name('t')
+
+
+    def spectral_decomposition(self):
+        """
+        Return the spectral decomposition of this operator as a list of
+        (eigenvalue, orthogonal projector) pairs.
+
+        This is the unique spectral decomposition, up to the order of
+        the projection operators, with distinct eigenvalues. So, the
+        projections are generally onto subspaces of dimension greater
+        than one.
+
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import RealSymmetricEJA
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(4)
+            sage: x = sum(J.gens())
+            sage: A = x.subalgebra_generated_by(orthonormalize_basis=True)
+            sage: L0x = A(x).operator()
+            sage: sd = L0x.spectral_decomposition()
+            sage: l0 = sd[0][0]
+            sage: l1 = sd[1][0]
+            sage: P0 = sd[0][1]
+            sage: P1 = sd[1][1]
+            sage: P0*l0 + P1*l1 == L0x
+            True
+            sage: P0 + P1 == P0^0 # the identity
+            True
+            sage: P0^2 == P0
+            True
+            sage: P1^2 == P1
+            True
+            sage: P0*P1 == A.zero().operator()
+            True
+            sage: P1*P0 == A.zero().operator()
+            True
+
+        """
+        if not self.matrix().is_symmetric():
+            raise ValueError('algebra basis is not orthonormal')
+
+        D,P = self.matrix().jordan_form(subdivide=False,transformation=True)
+        eigenvalues = D.diagonal()
+        us = P.columns()
+        projectors = []
+        for i in range(len(us)):
+            # they won't be normalized, but they have to be
+            # for the spectral theorem to work.
+            us[i] = us[i]/us[i].norm()
+            mat = us[i].column()*us[i].row()
+            Pi = FiniteDimensionalEuclideanJordanAlgebraOperator(
+                   self.domain(),
+                   self.codomain(),
+                   mat)
+            projectors.append(Pi)
+        return list(zip(eigenvalues, projectors))