From f0cabe7c6e37781e4f92c9ba0e0c7413a5f6b939 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Sat, 13 Mar 2021 20:46:54 -0500 Subject: [PATCH] eja: speed up det() for Cartesian product elements. --- mjo/eja/TODO | 6 +----- mjo/eja/eja_algebra.py | 4 +++- mjo/eja/eja_element.py | 43 ++++++++++++++++++++++++++++++------------ 3 files changed, 35 insertions(+), 18 deletions(-) diff --git a/mjo/eja/TODO b/mjo/eja/TODO index 82c49d5..9cee421 100644 --- a/mjo/eja/TODO +++ b/mjo/eja/TODO @@ -7,9 +7,5 @@ the nontrivial factor. On the other hand, it's nice that we can test out some alternate code paths... -4. Conjecture: if x = (x1,x2), then det(x) = det(x1)det(x2). This - should be used to fix the fact that det(x) is monstrously slow in - Cartesian product algebras, and thus randomly in the doctests. - -5. Add dimension bounds on any tests over AA that compute element +4. Add dimension bounds on any tests over AA that compute element subalgebras. diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index 907f40d..2bad32c 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -166,7 +166,8 @@ from sage.modules.free_module import FreeModule, VectorSpace from sage.rings.all import (ZZ, QQ, AA, QQbar, RR, RLF, CLF, PolynomialRing, QuadraticField) -from mjo.eja.eja_element import FiniteDimensionalEJAElement +from mjo.eja.eja_element import (CartesianProductEJAElement, + FiniteDimensionalEJAElement) from mjo.eja.eja_operator import FiniteDimensionalEJAOperator from mjo.eja.eja_utils import _all2list @@ -3089,6 +3090,7 @@ class CartesianProductEJA(FiniteDimensionalEJA): sage: actual == expected # long time True """ + Element = CartesianProductEJAElement def __init__(self, factors, **kwargs): m = len(factors) if m == 0: diff --git a/mjo/eja/eja_element.py b/mjo/eja/eja_element.py index 693e8e2..47f6ff0 100644 --- a/mjo/eja/eja_element.py +++ b/mjo/eja/eja_element.py @@ -6,6 +6,7 @@ from sage.modules.with_basis.indexed_element import IndexedFreeModuleElement from mjo.eja.eja_operator import FiniteDimensionalEJAOperator from mjo.eja.eja_utils import _scale + class FiniteDimensionalEJAElement(IndexedFreeModuleElement): """ An element of a Euclidean Jordan algebra. @@ -1121,18 +1122,10 @@ class FiniteDimensionalEJAElement(IndexedFreeModuleElement): B = self.parent().matrix_basis() W = self.parent().matrix_space() - if hasattr(W, 'cartesian_factors'): - # 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 - # _scale() function is slow. - pairs = zip(B, self.to_vector()) - return W.sum( _scale(b, alpha) for (b,alpha) in pairs ) - else: - # This is just a manual "from_vector()", but of course - # matrix spaces aren't vector spaces in sage, so they - # don't have a from_vector() method. - return W.linear_combination( zip(B, self.to_vector()) ) + # This is just a manual "from_vector()", but of course + # matrix spaces aren't vector spaces in sage, so they + # don't have a from_vector() method. + return W.linear_combination( zip(B, self.to_vector()) ) @@ -1650,3 +1643,29 @@ 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() ) + + def to_matrix(self): + # An override is necessary to call our custom _scale(). + 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. This is hidden behind an "if" because the + # _scale() function is slow. + pairs = zip(B, self.to_vector()) + return W.sum( _scale(b, alpha) for (b,alpha) in pairs ) -- 2.43.2