X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Feja%2Feja_algebra.py;h=d06efb59b14c1c7ecf1fb3fda375168ded981e6e;hb=634dab08d886610e41dfd363a0a608a9405abe8e;hp=fb840edb93ed564c7372eacac2e3b9913a4a2350;hpb=f9e1e977040c2c3b1019a0bcc5776384a13e96c4;p=sage.d.git diff --git a/mjo/eja/eja_algebra.py b/mjo/eja/eja_algebra.py index fb840ed..d06efb5 100644 --- a/mjo/eja/eja_algebra.py +++ b/mjo/eja/eja_algebra.py @@ -5,9 +5,10 @@ are used in optimization, and have some additional nice methods beyond what can be supported in a general Jordan Algebra. """ -from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra import FiniteDimensionalAlgebra +#from sage.algebras.finite_dimensional_algebras.finite_dimensional_algebra import FiniteDimensionalAlgebra from sage.algebras.quatalg.quaternion_algebra import QuaternionAlgebra from sage.categories.finite_dimensional_algebras_with_basis import FiniteDimensionalAlgebrasWithBasis +from sage.combinat.free_module import CombinatorialFreeModule from sage.matrix.constructor import matrix from sage.misc.cachefunc import cached_method from sage.misc.prandom import choice @@ -20,50 +21,14 @@ from sage.structure.element import is_Matrix from sage.structure.category_object import normalize_names from mjo.eja.eja_element import FiniteDimensionalEuclideanJordanAlgebraElement -from mjo.eja.eja_utils import _vec2mat, _mat2vec - -class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): - @staticmethod - def __classcall_private__(cls, - field, - mult_table, - rank, - names='e', - assume_associative=False, - category=None, - natural_basis=None): - n = len(mult_table) - mult_table = [b.base_extend(field) for b in mult_table] - for b in mult_table: - b.set_immutable() - if not (is_Matrix(b) and b.dimensions() == (n, n)): - raise ValueError("input is not a multiplication table") - mult_table = tuple(mult_table) - - cat = FiniteDimensionalAlgebrasWithBasis(field) - cat.or_subcategory(category) - if assume_associative: - cat = cat.Associative() - - names = normalize_names(n, names) - - fda = super(FiniteDimensionalEuclideanJordanAlgebra, cls) - return fda.__classcall__(cls, - field, - mult_table, - rank, - assume_associative=assume_associative, - names=names, - category=cat, - natural_basis=natural_basis) - +from mjo.eja.eja_utils import _mat2vec +class FiniteDimensionalEuclideanJordanAlgebra(CombinatorialFreeModule): def __init__(self, field, mult_table, rank, - names='e', - assume_associative=False, + prefix='e', category=None, natural_basis=None): """ @@ -86,11 +51,76 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): self._rank = rank self._natural_basis = natural_basis self._multiplication_table = mult_table + if category is None: + category = FiniteDimensionalAlgebrasWithBasis(field).Unital() fda = super(FiniteDimensionalEuclideanJordanAlgebra, self) fda.__init__(field, - mult_table, - names=names, + range(len(mult_table)), + prefix=prefix, category=category) + self.print_options(bracket='') + + + def _element_constructor_(self, elt): + """ + Construct an element of this algebra from its natural + representation. + + This gets called only after the parent element _call_ method + fails to find a coercion for the argument. + + SETUP:: + + sage: from mjo.eja.eja_algebra import (JordanSpinEJA, + ....: RealCartesianProductEJA, + ....: RealSymmetricEJA) + + EXAMPLES: + + The identity in `S^n` is converted to the identity in the EJA:: + + sage: J = RealSymmetricEJA(3) + sage: I = matrix.identity(QQ,3) + sage: J(I) == J.one() + True + + This skew-symmetric matrix can't be represented in the EJA:: + + sage: J = RealSymmetricEJA(3) + sage: A = matrix(QQ,3, lambda i,j: i-j) + sage: J(A) + Traceback (most recent call last): + ... + ArithmeticError: vector is not in free module + + TESTS: + + Ensure that we can convert any element of the two non-matrix + simple algebras (whose natural representations are their usual + vector representations) back and forth faithfully:: + + sage: set_random_seed() + sage: J = RealCartesianProductEJA(5) + sage: x = J.random_element() + sage: J(x.to_vector().column()) == x + True + sage: J = JordanSpinEJA(5) + sage: x = J.random_element() + sage: J(x.to_vector().column()) == x + True + + """ + natural_basis = self.natural_basis() + if elt not in natural_basis[0].matrix_space(): + raise ValueError("not a naturally-represented algebra element") + + # Thanks for nothing! Matrix spaces aren't vector + # spaces in Sage, so we have to figure out its + # natural-basis coordinates ourselves. + V = VectorSpace(elt.base_ring(), elt.nrows()*elt.ncols()) + W = V.span_of_basis( _mat2vec(s) for s in natural_basis ) + coords = W.coordinate_vector(_mat2vec(elt)) + return self.from_vector(coords) def _repr_(self): @@ -111,9 +141,15 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): Euclidean Jordan algebra of degree 3 over Real Double Field """ + # TODO: change this to say "dimension" and fix all the tests. fmt = "Euclidean Jordan algebra of degree {} over {}" - return fmt.format(self.degree(), self.base_ring()) + return fmt.format(self.dimension(), self.base_ring()) + def product_on_basis(self, i, j): + ei = self.basis()[i] + ej = self.basis()[j] + Lei = self._multiplication_table[i] + return self.from_vector(Lei*ej.to_vector()) def _a_regular_element(self): """ @@ -153,7 +189,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): """ z = self._a_regular_element() V = self.vector_space() - V1 = V.span_of_basis( (z**k).vector() for k in range(self.rank()) ) + V1 = V.span_of_basis( (z**k).to_vector() for k in range(self.rank()) ) b = (V1.basis() + V1.complement().basis()) return V.span_of_basis(b) @@ -205,11 +241,12 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): n = self.dimension() # Construct a new algebra over a multivariate polynomial ring... - names = ['X' + str(i) for i in range(1,n+1)] + names = tuple('X' + str(i) for i in range(1,n+1)) R = PolynomialRing(self.base_ring(), names) - J = FiniteDimensionalEuclideanJordanAlgebra(R, - self._multiplication_table, - r) + J = FiniteDimensionalEuclideanJordanAlgebra( + R, + tuple(self._multiplication_table), + r) idmat = matrix.identity(J.base_ring(), n) @@ -231,11 +268,17 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): # We want the middle equivalent thing in our matrix, but use # the first equivalent thing instead so that we can pass in # standard coordinates. - x = J(W(R.gens())) - l1 = [matrix.column(W.coordinates((x**k).vector())) for k in range(r)] + x = J.from_vector(W(R.gens())) + + # Handle the zeroth power separately, because computing + # the unit element in J is mathematically suspect. + x0 = W.coordinate_vector(self.one().to_vector()) + l1 = [ x0.column() ] + l1 += [ W.coordinate_vector((x**k).to_vector()).column() + for k in range(1,r) ] l2 = [idmat.column(k-1).column() for k in range(r+1, n+1)] A_of_x = matrix.block(R, 1, n, (l1 + l2)) - xr = W.coordinates((x**r).vector()) + xr = W.coordinate_vector((x**r).to_vector()) return (A_of_x, x, xr, A_of_x.det()) @@ -265,7 +308,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): sage: J = JordanSpinEJA(3) sage: p = J.characteristic_polynomial(); p X1^2 - X2^2 - X3^2 + (-2*t)*X1 + t^2 - sage: xvec = J.one().vector() + sage: xvec = J.one().to_vector() sage: p(*xvec) t^2 - 2*t + 1 @@ -353,7 +396,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): sage: J = RealSymmetricEJA(2) sage: J.basis() - Family (e0, e1, e2) + Finite family {0: e0, 1: e1, 2: e2} sage: J.natural_basis() ( [1 0] [0 1] [0 0] @@ -364,7 +407,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): sage: J = JordanSpinEJA(2) sage: J.basis() - Family (e0, e1) + Finite family {0: e0, 1: e1} sage: J.natural_basis() ( [1] [0] @@ -373,11 +416,72 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): """ if self._natural_basis is None: - return tuple( b.vector().column() for b in self.basis() ) + return tuple( b.to_vector().column() for b in self.basis() ) else: return self._natural_basis + @cached_method + def one(self): + """ + Return the unit element of this algebra. + + SETUP:: + + sage: from mjo.eja.eja_algebra import (RealCartesianProductEJA, + ....: random_eja) + + EXAMPLES:: + + sage: J = RealCartesianProductEJA(5) + sage: J.one() + e0 + e1 + e2 + e3 + e4 + + TESTS:: + + The identity element acts like the identity:: + + sage: set_random_seed() + sage: J = random_eja() + sage: x = J.random_element() + sage: J.one()*x == x and x*J.one() == x + True + + The matrix of the unit element's operator is the identity:: + + sage: set_random_seed() + sage: J = random_eja() + sage: actual = J.one().operator().matrix() + sage: expected = matrix.identity(J.base_ring(), J.dimension()) + sage: actual == expected + True + + """ + # We can brute-force compute the matrices of the operators + # that correspond to the basis elements of this algebra. + # If some linear combination of those basis elements is the + # algebra identity, then the same linear combination of + # their matrices has to be the identity matrix. + # + # Of course, matrices aren't vectors in sage, so we have to + # appeal to the "long vectors" isometry. + oper_vecs = [ _mat2vec(g.operator().matrix()) for g in self.gens() ] + + # Now we use basis linear algebra to find the coefficients, + # of the matrices-as-vectors-linear-combination, which should + # work for the original algebra basis too. + A = matrix.column(self.base_ring(), oper_vecs) + + # We used the isometry on the left-hand side already, but we + # still need to do it for the right-hand side. Recall that we + # wanted something that summed to the identity matrix. + b = _mat2vec( matrix.identity(self.base_ring(), self.dimension()) ) + + # Now if there's an identity element in the algebra, this should work. + coeffs = A.solve_right(b) + return self.linear_combination(zip(self.gens(), coeffs)) + + def rank(self): """ Return the rank of this EJA. @@ -453,7 +557,7 @@ class FiniteDimensionalEuclideanJordanAlgebra(FiniteDimensionalAlgebra): Vector space of dimension 3 over Rational Field """ - return self.zero().vector().parent().ambient_vector_space() + return self.zero().to_vector().parent().ambient_vector_space() Element = FiniteDimensionalEuclideanJordanAlgebraElement @@ -492,18 +596,17 @@ class RealCartesianProductEJA(FiniteDimensionalEuclideanJordanAlgebra): e2 """ - @staticmethod - def __classcall_private__(cls, n, field=QQ): - # The FiniteDimensionalAlgebra constructor takes a list of - # matrices, the ith representing right multiplication by the ith - # basis element in the vector space. So if e_1 = (1,0,0), then - # right (Hadamard) multiplication of x by e_1 picks out the first + def __init__(self, n, field=QQ): + # The superclass constructor takes a list of matrices, the ith + # representing right multiplication by the ith basis element + # in the vector space. So if e_1 = (1,0,0), then right + # (Hadamard) multiplication of x by e_1 picks out the first # component of x; and likewise for the ith basis element e_i. Qs = [ matrix(field, n, n, lambda k,j: 1*(k == j == i)) for i in xrange(n) ] - fdeja = super(RealCartesianProductEJA, cls) - return fdeja.__classcall_private__(cls, field, Qs, rank=n) + fdeja = super(RealCartesianProductEJA, self) + return fdeja.__init__(field, Qs, rank=n) def inner_product(self, x, y): return _usual_ip(x,y) @@ -688,30 +791,21 @@ def _multiplication_table_from_matrix_basis(basis): dimension = basis[0].nrows() V = VectorSpace(field, dimension**2) - W = V.span( _mat2vec(s) for s in basis ) - - # Taking the span above reorders our basis (thanks, jerk!) so we - # need to put our "matrix basis" in the same order as the - # (reordered) vector basis. - S = tuple( _vec2mat(b) for b in W.basis() ) + W = V.span_of_basis( _mat2vec(s) for s in basis ) Qs = [] - for s in S: + for s in basis: # Brute force the multiplication-by-s matrix by looping # through all elements of the basis and doing the computation - # to find out what the corresponding row should be. BEWARE: - # these multiplication tables won't be symmetric! It therefore - # becomes REALLY IMPORTANT that the underlying algebra - # constructor uses ROW vectors and not COLUMN vectors. That's - # why we're computing rows here and not columns. - Q_rows = [] - for t in S: - this_row = _mat2vec((s*t + t*s)/2) - Q_rows.append(W.coordinates(this_row)) - Q = matrix(field, W.dimension(), Q_rows) + # to find out what the corresponding row should be. + Q_cols = [] + for t in basis: + this_col = _mat2vec((s*t + t*s)/2) + Q_cols.append(W.coordinates(this_col)) + Q = matrix.column(field, W.dimension(), Q_cols) Qs.append(Q) - return (Qs, S) + return Qs def _embed_complex_matrix(M): @@ -942,7 +1036,7 @@ def _unembed_quaternion_matrix(M): # The usual inner product on R^n. def _usual_ip(x,y): - return x.vector().inner_product(y.vector()) + return x.to_vector().inner_product(y.to_vector()) # The inner product used for the real symmetric simple EJA. # We keep it as a separate function because e.g. the complex @@ -976,12 +1070,12 @@ class RealSymmetricEJA(FiniteDimensionalEuclideanJordanAlgebra): TESTS: - The degree of this algebra is `(n^2 + n) / 2`:: + The dimension of this algebra is `(n^2 + n) / 2`:: sage: set_random_seed() sage: n = ZZ.random_element(1,5) sage: J = RealSymmetricEJA(n) - sage: J.degree() == (n^2 + n)/2 + sage: J.dimension() == (n^2 + n)/2 True The Jordan multiplication is what we think it is:: @@ -1001,17 +1095,15 @@ class RealSymmetricEJA(FiniteDimensionalEuclideanJordanAlgebra): True """ - @staticmethod - def __classcall_private__(cls, n, field=QQ): + def __init__(self, n, field=QQ): S = _real_symmetric_basis(n, field=field) - (Qs, T) = _multiplication_table_from_matrix_basis(S) + Qs = _multiplication_table_from_matrix_basis(S) - fdeja = super(RealSymmetricEJA, cls) - return fdeja.__classcall_private__(cls, - field, - Qs, - rank=n, - natural_basis=T) + fdeja = super(RealSymmetricEJA, self) + return fdeja.__init__(field, + Qs, + rank=n, + natural_basis=S) def inner_product(self, x, y): return _matrix_ip(x,y) @@ -1030,12 +1122,12 @@ class ComplexHermitianEJA(FiniteDimensionalEuclideanJordanAlgebra): TESTS: - The degree of this algebra is `n^2`:: + The dimension of this algebra is `n^2`:: sage: set_random_seed() sage: n = ZZ.random_element(1,5) sage: J = ComplexHermitianEJA(n) - sage: J.degree() == n^2 + sage: J.dimension() == n^2 True The Jordan multiplication is what we think it is:: @@ -1055,17 +1147,16 @@ class ComplexHermitianEJA(FiniteDimensionalEuclideanJordanAlgebra): True """ - @staticmethod - def __classcall_private__(cls, n, field=QQ): + def __init__(self, n, field=QQ): S = _complex_hermitian_basis(n) - (Qs, T) = _multiplication_table_from_matrix_basis(S) + Qs = _multiplication_table_from_matrix_basis(S) + + fdeja = super(ComplexHermitianEJA, self) + return fdeja.__init__(field, + Qs, + rank=n, + natural_basis=S) - fdeja = super(ComplexHermitianEJA, cls) - return fdeja.__classcall_private__(cls, - field, - Qs, - rank=n, - natural_basis=T) def inner_product(self, x, y): # Since a+bi on the diagonal is represented as @@ -1091,12 +1182,12 @@ class QuaternionHermitianEJA(FiniteDimensionalEuclideanJordanAlgebra): TESTS: - The degree of this algebra is `n^2`:: + The dimension of this algebra is `n^2`:: sage: set_random_seed() sage: n = ZZ.random_element(1,5) sage: J = QuaternionHermitianEJA(n) - sage: J.degree() == 2*(n^2) - n + sage: J.dimension() == 2*(n^2) - n True The Jordan multiplication is what we think it is:: @@ -1116,17 +1207,15 @@ class QuaternionHermitianEJA(FiniteDimensionalEuclideanJordanAlgebra): True """ - @staticmethod - def __classcall_private__(cls, n, field=QQ): + def __init__(self, n, field=QQ): S = _quaternion_hermitian_basis(n) - (Qs, T) = _multiplication_table_from_matrix_basis(S) + Qs = _multiplication_table_from_matrix_basis(S) - fdeja = super(QuaternionHermitianEJA, cls) - return fdeja.__classcall_private__(cls, - field, - Qs, - rank=n, - natural_basis=T) + fdeja = super(QuaternionHermitianEJA, self) + return fdeja.__init__(field, + Qs, + rank=n, + natural_basis=S) def inner_product(self, x, y): # Since a+bi+cj+dk on the diagonal is represented as @@ -1174,8 +1263,7 @@ class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra): 0 """ - @staticmethod - def __classcall_private__(cls, n, field=QQ): + def __init__(self, n, field=QQ): Qs = [] id_matrix = matrix.identity(field, n) for i in xrange(n): @@ -1192,8 +1280,8 @@ class JordanSpinEJA(FiniteDimensionalEuclideanJordanAlgebra): # The rank of the spin algebra is two, unless we're in a # one-dimensional ambient space (because the rank is bounded by # the ambient dimension). - fdeja = super(JordanSpinEJA, cls) - return fdeja.__classcall_private__(cls, field, Qs, rank=min(n,2)) + fdeja = super(JordanSpinEJA, self) + return fdeja.__init__(field, Qs, rank=min(n,2)) def inner_product(self, x, y): return _usual_ip(x,y)