X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_algebra.py;h=02cf32c7d27ff93728ad548a964ff6bee5f277d4;hb=f1752f438aa849da1e909c67cac2cd7ac670e86f;hp=80147dfa62ca0ab408d6da2f5a7ee047eb788f05;hpb=1bae84b14cb9352727a670be01b5ef06dbcc7328;p=sage.d.git diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index 80147df..02cf32c 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -5,7 +5,7 @@ are used in optimization, and have some additional nice methods beyond what can be supported in a general Jordan Algebra. """ -from itertools import repeat +from itertools import izip, repeat from sage.algebras.quatalg.quaternion_algebra import QuaternionAlgebra from sage.categories.magmatic_algebras import MagmaticAlgebras @@ -61,9 +61,6 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule): self._rank = rank self._natural_basis = natural_basis - # TODO: HACK for the charpoly.. needs redesign badly. - self._basis_normalizers = None - if category is None: category = MagmaticAlgebras(field).FiniteDimensional() category = category.WithBasis().Unital() @@ -259,19 +256,6 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule): store the trace/determinant (a_{r-1} and a_{0} respectively) separate from the entire characteristic polynomial. """ - if self._basis_normalizers is not None: - # Must be a matrix class? - # WARNING/TODO: this whole mess is mis-designed. - n = self.natural_basis_space().nrows() - field = self.base_ring().base_ring() # yeeeeaaaahhh - J = self.__class__(n, field, False) - (_,x,_,_) = J._charpoly_matrix_system() - p = J._charpoly_coeff(i) - # p might be missing some vars, have to substitute "optionally" - pairs = zip(x.base_ring().gens(), self._basis_normalizers) - substitutions = { v: v*c for (v,c) in pairs } - return p.subs(substitutions) - (A_of_x, x, xr, detA) = self._charpoly_matrix_system() R = A_of_x.base_ring() if i >= self.rank(): @@ -425,7 +409,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule): # assign a[r] goes out-of-bounds. a.append(1) # corresponds to x^r - return sum( a[k]*(t**k) for k in range(len(a)) ) + return sum( a[k]*(t**k) for k in xrange(len(a)) ) def inner_product(self, x, y): @@ -507,7 +491,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule): """ M = list(self._multiplication_table) # copy - for i in range(len(M)): + for i in xrange(len(M)): # M had better be "square" M[i] = [self.monomial(i)] + M[i] M = [["*"] + list(self.gens())] + M @@ -815,8 +799,8 @@ class RealCartesianProductEJA(FiniteDimensionalEuclideanJordanAlgebra): """ def __init__(self, n, field=QQ, **kwargs): V = VectorSpace(field, n) - mult_table = [ [ V.gen(i)*(i == j) for j in range(n) ] - for i in range(n) ] + mult_table = [ [ V.gen(i)*(i == j) for j in xrange(n) ] + for i in xrange(n) ] fdeja = super(RealCartesianProductEJA, self) return fdeja.__init__(field, mult_table, rank=n, **kwargs) @@ -901,14 +885,20 @@ class MatrixEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra): # field can have dimension 4 (quaternions) too. return 2 - @classmethod - def _denormalized_basis(cls, n, field): - raise NotImplementedError + def __init__(self, field, basis, rank, normalize_basis=True, **kwargs): + """ + Compared to the superclass constructor, we take a basis instead of + a multiplication table because the latter can be computed in terms + of the former when the product is known (like it is here). + """ + # Used in this class's fast _charpoly_coeff() override. + self._basis_normalizers = None - def __init__(self, n, field=QQ, normalize_basis=True, **kwargs): - S = self._denormalized_basis(n, field) + # We're going to loop through this a few times, so now's a good + # time to ensure that it isn't a generator expression. + basis = tuple(basis) - if n > 1 and normalize_basis: + if rank > 1 and normalize_basis: # We'll need sqrt(2) to normalize the basis, and this # winds up in the multiplication table, so the whole # algebra needs to be over the field extension. @@ -917,21 +907,47 @@ class MatrixEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra): p = z**2 - 2 if p.is_irreducible(): field = NumberField(p, 'sqrt2', embedding=RLF(2).sqrt()) - S = [ s.change_ring(field) for s in S ] + basis = tuple( s.change_ring(field) for s in basis ) self._basis_normalizers = tuple( - ~(self.natural_inner_product(s,s).sqrt()) for s in S ) - S = tuple( s*c for (s,c) in zip(S,self._basis_normalizers) ) + ~(self.natural_inner_product(s,s).sqrt()) for s in basis ) + basis = tuple(s*c for (s,c) in izip(basis,self._basis_normalizers)) - Qs = self.multiplication_table_from_matrix_basis(S) + Qs = self.multiplication_table_from_matrix_basis(basis) fdeja = super(MatrixEuclideanJordanAlgebra, self) return fdeja.__init__(field, Qs, - rank=n, - natural_basis=S, + rank=rank, + natural_basis=basis, **kwargs) + @cached_method + def _charpoly_coeff(self, i): + """ + Override the parent method with something that tries to compute + over a faster (non-extension) field. + """ + if self._basis_normalizers is None: + # We didn't normalize, so assume that the basis we started + # with had entries in a nice field. + return super(MatrixEuclideanJordanAlgebra, self)._charpoly_coeff(i) + else: + basis = ( (b/n) for (b,n) in izip(self.natural_basis(), + self._basis_normalizers) ) + field = self.base_ring().base_ring() # yeeeaahhhhhhh + J = MatrixEuclideanJordanAlgebra(field, + basis, + self.rank(), + normalize_basis=False) + (_,x,_,_) = J._charpoly_matrix_system() + p = J._charpoly_coeff(i) + # p might be missing some vars, have to substitute "optionally" + pairs = izip(x.base_ring().gens(), self._basis_normalizers) + substitutions = { v: v*c for (v,c) in pairs } + return p.subs(substitutions) + + @staticmethod def multiplication_table_from_matrix_basis(basis): """ @@ -953,9 +969,9 @@ class MatrixEuclideanJordanAlgebra(FiniteDimensionalEuclideanJordanAlgebra): V = VectorSpace(field, dimension**2) W = V.span_of_basis( _mat2vec(s) for s in basis ) n = len(basis) - mult_table = [[W.zero() for j in range(n)] for i in range(n)] - for i in range(n): - for j in range(n): + mult_table = [[W.zero() for j in xrange(n)] for i in xrange(n)] + for i in xrange(n): + for j in xrange(n): mat_entry = (basis[i]*basis[j] + basis[j]*basis[i])/2 mult_table[i][j] = W.coordinate_vector(_mat2vec(mat_entry)) @@ -1011,24 +1027,16 @@ class RealMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra): @staticmethod def real_embed(M): """ - Embed the matrix ``M`` into a space of real matrices. - - The matrix ``M`` can have entries in any field at the moment: - the real numbers, complex numbers, or quaternions. And although - they are not a field, we can probably support octonions at some - point, too. This function returns a real matrix that "acts like" - the original with respect to matrix multiplication; i.e. - - real_embed(M*N) = real_embed(M)*real_embed(N) - + The identity function, for embedding real matrices into real + matrices. """ return M - @staticmethod def real_unembed(M): """ - The inverse of :meth:`real_embed`. + The identity function, for unembedding real matrices from real + matrices. """ return M @@ -1133,7 +1141,7 @@ class RealSymmetricEJA(RealMatrixEuclideanJordanAlgebra): else: Sij = Eij + Eij.transpose() S.append(Sij) - return tuple(S) + return S @staticmethod @@ -1141,6 +1149,10 @@ class RealSymmetricEJA(RealMatrixEuclideanJordanAlgebra): return 4 # Dimension 10 + def __init__(self, n, field=QQ, **kwargs): + basis = self._denormalized_basis(n, field) + super(RealSymmetricEJA, self).__init__(field, basis, n, **kwargs) + class ComplexMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra): @staticmethod @@ -1353,6 +1365,7 @@ class ComplexHermitianEJA(ComplexMatrixEuclideanJordanAlgebra): True """ + @classmethod def _denormalized_basis(cls, n, field): """ @@ -1404,9 +1417,13 @@ class ComplexHermitianEJA(ComplexMatrixEuclideanJordanAlgebra): # Since we embedded these, we can drop back to the "field" that we # started with instead of the complex extension "F". - return tuple( s.change_ring(field) for s in S ) + return ( s.change_ring(field) for s in S ) + def __init__(self, n, field=QQ, **kwargs): + basis = self._denormalized_basis(n,field) + super(ComplexHermitianEJA,self).__init__(field, basis, n, **kwargs) + class QuaternionMatrixEuclideanJordanAlgebra(MatrixEuclideanJordanAlgebra): @staticmethod @@ -1682,8 +1699,12 @@ class QuaternionHermitianEJA(QuaternionMatrixEuclideanJordanAlgebra): S.append(Sij_J) Sij_K = cls.real_embed(K*Eij - K*Eij.transpose()) S.append(Sij_K) - return tuple(S) + return S + + def __init__(self, n, field=QQ, **kwargs): + basis = self._denormalized_basis(n,field) + super(QuaternionHermitianEJA,self).__init__(field, basis, n, **kwargs) class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra): @@ -1726,9 +1747,9 @@ class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra): """ def __init__(self, n, field=QQ, **kwargs): V = VectorSpace(field, n) - mult_table = [[V.zero() for j in range(n)] for i in range(n)] - for i in range(n): - for j in range(n): + mult_table = [[V.zero() for j in xrange(n)] for i in xrange(n)] + for i in xrange(n): + for j in xrange(n): x = V.gen(i) y = V.gen(j) x0 = x[0]