]> gitweb.michael.orlitzky.com - sage.d.git/commitdiff
eja: get rid of the old rational basis constructor.
authorMichael Orlitzky <michael@orlitzky.com>
Mon, 30 Nov 2020 15:07:30 +0000 (10:07 -0500)
committerMichael Orlitzky <michael@orlitzky.com>
Mon, 30 Nov 2020 15:07:30 +0000 (10:07 -0500)
mjo/eja/eja_algebra.py

index 14666c2f2a4d47b2f1a5cde1448dcb3b778aa72a..b83f4a5c322e13b93303378a2aeab5377ec307fe 100644 (file)
@@ -214,6 +214,8 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule):
                 elt = self.from_vector(multiplication_table[i][j])
                 self._multiplication_table[i][j] = elt
 
+        self._multiplication_table = tuple(map(tuple, self._multiplication_table))
+
         # Save our inner product as a matrix, since the efficiency of
         # matrix multiplication will usually outweigh the fact that we
         # have to store a redundant upper- or lower-triangular part.
@@ -221,7 +223,8 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule):
         # in fact) in case some e.g. matrix multiplication routine can
         # take advantage of it.
         self._inner_product_matrix = matrix(field, inner_product_table)
-        self._inner_product_matrix._cache = {'hermitian': False}
+        self._inner_product_matrix._cache = {'hermitian': True}
+        self._inner_product_matrix.set_immutable()
 
         if check_axioms:
             if not self._is_jordanian():
@@ -1176,6 +1179,50 @@ class RationalBasisEuclideanJordanAlgebraNg(FiniteDimensionalEuclideanJordanAlge
         self._deortho_multiplication_table = None
         self._deortho_inner_product_table = None
 
+        if orthonormalize:
+            # Compute the deorthonormalized tables before we orthonormalize
+            # the given basis.
+            W = V.span_of_basis( vector_basis )
+
+            # TODO: use symmetry
+            self._deortho_multiplication_table = [ [0 for j in range(n)]
+                                                   for i in range(n) ]
+            self._deortho_inner_product_table = [ [0 for j in range(n)]
+                                                  for i in range(n) ]
+
+            # Note: the Jordan and inner-products are defined in terms
+            # of the ambient basis. It's important that their arguments
+            # are in ambient coordinates as well.
+            for i in range(n):
+                for j in range(i+1):
+                    # given basis w.r.t. ambient coords
+                    q_i = vector_basis[i]
+                    q_j = vector_basis[j]
+
+                    if basis_is_matrices:
+                        q_i = _vec2mat(q_i)
+                        q_j = _vec2mat(q_j)
+
+                    elt = jordan_product(q_i, q_j)
+                    ip = inner_product(q_i, q_j)
+
+                    if basis_is_matrices:
+                        # do another mat2vec because the multiplication
+                        # table is in terms of vectors
+                        elt = _mat2vec(elt)
+
+                    # TODO: use symmetry
+                    elt = W.coordinate_vector(elt)
+                    self._deortho_multiplication_table[i][j] = elt
+                    self._deortho_multiplication_table[j][i] = elt
+                    self._deortho_inner_product_table[i][j] = ip
+                    self._deortho_inner_product_table[j][i] = ip
+
+        if self._deortho_multiplication_table is not None:
+            self._deortho_multiplication_table = tuple(map(tuple, self._deortho_multiplication_table))
+        if self._deortho_inner_product_table is not None:
+            self._deortho_inner_product_table = tuple(map(tuple, self._deortho_inner_product_table))
+
         if orthonormalize:
             from mjo.eja.eja_utils import gram_schmidt
             vector_basis = gram_schmidt(vector_basis, inner_product)
@@ -1190,8 +1237,9 @@ class RationalBasisEuclideanJordanAlgebraNg(FiniteDimensionalEuclideanJordanAlge
 
         W = V.span_of_basis( vector_basis )
 
-        mult_table = [ [0 for i in range(n)] for j in range(n) ]
-        ip_table = [ [0 for i in range(n)] for j in range(n) ]
+        # TODO: use symmetry
+        mult_table = [ [0 for j in range(n)] for i in range(n) ]
+        ip_table = [ [0 for j in range(n)] for i in range(n) ]
 
         # Note: the Jordan and inner- products are defined in terms
         # of the ambient basis. It's important that their arguments
@@ -1214,6 +1262,7 @@ class RationalBasisEuclideanJordanAlgebraNg(FiniteDimensionalEuclideanJordanAlge
                     # table is in terms of vectors
                     elt = _mat2vec(elt)
 
+                # TODO: use symmetry
                 elt = W.coordinate_vector(elt)
                 mult_table[i][j] = elt
                 mult_table[j][i] = elt
@@ -1235,29 +1284,29 @@ class RationalBasisEuclideanJordanAlgebraNg(FiniteDimensionalEuclideanJordanAlge
                          check_field,
                          check_axioms)
 
-class RationalBasisEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra):
-    r"""
-    Algebras whose basis consists of vectors with rational
-    entries. Equivalently, algebras whose multiplication tables
-    contain only rational coefficients.
-
-    When an EJA has a basis that can be made rational, we can speed up
-    the computation of its characteristic polynomial by doing it over
-    ``QQ``. All of the named EJA constructors that we provide fall
-    into this category.
-    """
     @cached_method
     def _charpoly_coefficients(self):
         r"""
-        Override the parent method with something that tries to compute
-        over a faster (non-extension) field.
-
         SETUP::
 
-            sage: from mjo.eja.eja_algebra import JordanSpinEJA
+            sage: from mjo.eja.eja_algebra import (BilinearFormEJA,
+            ....:                                  JordanSpinEJA)
 
         EXAMPLES:
 
+        The returned coefficients should be the same as if we'd never
+        orthonormalized the basis to begin with::
+
+            sage: B = matrix(QQ, [[1,   0,   0],
+            ....:                 [0,  25, -32],
+            ....:                 [0, -32,  41] ])
+            sage: J1 = BilinearFormEJA(B)
+            sage: J2 = BilinearFormEJA(B,QQ,orthonormalize=False)
+            sage: J1._charpoly_coefficients()
+            (X1^2 - 25*X2^2 + 64*X2*X3 - 41*X3^2, -2*X1)
+            sage: J2._charpoly_coefficients()
+            (X1^2 - 25*X2^2 + 64*X2*X3 - 41*X3^2, -2*X1)
+
         The base ring of the resulting polynomial coefficients is what
         it should be, and not the rationals (unless the algebra was
         already over the rationals)::
@@ -1270,30 +1319,21 @@ class RationalBasisEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebr
             Algebraic Real Field
             sage: a0.base_ring()
             Algebraic Real Field
-
         """
         if self.base_ring() is QQ:
             # There's no need to construct *another* algebra over the
             # rationals if this one is already over the rationals.
-            superclass = super(RationalBasisEuclideanJordanAlgebra, self)
+            superclass = super(RationalBasisEuclideanJordanAlgebraNg, self)
             return superclass._charpoly_coefficients()
 
-        mult_table = tuple(
-            tuple(map(lambda x: x.to_vector(), ls))
-            for ls in self._multiplication_table
-        )
-
         # Do the computation over the rationals. The answer will be
-        # the same, because our basis coordinates are (essentially)
-        # rational.
+        # the same, because all we've done is a change of basis.
         J = FiniteDimensionalEuclideanJordanAlgebra(QQ,
-                                                    mult_table,
-                                                    check_field=False,
-                                                    check_axioms=False)
+                                                    self._deortho_multiplication_table,
+                                                    self._deortho_inner_product_table)
         a = J._charpoly_coefficients()
         return tuple(map(lambda x: x.change_ring(self.base_ring()), a))
 
-
 class ConcreteEuclideanJordanAlgebra:
     r"""
     A class for the Euclidean Jordan algebras that we know by name.