]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/eja/eja_element.py
eja: add subalgebra() method.
[sage.d.git] / mjo / eja / eja_element.py
index d3e9a33ceba6e1fc4e58b417f3b953e2e2d1c3d7..42e5782c4539c97fcd71805b50c81cd4a96877a5 100644 (file)
@@ -1,4 +1,5 @@
 from sage.matrix.constructor import matrix
+from sage.misc.cachefunc import cached_method
 from sage.modules.free_module import VectorSpace
 from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement
 
@@ -438,6 +439,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         return ((-1)**r)*p(*self.to_vector())
 
 
+    @cached_method
     def inverse(self):
         """
         Return the Jordan-multiplicative inverse of this element.
@@ -482,7 +484,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             sage: JordanSpinEJA(3).zero().inverse()
             Traceback (most recent call last):
             ...
-            ValueError: element is not invertible
+            ZeroDivisionError: element is not invertible
 
         TESTS:
 
@@ -534,40 +536,38 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             sage: slow == fast                                   # long time
             True
         """
-        if not self.is_invertible():
-            raise ValueError("element is not invertible")
-
+        not_invertible_msg = "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!
+            if self.det().is_zero():
+                raise ZeroDivisionError(not_invertible_msg)
             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)
+        try:
+            inv = (~self.quadratic_representation())(self)
+            self.is_invertible.set_cache(True)
+            return inv
+        except ZeroDivisionError:
+            self.is_invertible.set_cache(False)
+            raise ZeroDivisionError(not_invertible_msg)
 
 
+    @cached_method
     def is_invertible(self):
         """
         Return whether or not this element is invertible.
 
         ALGORITHM:
 
-        The usual way to do this is to check if the determinant is
-        zero, but we need the characteristic polynomial for the
-        determinant. The minimal polynomial is a lot easier to get,
-        so we use Corollary 2 in Chapter V of Koecher to check
-        whether or not the parent 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.
+        If computing my determinant will be fast, we do so and compare
+        with zero (Proposition II.2.4 in Faraut and
+        Koranyi). Otherwise, Proposition II.3.2 in Faraut and Koranyi
+        reduces the problem to the invertibility of my quadratic
+        representation.
 
         SETUP::
 
@@ -600,7 +600,6 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             sage: fast = x.is_invertible()                       # long time
             sage: slow == fast                                   # long time
             True
-
         """
         if self.is_zero():
             if self.parent().is_trivial():
@@ -609,15 +608,17 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
                 return False
 
         if self.parent()._charpoly_coefficients.is_in_cache():
-            # The determinant will be quicker than computing the minimal
-            # polynomial from scratch, most likely.
+            # The determinant will be quicker than inverting the
+            # quadratic representation, 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()
-        p = self.minimal_polynomial()
-        return not (p(zero) == zero)
+        # The easiest way to determine if I'm invertible is to try.
+        try:
+            inv = (~self.quadratic_representation())(self)
+            self.inverse.set_cache(inv)
+            return True
+        except ZeroDivisionError:
+            return False
 
 
     def is_primitive_idempotent(self):
@@ -1116,6 +1117,7 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
         return W.linear_combination( zip(B, self.to_vector()) )
 
 
+
     def norm(self):
         """
         The norm of this element with respect to :meth:`inner_product`.
@@ -1411,8 +1413,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             True
 
         """
-        from mjo.eja.eja_element_subalgebra import FiniteDimensionalEJAElementSubalgebra
-        return FiniteDimensionalEJAElementSubalgebra(self, **kwargs)
+        powers = tuple( self**k for k in range(self.degree()) )
+        A = self.parent().subalgebra(powers, associative=True, **kwargs)
+        A.one.set_cache(A(self.parent().one()))
+        return A
 
 
     def subalgebra_idempotent(self):
@@ -1519,6 +1523,15 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
             sage: J.random_element().trace() in RLF
             True
 
+        The trace is linear::
+
+            sage: set_random_seed()
+            sage: J = random_eja()
+            sage: x,y = J.random_elements(2)
+            sage: alpha = J.base_ring().random_element()
+            sage: (alpha*x + y).trace() == alpha*x.trace() + y.trace()
+            True
+
         """
         P = self.parent()
         r = P.rank()
@@ -1602,3 +1615,38 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement):
 
         """
         return self.trace_inner_product(self).sqrt()
+
+
+
+class CartesianProductEJAElement(FiniteDimensionalEJAElement):
+
+    def to_matrix(self):
+        r"""
+        SETUP::
+
+            sage: from mjo.eja.eja_algebra import (HadamardEJA,
+            ....:                                  RealSymmetricEJA)
+
+        EXAMPLES::
+
+            sage: J1 = HadamardEJA(1)
+            sage: J2 = RealSymmetricEJA(2)
+            sage: J = cartesian_product([J1,J2])
+            sage: x = sum(J.gens())
+            sage: x.to_matrix()[0]
+            [1]
+            sage: x.to_matrix()[1]
+            [                  1 0.7071067811865475?]
+            [0.7071067811865475?                   1]
+
+        """
+        B = self.parent().matrix_basis()
+        W = self.parent().matrix_space()
+
+        # Aaaaand linear combinations don't work in Cartesian
+        # product spaces, even though they provide a method
+        # with that name.
+        pairs = zip(B, self.to_vector())
+        return sum( ( W(tuple(alpha*b_i for b_i in b))
+                      for (b,alpha) in pairs ),
+                    W.zero())