from sage.matrix.constructor import matrix from sage.categories.all import FreeModules from sage.categories.map import Map class FiniteDimensionalEuclideanJordanAlgebraOperator(Map): def __init__(self, domain_eja, codomain_eja, mat): # if not ( # isinstance(domain_eja, FiniteDimensionalEuclideanJordanAlgebra) and # isinstance(codomain_eja, FiniteDimensionalEuclideanJordanAlgebra) ): # raise ValueError('(co)domains must be finite-dimensional Euclidean ' # 'Jordan algebras') F = domain_eja.base_ring() if not (F == codomain_eja.base_ring()): raise ValueError("domain and codomain must have the same base ring") if not (F == mat.base_ring()): raise ValueError("domain and matrix must have the same base ring") # We need to supply something here to avoid getting the # default Homset of the parent FiniteDimensionalAlgebra class, # which messes up e.g. equality testing. We use FreeModules(F) # instead of VectorSpaces(F) because our characteristic polynomial # algorithm will need to F to be a polynomial ring at some point. # When F is a field, FreeModules(F) returns VectorSpaces(F) anyway. parent = domain_eja.Hom(codomain_eja, FreeModules(F)) # The Map initializer will set our parent to a homset, which # is explicitly NOT what we want, because these ain't algebra # homomorphisms. super(FiniteDimensionalEuclideanJordanAlgebraOperator,self).__init__(parent) # Keep a matrix around to do all of the real work. It would # be nice if we could use a VectorSpaceMorphism instead, but # those use row vectors that we don't want to accidentally # expose to our users. self._matrix = mat def _call_(self, x): """ Allow this operator to be called only on elements of an EJA. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import JordanSpinEJA EXAMPLES:: sage: J = JordanSpinEJA(3) sage: x = J.linear_combination(zip(J.gens(),range(len(J.gens())))) sage: id = identity_matrix(J.base_ring(), J.dimension()) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: f(x) == x True """ return self.codomain().from_vector(self.matrix()*x.to_vector()) def _add_(self, other): """ Add the ``other`` EJA operator to this one. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import ( ....: JordanSpinEJA, ....: RealSymmetricEJA ) EXAMPLES: When we add two EJA operators, we get another one back:: sage: J = RealSymmetricEJA(2) sage: id = identity_matrix(J.base_ring(), J.dimension()) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: g = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: f + g Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [2 0 0] [0 2 0] [0 0 2] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... If you try to add two identical vector space operators but on different EJAs, that should blow up:: sage: J1 = RealSymmetricEJA(2) sage: id1 = identity_matrix(J1.base_ring(), 3) sage: J2 = JordanSpinEJA(3) sage: id2 = identity_matrix(J2.base_ring(), 3) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J1,J1,id1) sage: g = FiniteDimensionalEuclideanJordanAlgebraOperator(J2,J2,id2) sage: f + g Traceback (most recent call last): ... TypeError: unsupported operand parent(s) for +: ... """ return FiniteDimensionalEuclideanJordanAlgebraOperator( self.domain(), self.codomain(), self.matrix() + other.matrix()) def _composition_(self, other, homset): """ Compose two EJA operators to get another one (and NOT a formal composite object) back. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import ( ....: JordanSpinEJA, ....: HadamardEJA, ....: RealSymmetricEJA) EXAMPLES:: sage: J1 = JordanSpinEJA(3) sage: J2 = HadamardEJA(2) sage: J3 = RealSymmetricEJA(1) sage: mat1 = matrix(QQ, [[1,2,3], ....: [4,5,6]]) sage: mat2 = matrix(QQ, [[7,8]]) sage: g = FiniteDimensionalEuclideanJordanAlgebraOperator(J1, ....: J2, ....: mat1) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J2, ....: J3, ....: mat2) sage: f*g Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [39 54 69] Domain: Euclidean Jordan algebra of dimension 3 over Rational Field Codomain: Euclidean Jordan algebra of dimension 1 over Rational Field """ return FiniteDimensionalEuclideanJordanAlgebraOperator( other.domain(), self.codomain(), self.matrix()*other.matrix()) def __eq__(self, other): if self.domain() != other.domain(): return False if self.codomain() != other.codomain(): return False if self.matrix() != other.matrix(): return False return True def __invert__(self): """ Invert this EJA operator. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES:: sage: J = RealSymmetricEJA(2) sage: id = identity_matrix(J.base_ring(), J.dimension()) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: ~f Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [1 0 0] [0 1 0] [0 0 1] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... """ return FiniteDimensionalEuclideanJordanAlgebraOperator( self.codomain(), self.domain(), ~self.matrix()) def __mul__(self, other): """ Compose two EJA operators, or scale myself by an element of the ambient vector space. We need to override the real ``__mul__`` function to prevent the coercion framework from throwing an error when it fails to convert a base ring element into a morphism. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES: We can scale an operator on a rational algebra by a rational number:: sage: J = RealSymmetricEJA(2) sage: e0,e1,e2 = J.gens() sage: x = 2*e0 + 4*e1 + 16*e2 sage: x.operator() Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [ 2 2 0] [ 2 9 2] [ 0 2 16] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... sage: x.operator()*(1/2) Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [ 1 1 0] [ 1 9/2 1] [ 0 1 8] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... """ try: if other in self.codomain().base_ring(): return FiniteDimensionalEuclideanJordanAlgebraOperator( self.domain(), self.codomain(), self.matrix()*other) except NotImplementedError: # This can happen with certain arguments if the base_ring() # is weird and doesn't know how to test membership. pass # This should eventually delegate to _composition_ after performing # some sanity checks for us. mor = super(FiniteDimensionalEuclideanJordanAlgebraOperator,self) return mor.__mul__(other) def _neg_(self): """ Negate this EJA operator. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES:: sage: J = RealSymmetricEJA(2) sage: id = identity_matrix(J.base_ring(), J.dimension()) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: -f Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [-1 0 0] [ 0 -1 0] [ 0 0 -1] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... """ return FiniteDimensionalEuclideanJordanAlgebraOperator( self.domain(), self.codomain(), -self.matrix()) def __pow__(self, n): """ Raise this EJA operator to the power ``n``. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA TESTS: Ensure that we get back another EJA operator that can be added, subtracted, et cetera:: sage: J = RealSymmetricEJA(2) sage: id = identity_matrix(J.base_ring(), J.dimension()) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: f^0 + f^1 + f^2 Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [3 0 0] [0 3 0] [0 0 3] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... """ if (n == 1): return self elif (n == 0): # Raising a vector space morphism to the zero power gives # you back a special IdentityMorphism that is useless to us. rows = self.codomain().dimension() cols = self.domain().dimension() mat = matrix.identity(self.base_ring(), rows, cols) else: mat = self.matrix()**n return FiniteDimensionalEuclideanJordanAlgebraOperator( self.domain(), self.codomain(), mat) def _repr_(self): r""" A text representation of this linear operator on a Euclidean Jordan Algebra. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import JordanSpinEJA EXAMPLES:: sage: J = JordanSpinEJA(2) sage: id = identity_matrix(J.base_ring(), J.dimension()) sage: FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [1 0] [0 1] Domain: Euclidean Jordan algebra of dimension 2 over Rational Field Codomain: Euclidean Jordan algebra of dimension 2 over Rational Field """ msg = ("Linear operator between finite-dimensional Euclidean Jordan " "algebras represented by the matrix:\n", "{!r}\n", "Domain: {}\n", "Codomain: {}") return ''.join(msg).format(self.matrix(), self.domain(), self.codomain()) def _sub_(self, other): """ Subtract ``other`` from this EJA operator. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES:: sage: J = RealSymmetricEJA(2) sage: id = identity_matrix(J.base_ring(),J.dimension()) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,id) sage: f - (f*2) Linear operator between finite-dimensional Euclidean Jordan algebras represented by the matrix: [-1 0 0] [ 0 -1 0] [ 0 0 -1] Domain: Euclidean Jordan algebra of dimension 3 over... Codomain: Euclidean Jordan algebra of dimension 3 over... """ return (self + (-other)) def inverse(self): """ Return the inverse of this operator, if it exists. The reason this method is not simply an alias for the built-in :meth:`__invert__` is that the built-in inversion is a bit magic since it's intended to be a unary operator. If we alias ``inverse`` to ``__invert__``, then we wind up having to call e.g. ``A.inverse`` without parentheses. SETUP:: sage: from mjo.eja.eja_algebra import RealSymmetricEJA, random_eja EXAMPLES:: sage: J = RealSymmetricEJA(2) sage: x = sum(J.gens()) sage: x.operator().inverse().matrix() [3/2 -1 1/2] [ -1 2 -1] [1/2 -1 3/2] sage: x.operator().matrix().inverse() [3/2 -1 1/2] [ -1 2 -1] [1/2 -1 3/2] TESTS: The identity operator is its own inverse:: sage: set_random_seed() sage: J = random_eja() sage: idJ = J.one().operator() sage: idJ.inverse() == idJ True The inverse of the inverse is the operator we started with:: sage: set_random_seed() sage: x = random_eja().random_element() sage: L = x.operator() sage: not L.is_invertible() or (L.inverse().inverse() == L) True """ return ~self def is_invertible(self): """ Return whether or not this operator is invertible. SETUP:: sage: from mjo.eja.eja_algebra import (RealSymmetricEJA, ....: TrivialEJA, ....: random_eja) EXAMPLES:: sage: J = RealSymmetricEJA(2) sage: x = sum(J.gens()) sage: x.operator().matrix() [ 1 1/2 0] [1/2 1 1/2] [ 0 1/2 1] sage: x.operator().matrix().is_invertible() True sage: x.operator().is_invertible() True The zero operator is invertible in a trivial algebra:: sage: J = TrivialEJA() sage: J.zero().operator().is_invertible() True TESTS: The identity operator is always invertible:: sage: set_random_seed() sage: J = random_eja() sage: J.one().operator().is_invertible() True The zero operator is never invertible in a nontrivial algebra:: sage: set_random_seed() sage: J = random_eja() sage: not J.is_trivial() and J.zero().operator().is_invertible() False """ return self.matrix().is_invertible() def matrix(self): """ Return the matrix representation of this operator with respect to the default bases of its (co)domain. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES:: sage: J = RealSymmetricEJA(2) sage: mat = matrix(J.base_ring(), J.dimension(), range(9)) sage: f = FiniteDimensionalEuclideanJordanAlgebraOperator(J,J,mat) sage: f.matrix() [0 1 2] [3 4 5] [6 7 8] """ return self._matrix def minimal_polynomial(self): """ Return the minimal polynomial of this linear operator, in the variable ``t``. SETUP:: sage: from mjo.eja.eja_operator import FiniteDimensionalEuclideanJordanAlgebraOperator sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES:: sage: J = RealSymmetricEJA(3) sage: J.one().operator().minimal_polynomial() t - 1 """ # The matrix method returns a polynomial in 'x' but want one in 't'. return self.matrix().minimal_polynomial().change_variable_name('t') def spectral_decomposition(self): """ Return the spectral decomposition of this operator as a list of (eigenvalue, orthogonal projector) pairs. This is the unique spectral decomposition, up to the order of the projection operators, with distinct eigenvalues. So, the projections are generally onto subspaces of dimension greater than one. SETUP:: sage: from mjo.eja.eja_algebra import RealSymmetricEJA EXAMPLES:: sage: J = RealSymmetricEJA(4,AA) sage: x = sum(J.gens()) sage: A = x.subalgebra_generated_by(orthonormalize_basis=True) sage: L0x = A(x).operator() sage: sd = L0x.spectral_decomposition() sage: l0 = sd[0][0] sage: l1 = sd[1][0] sage: P0 = sd[0][1] sage: P1 = sd[1][1] sage: P0*l0 + P1*l1 == L0x True sage: P0 + P1 == P0^0 # the identity True sage: P0^2 == P0 True sage: P1^2 == P1 True sage: P0*P1 == A.zero().operator() True sage: P1*P0 == A.zero().operator() True """ if not self.matrix().is_symmetric(): raise ValueError('algebra basis is not orthonormal') D,P = self.matrix().jordan_form(subdivide=False,transformation=True) eigenvalues = D.diagonal() us = P.columns() projectors = [] for i in range(len(us)): # they won't be normalized, but they have to be # for the spectral theorem to work. us[i] = us[i]/us[i].norm() mat = us[i].column()*us[i].row() Pi = FiniteDimensionalEuclideanJordanAlgebraOperator( self.domain(), self.codomain(), mat) projectors.append(Pi) return list(zip(eigenvalues, projectors))