From e78216245caff4f11de66433f28fd995f4670b78 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Fri, 20 Nov 2020 21:32:57 -0500 Subject: [PATCH] eja: echelonize a matrix manually to ensure a decent algorithm is used. --- mjo/eja/eja_element_subalgebra.py | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/mjo/eja/eja_element_subalgebra.py b/mjo/eja/eja_element_subalgebra.py index a26381b..d2b54fe 100644 --- a/mjo/eja/eja_element_subalgebra.py +++ b/mjo/eja/eja_element_subalgebra.py @@ -1,4 +1,5 @@ from sage.matrix.constructor import matrix +from sage.rings.all import QQ from mjo.eja.eja_subalgebra import FiniteDimensionalEuclideanJordanSubalgebra @@ -18,6 +19,19 @@ class FiniteDimensionalEuclideanJordanElementSubalgebra(FiniteDimensionalEuclide P = matrix(field, power_vectors) if orthonormalize_basis == False: + # Echelonize the matrix ourselves, because otherwise the + # call to P.pivot_rows() below can choose a non-optimal + # row-reduction algorithm. In particular, scaling can + # help over AA because it avoids the RecursionError that + # gets thrown when we have to look too hard for a root. + # + # Beware: QQ supports an entirely different set of "algorithm" + # keywords than do AA and RR. + algo = None + if field is not QQ: + algo = "scaled_partial_pivoting" + P.echelonize(algorithm=algo) + # In this case, we just need to figure out which elements # of the "powers" list are redundant... First compute the # vector subspace spanned by the powers of the given -- 2.44.2