From a518a81a52fa629c69ab67e2c299f063ada75f00 Mon Sep 17 00:00:00 2001
From: Michael Orlitzky <michael@orlitzky.com>
Date: Thu, 24 Dec 2020 18:36:14 -0500
Subject: [PATCH] eja: cache inverse() and is_invertible() more intelligently.

---
 mjo/eja/eja_element.py | 55 +++++++++++++++++++++---------------------
 1 file changed, 28 insertions(+), 27 deletions(-)

diff --git a/mjo/eja/eja_element.py b/mjo/eja/eja_element.py
index e30dbb1..66138b2 100644
--- a/mjo/eja/eja_element.py
+++ b/mjo/eja/eja_element.py
@@ -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):
-- 
2.43.2