]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/euclidean_jordan_algebra.py
eja: fix an erroneous test case.
[sage.d.git] / mjo / eja / euclidean_jordan_algebra.py
index 920ecc988d21c8ce698fdd4dcf0aefceffdb8236..30d04f4a352b7d553470bd456d18d224e2ec5428 100644 (file)
@@ -5,12 +5,405 @@ are used in optimization, and have some additional nice methods beyond
 what can be supported in a general Jordan Algebra.
 """
 
-from sage.categories.magmatic_algebras import MagmaticAlgebras
+from sage.categories.finite_dimensional_algebras_with_basis import FiniteDimensionalAlgebrasWithBasis
+from sage.categories.morphism import SetMorphism
 from sage.structure.element import is_Matrix
 from sage.structure.category_object import normalize_names
 
 from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra import FiniteDimensionalAlgebra
 from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra_element import FiniteDimensionalAlgebraElement
+from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra_morphism import FiniteDimensionalAlgebraMorphism, FiniteDimensionalAlgebraHomset
+
+
+class FiniteDimensionalEuclideanJordanAlgebraHomset(FiniteDimensionalAlgebraHomset):
+
+    def has_coerce_map_from(self, S):
+        """
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: H = J.Hom(J)
+            sage: H.has_coerce_map_from(QQ)
+            True
+
+        """
+        try:
+            # The Homset classes override has_coerce_map_from() with
+            # something that crashes when it's given e.g. QQ.
+            if S.is_subring(self.codomain().base_ring()):
+                return True
+        except:
+            pclass = super(FiniteDimensionalEuclideanJordanAlgebraHomset, self)
+            return pclass.has_coerce_map_from(S)
+
+
+    def _coerce_map_from_(self, S):
+        """
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: H = J.Hom(J)
+            sage: H.coerce(2)
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [2 0 0]
+            [0 2 0]
+            [0 0 2]
+
+        """
+        C = self.codomain()
+        R = C.base_ring()
+        if S.is_subring(R):
+            h = S.hom(self.codomain())
+            return SetMorphism(Hom(S,C), lambda x: h(x).operator())
+
+
+    def __call__(self, x):
+        """
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: H = J.Hom(J)
+            sage: H(2)
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [2 0 0]
+            [0 2 0]
+            [0 0 2]
+
+        """
+        if x in self.base_ring():
+            cols = self.domain().dimension()
+            rows = self.codomain().dimension()
+            x = x*identity_matrix(self.codomain().base_ring(), rows, cols)
+        return FiniteDimensionalEuclideanJordanAlgebraMorphism(self, x)
+
+
+    def one(self):
+        """
+        Return the identity morphism, but as a member of the right
+        space (so that we can add it, multiply it, etc.)
+        """
+        cols = self.domain().dimension()
+        rows = self.codomain().dimension()
+        mat = identity_matrix(self.base_ring(), rows, cols)
+        return FiniteDimensionalEuclideanJordanAlgebraMorphism(self, mat)
+
+
+
+class FiniteDimensionalEuclideanJordanAlgebraMorphism(FiniteDimensionalAlgebraMorphism):
+    """
+    A linear map between two finite-dimensional EJAs.
+
+    This is a very thin wrapper around FiniteDimensionalAlgebraMorphism
+    that does only a few things:
+
+      1. Avoids the ``unitary`` and ``check`` arguments to the constructor
+         that will always be ``False``. This is necessary because these
+         are homomorphisms with respect to ADDITION, but the SageMath
+         machinery wants to check that they're homomorphisms with respect
+         to (Jordan) MULTIPLICATION. That obviously doesn't work.
+
+      2. Inputs and outputs the underlying matrix with respect to COLUMN
+         vectors, unlike the parent class.
+
+      3. Allows us to add, subtract, negate, multiply (compose), and
+         invert morphisms in the obvious way.
+
+    If this seems a bit heavyweight, it is. I would have been happy to
+    use a the ring morphism that underlies the finite-dimensional
+    algebra morphism, but they don't seem to be callable on elements of
+    our EJA, and you can't add/multiply/etc. them.
+    """
+    def _add_(self, other):
+        """
+        Add two EJA morphisms in the obvious way.
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(3)
+            sage: x = J.zero()
+            sage: y = J.one()
+            sage: x.operator() + y.operator()
+            Morphism from Euclidean Jordan algebra of degree 6 over Rational
+            Field to Euclidean Jordan algebra of degree 6 over Rational Field
+            given by matrix
+            [1 0 0 0 0 0]
+            [0 1 0 0 0 0]
+            [0 0 1 0 0 0]
+            [0 0 0 1 0 0]
+            [0 0 0 0 1 0]
+            [0 0 0 0 0 1]
+
+        TESTS::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: y = J.random_element()
+            sage: (x.operator() + y.operator()) in J.Hom(J)
+            True
+
+        """
+        P = self.parent()
+        if not other in P:
+            raise ValueError("summands must live in the same space")
+
+        return FiniteDimensionalEuclideanJordanAlgebraMorphism(
+                  P,
+                  self.matrix() + other.matrix() )
+
+
+    def __init__(self, parent, f):
+        FiniteDimensionalAlgebraMorphism.__init__(self,
+                                                  parent,
+                                                  f.transpose(),
+                                                  unitary=False,
+                                                  check=False)
+
+
+    def __invert__(self):
+        """
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: x = J.linear_combination(zip(range(len(J.gens())), J.gens()))
+            sage: x.is_invertible()
+            True
+            sage: ~x.operator()
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [-3/2    2 -1/2]
+            [   1    0    0]
+            [-1/2    0  1/2]
+            sage: x.operator_matrix().inverse()
+            [-3/2    2 -1/2]
+            [   1    0    0]
+            [-1/2    0  1/2]
+
+        TESTS::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: not x.is_invertible() or (
+            ....:   (~x.operator()).matrix() == x.operator_matrix().inverse() )
+            True
+
+        """
+        A = self.matrix()
+        if not A.is_invertible():
+            raise ValueError("morphism is not invertible")
+
+        P = self.parent()
+        return FiniteDimensionalEuclideanJordanAlgebraMorphism(self.parent(),
+                                                                A.inverse())
+
+    def _lmul_(self, right):
+        """
+        Compose two EJA morphisms using multiplicative notation.
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: x = J.zero()
+            sage: y = J.one()
+            sage: x.operator() * y.operator()
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [0 0 0]
+            [0 0 0]
+            [0 0 0]
+
+        ::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: x = J.linear_combination(zip(range(len(J.gens())), J.gens()))
+            sage: x.operator()
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [  0   1   0]
+            [1/2   1 1/2]
+            [  0   1   2]
+            sage: 2*x.operator()
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [0 2 0]
+            [1 2 1]
+            [0 2 4]
+            sage: x.operator()*2
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [0 2 0]
+            [1 2 1]
+            [0 2 4]
+
+        TESTS::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: y = J.random_element()
+            sage: (x.operator() * y.operator()) in J.Hom(J)
+            True
+
+        """
+        try:
+            # I think the morphism classes break the coercion framework
+            # somewhere along the way, so we have to do this ourselves.
+            right = self.parent().coerce(right)
+        except:
+            pass
+
+        if not right.codomain() is self.domain():
+            raise ValueError("(co)domains must agree for composition")
+
+        return FiniteDimensionalEuclideanJordanAlgebraMorphism(
+                 self.parent(),
+                 self.matrix()*right.matrix() )
+
+    __mul__ = _lmul_
+
+
+    def __pow__(self, n):
+        """
+
+        TESTS::
+
+            sage: J = JordanSpinEJA(4)
+            sage: e0,e1,e2,e3 = J.gens()
+            sage: x = -5/2*e0 + 1/2*e2 + 20*e3
+            sage: Qx = x.quadratic_representation()
+            sage: Qx^0
+            Morphism from Euclidean Jordan algebra of degree 4 over Rational
+            Field to Euclidean Jordan algebra of degree 4 over Rational Field
+            given by matrix
+            [1 0 0 0]
+            [0 1 0 0]
+            [0 0 1 0]
+            [0 0 0 1]
+            sage: (x^0).quadratic_representation() == Qx^0
+            True
+
+        """
+        if n == 0:
+            # We get back the stupid identity morphism which doesn't
+            # live in the right space.
+            return self.parent().one()
+        elif n == 1:
+            return self
+        else:
+            return FiniteDimensionalAlgebraMorphism.__pow__(self,n)
+
+
+    def _neg_(self):
+        """
+        Negate this morphism.
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: x = J.one()
+            sage: -x.operator()
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational Field
+            given by matrix
+            [-1  0  0]
+            [ 0 -1  0]
+            [ 0  0 -1]
+
+        TESTS::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: -x.operator() in J.Hom(J)
+            True
+
+        """
+        return FiniteDimensionalEuclideanJordanAlgebraMorphism(
+                  self.parent(),
+                  -self.matrix() )
+
+
+    def _repr_(self):
+        """
+        We override only the representation that is shown to the user,
+        because we want the matrix to be with respect to COLUMN vectors.
+
+        TESTS:
+
+        Ensure that we see the transpose of the underlying matrix object:
+
+            sage: J = RealSymmetricEJA(3)
+            sage: x = J.linear_combination(zip(range(len(J.gens())), J.gens()))
+            sage: L = x.operator()
+            sage: L
+            Morphism from Euclidean Jordan algebra of degree 6 over Rational
+            Field to Euclidean Jordan algebra of degree 6 over Rational Field
+            given by matrix
+            [  0   1   2   0   0   0]
+            [1/2 3/2   2 1/2   1   0]
+            [  1   2 5/2   0 1/2   1]
+            [  0   1   0   3   4   0]
+            [  0   1 1/2   2   4   2]
+            [  0   0   2   0   4   5]
+            sage: L._matrix
+            [  0 1/2   1   0   0   0]
+            [  1 3/2   2   1   1   0]
+            [  2   2 5/2   0 1/2   2]
+            [  0 1/2   0   3   2   0]
+            [  0   1 1/2   4   4   4]
+            [  0   0   1   0   2   5]
+
+        """
+        return "Morphism from {} to {} given by matrix\n{}".format(
+            self.domain(), self.codomain(), self.matrix())
+
+
+    def __sub__(self, other):
+        """
+        Subtract one morphism from another using addition and negation.
+
+        EXAMPLES::
+
+            sage: J = RealSymmetricEJA(2)
+            sage: L1 = J.one().operator()
+            sage: L1 - L1
+            Morphism from Euclidean Jordan algebra of degree 3 over Rational
+            Field to Euclidean Jordan algebra of degree 3 over Rational
+            Field given by matrix
+            [0 0 0]
+            [0 0 0]
+            [0 0 0]
+
+        TESTS::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x = J.random_element()
+            sage: y = J.random_element()
+            sage: x.operator() - y.operator() in J.Hom(J)
+            True
+
+        """
+        return self + (-other)
+
+
+    def matrix(self):
+        """
+        Return the matrix of this morphism with respect to a left-action
+        on column vectors.
+        """
+        return FiniteDimensionalAlgebraMorphism.matrix(self).transpose()
+
 
 class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
     @staticmethod
@@ -30,7 +423,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 raise ValueError("input is not a multiplication table")
         mult_table = tuple(mult_table)
 
-        cat = MagmaticAlgebras(field).FiniteDimensional().WithBasis()
+        cat = FiniteDimensionalAlgebrasWithBasis(field)
         cat.or_subcategory(category)
         if assume_associative:
             cat = cat.Associative()
@@ -48,6 +441,15 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                                  natural_basis=natural_basis)
 
 
+    def _Hom_(self, B, cat):
+        """
+        Construct a homset of ``self`` and ``B``.
+        """
+        return FiniteDimensionalEuclideanJordanAlgebraHomset(self,
+                                                             B,
+                                                             category=cat)
+
+
     def __init__(self,
                  field,
                  mult_table,
@@ -340,6 +742,16 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
         An element of a Euclidean Jordan algebra.
         """
 
+        def __dir__(self):
+            """
+            Oh man, I should not be doing this. This hides the "disabled"
+            methods ``left_matrix`` and ``matrix`` from introspection;
+            in particular it removes them from tab-completion.
+            """
+            return filter(lambda s: s not in ['left_matrix', 'matrix'],
+                          dir(self.__class__) )
+
+
         def __init__(self, A, elt=None):
             """
             EXAMPLES:
@@ -658,8 +1070,10 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
             """
             Return the Jordan-multiplicative inverse of this element.
 
-            We can't use the superclass method because it relies on the
-            algebra being associative.
+            ALGORITHM:
+
+            We appeal to the quadratic representation as in Koecher's
+            Theorem 12 in Chapter III, Section 5.
 
             EXAMPLES:
 
@@ -670,12 +1084,12 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 sage: n = ZZ.random_element(1,10)
                 sage: J = JordanSpinEJA(n)
                 sage: x = J.random_element()
-                sage: while x.is_zero():
+                sage: while not x.is_invertible():
                 ....:     x = J.random_element()
                 sage: x_vec = x.vector()
                 sage: x0 = x_vec[0]
                 sage: x_bar = x_vec[1:]
-                sage: coeff = 1/(x0^2 - x_bar.inner_product(x_bar))
+                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: x.inverse() == J(x_inverse)
@@ -698,34 +1112,27 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 sage: (not x.is_invertible()) or (x.inverse()*x == J.one())
                 True
 
-            """
-            if not self.is_invertible():
-                raise ValueError("element not invertible")
+            The inverse of the inverse is what we started with::
 
-            if self.parent().is_associative():
-                elt = FiniteDimensionalAlgebraElement(self.parent(), self)
-                # elt is in the right coordinates, but has the wrong class.
-                return self.parent()(elt.inverse().vector())
+                sage: set_random_seed()
+                sage: J = random_eja()
+                sage: x = J.random_element()
+                sage: (not x.is_invertible()) or (x.inverse().inverse() == x)
+                True
 
-            # We do this a little different than the usual recursive
-            # call to a finite-dimensional algebra element, because we
-            # wind up with an inverse that lives in the subalgebra and
-            # we need information about the parent to convert it back.
-            V = self.span_of_powers()
-            assoc_subalg = self.subalgebra_generated_by()
-            # Mis-design warning: the basis used for span_of_powers()
-            # and subalgebra_generated_by() must be the same, and in
-            # the same order!
-            elt = assoc_subalg(V.coordinates(self.vector()))
+            The zero element is never invertible::
 
-            # This will be in the subalgebra's coordinates...
-            fda_elt = FiniteDimensionalAlgebraElement(assoc_subalg, elt)
-            subalg_inverse = fda_elt.inverse()
+                sage: set_random_seed()
+                sage: J = random_eja().zero().inverse()
+                Traceback (most recent call last):
+                ...
+                ValueError: element is not invertible
+
+            """
+            if not self.is_invertible():
+                raise ValueError("element is not invertible")
 
-            # So we have to convert back...
-            basis = [ self.parent(v) for v in V.basis() ]
-            pairs = zip(subalg_inverse.vector(), basis)
-            return self.parent().linear_combination(pairs)
+            return (~self.quadratic_representation())(self)
 
 
         def is_invertible(self):
@@ -866,6 +1273,16 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
             return self.span_of_powers().dimension()
 
 
+        def left_matrix(self):
+            """
+            Our parent class defines ``left_matrix`` and ``matrix``
+            methods whose names are misleading. We don't want them.
+            """
+            raise NotImplementedError("use operator_matrix() instead")
+
+        matrix = left_matrix
+
+
         def minimal_polynomial(self):
             """
             Return the minimal polynomial of this element,
@@ -987,14 +1404,37 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
             return W.linear_combination(zip(self.vector(), B))
 
 
+        def operator(self):
+            """
+            Return the left-multiplication-by-this-element
+            operator on the ambient algebra.
+
+            TESTS::
+
+                sage: set_random_seed()
+                sage: J = random_eja()
+                sage: x = J.random_element()
+                sage: y = J.random_element()
+                sage: x.operator()(y) == x*y
+                True
+                sage: y.operator()(x) == x*y
+                True
+
+            """
+            P = self.parent()
+            return FiniteDimensionalEuclideanJordanAlgebraMorphism(
+                     Hom(P,P),
+                     self.operator_matrix() )
+
+
+
         def operator_matrix(self):
             """
             Return the matrix that represents left- (or right-)
             multiplication by this element in the parent algebra.
 
-            We have to override this because the superclass method
-            returns a matrix that acts on row vectors (that is, on
-            the right).
+            We implement this ourselves to work around the fact that
+            our parent class represents everything with row vectors.
 
             EXAMPLES:
 
@@ -1077,7 +1517,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 sage: D = (x0^2 - x_bar.inner_product(x_bar))*D
                 sage: D = D + 2*x_bar.tensor_product(x_bar)
                 sage: Q = block_matrix(2,2,[A,B,C,D])
-                sage: Q == x.quadratic_representation()
+                sage: Q == x.quadratic_representation().matrix()
                 True
 
             Test all of the properties from Theorem 11.2 in Alizadeh::
@@ -1086,38 +1526,76 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
                 sage: J = random_eja()
                 sage: x = J.random_element()
                 sage: y = J.random_element()
+                sage: Lx = x.operator()
+                sage: Lxx = (x*x).operator()
+                sage: Qx = x.quadratic_representation()
+                sage: Qy = y.quadratic_representation()
+                sage: Qxy = x.quadratic_representation(y)
+                sage: Qex = J.one().quadratic_representation(x)
+                sage: n = ZZ.random_element(10)
+                sage: Qxn = (x^n).quadratic_representation()
 
             Property 1:
 
-                sage: actual = x.quadratic_representation(y)
-                sage: expected = ( (x+y).quadratic_representation()
-                ....:              -x.quadratic_representation()
-                ....:              -y.quadratic_representation() ) / 2
-                sage: actual == expected
+                sage: 2*Qxy == (x+y).quadratic_representation() - Qx - Qy
                 True
 
             Property 2:
 
                 sage: alpha = QQ.random_element()
-                sage: actual = (alpha*x).quadratic_representation()
-                sage: expected = (alpha^2)*x.quadratic_representation()
-                sage: actual == expected
+                sage: (alpha*x).quadratic_representation() == (alpha^2)*Qx
+                True
+
+            Property 3:
+
+                sage: not x.is_invertible() or ( Qx(x.inverse()) == x )
+                True
+
+                sage: not x.is_invertible() or (
+                ....:   ~Qx
+                ....:   ==
+                ....:   x.inverse().quadratic_representation() )
+                True
+
+                sage: Qxy(J.one()) == x*y
+                True
+
+            Property 4:
+
+                sage: not x.is_invertible() or (
+                ....:   x.quadratic_representation(x.inverse())*Qx
+                ....:   == Qx*x.quadratic_representation(x.inverse()) )
+                True
+
+                sage: not x.is_invertible() or (
+                ....:   x.quadratic_representation(x.inverse())*Qx
+                ....:   ==
+                ....:   2*x.operator()*Qex - Qx )
+                True
+
+                sage: 2*x.operator()*Qex - Qx == Lxx
                 True
 
             Property 5:
 
-                sage: Qy = y.quadratic_representation()
-                sage: actual = J(Qy*x.vector()).quadratic_representation()
-                sage: expected = Qy*x.quadratic_representation()*Qy
-                sage: actual == expected
+                sage: Qy(x).quadratic_representation() == Qy*Qx*Qy
                 True
 
             Property 6:
 
-                sage: k = ZZ.random_element(1,10)
-                sage: actual = (x^k).quadratic_representation()
-                sage: expected = (x.quadratic_representation())^k
-                sage: actual == expected
+                sage: Qxn == (Qx)^n
+                True
+
+            Property 7:
+
+                sage: not x.is_invertible() or (
+                ....:   Qx*x.inverse().operator() == Lx )
+                True
+
+            Property 8:
+
+                sage: not x.operator_commutes_with(y) or (
+                ....:   Qx(y)^n == Qxn(y^n) )
                 True
 
             """
@@ -1126,9 +1604,9 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra):
             elif not other in self.parent():
                 raise TypeError("'other' must live in the same algebra")
 
-            L = self.operator_matrix()
-            M = other.operator_matrix()
-            return ( L*M + M*L - (self*other).operator_matrix() )
+            L = self.operator()
+            M = other.operator()
+            return ( L*M + M*L - (self*other).operator() )
 
 
         def span_of_powers(self):