Representations of the Symmetric Group¶
Todo
construct the product of two irreducible representations.
implement Induction/Restriction of representations.
Warning
This code uses a different convention than in Sagan’s book “The Symmetric Group”
- class sage.combinat.symmetric_group_representations.GarsiaProcesiModule(SGA, shape)[source]¶
Bases:
UniqueRepresentation
,QuotientRing_generic
,SymmetricGroupRepresentation
A Garsia-Procesi module.
Let
be a partition of and be a commutative ring. The Garsia-Procesi module is defined by , wherewith
being the -the elementary symmetric function and , is the Tanisaki ideal.If we consider
, then the Garsia-Procesi module has the following interpretation. Let denote the (complex type A) flag variety. Consider the Springer fiber associated to a nilpotent matrix with Jordan blocks sizes . Springer showed that the cohomology ring admits a graded -action that agrees with the induced representation of the sign representation of the Young subgroup . From work of De Concini and Procesi, this -representation is isomorphic to . Moreover, the graded Frobenius image is known to be a modified Hall-Littlewood polynomial.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(QQ, 7) sage: GP421 = SGA.garsia_procesi_module([4, 2, 1]) sage: GP421.dimension() 105 sage: v = GP421.an_element(); v -gp1 - gp2 - gp3 - gp4 - gp5 - gp6 sage: SGA.an_element() * v -6*gp1 - 6*gp2 - 6*gp3 - 6*gp4 - 6*gp5 - 5*gp6
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(QQ, Integer(7)) >>> GP421 = SGA.garsia_procesi_module([Integer(4), Integer(2), Integer(1)]) >>> GP421.dimension() 105 >>> v = GP421.an_element(); v -gp1 - gp2 - gp3 - gp4 - gp5 - gp6 >>> SGA.an_element() * v -6*gp1 - 6*gp2 - 6*gp3 - 6*gp4 - 6*gp5 - 5*gp6
We verify the result is a modified Hall-Littlewood polynomial by using the
Hall-Littlewood polynomials, replacing and multiplying by the smallest power of so the coefficients are again polynomials:sage: GP421.graded_frobenius_image() q^4*s[4, 2, 1] + q^3*s[4, 3] + q^3*s[5, 1, 1] + (q^3+q^2)*s[5, 2] + (q^2+q)*s[6, 1] + s[7] sage: R.<q> = QQ[] sage: Sym = SymmetricFunctions(R) sage: s = Sym.s() sage: Qp = Sym.hall_littlewood(q).Qp() sage: mHL = s(Qp[4,2,1]); mHL s[4, 2, 1] + q*s[4, 3] + q*s[5, 1, 1] + (q^2+q)*s[5, 2] + (q^3+q^2)*s[6, 1] + q^4*s[7] sage: mHL.map_coefficients(lambda c: R(q^4*c(q^-1))) q^4*s[4, 2, 1] + q^3*s[4, 3] + q^3*s[5, 1, 1] + (q^3+q^2)*s[5, 2] + (q^2+q)*s[6, 1] + s[7]
>>> from sage.all import * >>> GP421.graded_frobenius_image() q^4*s[4, 2, 1] + q^3*s[4, 3] + q^3*s[5, 1, 1] + (q^3+q^2)*s[5, 2] + (q^2+q)*s[6, 1] + s[7] >>> R = QQ['q']; (q,) = R._first_ngens(1) >>> Sym = SymmetricFunctions(R) >>> s = Sym.s() >>> Qp = Sym.hall_littlewood(q).Qp() >>> mHL = s(Qp[Integer(4),Integer(2),Integer(1)]); mHL s[4, 2, 1] + q*s[4, 3] + q*s[5, 1, 1] + (q^2+q)*s[5, 2] + (q^3+q^2)*s[6, 1] + q^4*s[7] >>> mHL.map_coefficients(lambda c: R(q**Integer(4)*c(q**-Integer(1)))) q^4*s[4, 2, 1] + q^3*s[4, 3] + q^3*s[5, 1, 1] + (q^3+q^2)*s[5, 2] + (q^2+q)*s[6, 1] + s[7]
We show that the maximal degree component corresponds to the Yamanouchi words of content
:sage: B = GP421.graded_decomposition(4).basis() sage: top_deg = [Word([i+1 for i in b.lift().lift().exponents()[0]]) for b in B] sage: yamanouchi = [P.to_packed_word() for P in OrderedSetPartitions(range(7), [4, 2, 1]) ....: if P.to_packed_word().reversal().is_yamanouchi()] sage: set(top_deg) == set(yamanouchi) True
>>> from sage.all import * >>> B = GP421.graded_decomposition(Integer(4)).basis() >>> top_deg = [Word([i+Integer(1) for i in b.lift().lift().exponents()[Integer(0)]]) for b in B] >>> yamanouchi = [P.to_packed_word() for P in OrderedSetPartitions(range(Integer(7)), [Integer(4), Integer(2), Integer(1)]) ... if P.to_packed_word().reversal().is_yamanouchi()] >>> set(top_deg) == set(yamanouchi) True
- class Element(parent, rep, reduce=True)[source]¶
Bases:
QuotientRingElement
- degree()[source]¶
Return the degree of
self
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(3), 4) sage: GP22 = SGA.garsia_procesi_module([2, 2]) sage: for b in GP22.basis(): ....: print(b, b.degree()) gp2*gp3 2 gp1*gp3 2 gp3 1 gp2 1 gp1 1 1 0 sage: v = sum(GP22.basis()) sage: v.degree() 2
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(3)), Integer(4)) >>> GP22 = SGA.garsia_procesi_module([Integer(2), Integer(2)]) >>> for b in GP22.basis(): ... print(b, b.degree()) gp2*gp3 2 gp1*gp3 2 gp3 1 gp2 1 gp1 1 1 0 >>> v = sum(GP22.basis()) >>> v.degree() 2
- homogeneous_degree()[source]¶
Return the (homogeneous) degree of
self
if homogeneous otherwise raise an error.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(2), 4) sage: GP31 = SGA.garsia_procesi_module([3, 1]) sage: for b in GP31.basis(): ....: print(b, b.homogeneous_degree()) gp3 1 gp2 1 gp1 1 1 0 sage: v = sum(GP31.basis()); v gp1 + gp2 + gp3 + 1 sage: v.homogeneous_degree() Traceback (most recent call last): ... ValueError: element is not homogeneous
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(2)), Integer(4)) >>> GP31 = SGA.garsia_procesi_module([Integer(3), Integer(1)]) >>> for b in GP31.basis(): ... print(b, b.homogeneous_degree()) gp3 1 gp2 1 gp1 1 1 0 >>> v = sum(GP31.basis()); v gp1 + gp2 + gp3 + 1 >>> v.homogeneous_degree() Traceback (most recent call last): ... ValueError: element is not homogeneous
- monomial_coefficients(copy=None)[source]¶
Return the monomial coefficients of
self
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(3), 4) sage: GP31 = SGA.garsia_procesi_module([3, 1]) sage: v = GP31.an_element(); v -gp1 - gp2 - gp3 sage: v.monomial_coefficients() {0: 2, 1: 2, 2: 2, 3: 0}
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(3)), Integer(4)) >>> GP31 = SGA.garsia_procesi_module([Integer(3), Integer(1)]) >>> v = GP31.an_element(); v -gp1 - gp2 - gp3 >>> v.monomial_coefficients() {0: 2, 1: 2, 2: 2, 3: 0}
- to_vector(order=None)[source]¶
Return
self
as a (dense) free module vector.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(3), 4) sage: GP22 = SGA.garsia_procesi_module([2, 2]) sage: v = GP22.an_element(); v -gp1 - gp2 - gp3 sage: v.to_vector() (0, 0, 2, 2, 2, 0)
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(3)), Integer(4)) >>> GP22 = SGA.garsia_procesi_module([Integer(2), Integer(2)]) >>> v = GP22.an_element(); v -gp1 - gp2 - gp3 >>> v.to_vector() (0, 0, 2, 2, 2, 0)
- basis()[source]¶
Return a basis of
self
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(QQ, 4) sage: GP = SGA.garsia_procesi_module([2, 2]) sage: GP.basis() Family (gp2*gp3, gp1*gp3, gp3, gp2, gp1, 1)
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(QQ, Integer(4)) >>> GP = SGA.garsia_procesi_module([Integer(2), Integer(2)]) >>> GP.basis() Family (gp2*gp3, gp1*gp3, gp3, gp2, gp1, 1)
- dimension()[source]¶
Return the dimension of
self
.The graded Frobenius character of the Garsia-Procesi module
is given by the modified Hall-Littlewood polynomial .EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(QQ, 5) sage: Sym = SymmetricFunctions(QQ) sage: s = Sym.s() sage: Qp = Sym.hall_littlewood(1).Qp() sage: for la in Partitions(5): ....: print(SGA.garsia_procesi_module(la).dimension(), ....: sum(c * StandardTableaux(la).cardinality() ....: for la, c in s(Qp[la]))) 1 1 5 5 10 10 20 20 30 30 60 60 120 120
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(QQ, Integer(5)) >>> Sym = SymmetricFunctions(QQ) >>> s = Sym.s() >>> Qp = Sym.hall_littlewood(Integer(1)).Qp() >>> for la in Partitions(Integer(5)): ... print(SGA.garsia_procesi_module(la).dimension(), ... sum(c * StandardTableaux(la).cardinality() ... for la, c in s(Qp[la]))) 1 1 5 5 10 10 20 20 30 30 60 60 120 120
- get_order()[source]¶
Return the order of the elements in the basis.
EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(QQ, 4) sage: GP = SGA.garsia_procesi_module([2, 2]) sage: GP.get_order() (0, 1, 2, 3, 4, 5)
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(QQ, Integer(4)) >>> GP = SGA.garsia_procesi_module([Integer(2), Integer(2)]) >>> GP.get_order() (0, 1, 2, 3, 4, 5)
- graded_brauer_character()[source]¶
Return the graded Brauer character of
self
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(2), 5) sage: GP311 = SGA.garsia_procesi_module([3, 1, 1]) sage: GP311.graded_brauer_character() (6*q^3 + 9*q^2 + 4*q + 1, q + 1, q^3 - q^2 - q + 1)
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(2)), Integer(5)) >>> GP311 = SGA.garsia_procesi_module([Integer(3), Integer(1), Integer(1)]) >>> GP311.graded_brauer_character() (6*q^3 + 9*q^2 + 4*q + 1, q + 1, q^3 - q^2 - q + 1)
- graded_character()[source]¶
Return the graded character of
self
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(QQ, 5) sage: GP = SGA.garsia_procesi_module([2, 2, 1]) sage: gchi = GP.graded_character(); gchi (5*q^4 + 11*q^3 + 9*q^2 + 4*q + 1, -q^4 + q^3 + 3*q^2 + 2*q + 1, q^4 - q^3 + q^2 + 1, -q^4 - q^3 + q + 1, -q^4 + q^3 - q + 1, q^4 - q^3 - q^2 + 1, q^3 - q^2 - q + 1) sage: R.<q> = QQ[] sage: gchi == sum(q^d * D.character() ....: for d, D in GP.graded_decomposition().items()) True
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(QQ, Integer(5)) >>> GP = SGA.garsia_procesi_module([Integer(2), Integer(2), Integer(1)]) >>> gchi = GP.graded_character(); gchi (5*q^4 + 11*q^3 + 9*q^2 + 4*q + 1, -q^4 + q^3 + 3*q^2 + 2*q + 1, q^4 - q^3 + q^2 + 1, -q^4 - q^3 + q + 1, -q^4 + q^3 - q + 1, q^4 - q^3 - q^2 + 1, q^3 - q^2 - q + 1) >>> R = QQ['q']; (q,) = R._first_ngens(1) >>> gchi == sum(q**d * D.character() ... for d, D in GP.graded_decomposition().items()) True
- graded_decomposition(k=None)[source]¶
Return the decomposition of
self
as a direct sum of representations given by a fixed grading.INPUT:
k
– (optional) integer; if given, return the -th graded part
EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(2), 5) sage: GP32 = SGA.garsia_procesi_module([3, 2]) sage: decomp = GP32.graded_decomposition(); decomp {0: Subrepresentation with basis {0} of Garsia-Procesi ..., 1: Subrepresentation with basis {0, 1, 2, 3} of Garsia-Procesi ..., 2: Subrepresentation with basis {0, 1, 2, 3, 4} of Garsia-Procesi ...} sage: decomp[2] is GP32.graded_decomposition(2) True sage: GP32.graded_decomposition(10) Subrepresentation with basis {} of Garsia-Procesi module of shape [3, 2] over Finite Field of size 2
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(2)), Integer(5)) >>> GP32 = SGA.garsia_procesi_module([Integer(3), Integer(2)]) >>> decomp = GP32.graded_decomposition(); decomp {0: Subrepresentation with basis {0} of Garsia-Procesi ..., 1: Subrepresentation with basis {0, 1, 2, 3} of Garsia-Procesi ..., 2: Subrepresentation with basis {0, 1, 2, 3, 4} of Garsia-Procesi ...} >>> decomp[Integer(2)] is GP32.graded_decomposition(Integer(2)) True >>> GP32.graded_decomposition(Integer(10)) Subrepresentation with basis {} of Garsia-Procesi module of shape [3, 2] over Finite Field of size 2
- graded_frobenius_image()[source]¶
Return the graded Frobenius image of
self
.The graded Frobenius image is the sum of the
frobenius_image()
of each graded component, which is known to result in the modified Hall-Littlewood polynomial .EXAMPLES:
We verify that the result is the modified Hall-Littlewood polynomial for
:sage: R.<q> = QQ[] sage: Sym = SymmetricFunctions(R) sage: s = Sym.s() sage: Qp = Sym.hall_littlewood(q).Qp() sage: SGA = SymmetricGroupAlgebra(QQ, 5) sage: for la in Partitions(5): ....: f = SGA.garsia_procesi_module(la).graded_frobenius_image() ....: d = f[la].degree() ....: assert f.map_coefficients(lambda c: R(c(~q)*q^d)) == s(Qp[la])
>>> from sage.all import * >>> R = QQ['q']; (q,) = R._first_ngens(1) >>> Sym = SymmetricFunctions(R) >>> s = Sym.s() >>> Qp = Sym.hall_littlewood(q).Qp() >>> SGA = SymmetricGroupAlgebra(QQ, Integer(5)) >>> for la in Partitions(Integer(5)): ... f = SGA.garsia_procesi_module(la).graded_frobenius_image() ... d = f[la].degree() ... assert f.map_coefficients(lambda c: R(c(~q)*q**d)) == s(Qp[la])
- graded_representation_matrix(elt, q=None)[source]¶
Return the matrix corresponding to the left action of the symmetric group (algebra) element
elt
onself
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(GF(3), 3) sage: GP = SGA.garsia_procesi_module([1, 1, 1]) sage: elt = SGA.an_element(); elt [1, 2, 3] + 2*[1, 3, 2] + [3, 1, 2] sage: X = GP.graded_representation_matrix(elt); X [ 0 0 0 0 0 0] [ 0 q^2 0 0 0 0] [ 0 q^2 0 0 0 0] [ 0 0 0 q 0 0] [ 0 0 0 q 0 0] [ 0 0 0 0 0 1] sage: X.parent() Full MatrixSpace of 6 by 6 dense matrices over Univariate Polynomial Ring in q over Finite Field of size 3 sage: R.<q> = GF(3)[] sage: t = R.quotient([q^2+2*q+1]).gen() sage: GP.graded_representation_matrix(elt, t) [ 0 0 0 0 0 0] [ 0 qbar + 2 0 0 0 0] [ 0 qbar + 2 0 0 0 0] [ 0 0 0 qbar 0 0] [ 0 0 0 qbar 0 0] [ 0 0 0 0 0 1]
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(GF(Integer(3)), Integer(3)) >>> GP = SGA.garsia_procesi_module([Integer(1), Integer(1), Integer(1)]) >>> elt = SGA.an_element(); elt [1, 2, 3] + 2*[1, 3, 2] + [3, 1, 2] >>> X = GP.graded_representation_matrix(elt); X [ 0 0 0 0 0 0] [ 0 q^2 0 0 0 0] [ 0 q^2 0 0 0 0] [ 0 0 0 q 0 0] [ 0 0 0 q 0 0] [ 0 0 0 0 0 1] >>> X.parent() Full MatrixSpace of 6 by 6 dense matrices over Univariate Polynomial Ring in q over Finite Field of size 3 >>> R = GF(Integer(3))['q']; (q,) = R._first_ngens(1) >>> t = R.quotient([q**Integer(2)+Integer(2)*q+Integer(1)]).gen() >>> GP.graded_representation_matrix(elt, t) [ 0 0 0 0 0 0] [ 0 qbar + 2 0 0 0 0] [ 0 qbar + 2 0 0 0 0] [ 0 0 0 qbar 0 0] [ 0 0 0 qbar 0 0] [ 0 0 0 0 0 1]
- one_basis()[source]¶
Return the index of the basis element
.EXAMPLES:
sage: SGA = SymmetricGroupAlgebra(QQ, 4) sage: GP = SGA.garsia_procesi_module([2, 2]) sage: GP.one_basis() 5
>>> from sage.all import * >>> SGA = SymmetricGroupAlgebra(QQ, Integer(4)) >>> GP = SGA.garsia_procesi_module([Integer(2), Integer(2)]) >>> GP.one_basis() 5
- class sage.combinat.symmetric_group_representations.SpechtRepresentation(parent, partition)[source]¶
Bases:
SymmetricGroupRepresentation_generic_class
- representation_matrix(permutation)[source]¶
Return the matrix representing the
permutation
in this irreducible representation.Note
This method caches the results.
EXAMPLES:
sage: spc = SymmetricGroupRepresentation([3,1], 'specht') sage: spc.representation_matrix(Permutation([2,1,3,4])) [ 0 -1 0] [-1 0 0] [ 0 0 1] sage: spc.representation_matrix(Permutation([3,2,1,4])) [0 0 1] [0 1 0] [1 0 0]
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(3),Integer(1)], 'specht') >>> spc.representation_matrix(Permutation([Integer(2),Integer(1),Integer(3),Integer(4)])) [ 0 -1 0] [-1 0 0] [ 0 0 1] >>> spc.representation_matrix(Permutation([Integer(3),Integer(2),Integer(1),Integer(4)])) [0 0 1] [0 1 0] [1 0 0]
- scalar_product(u, v)[source]¶
Return
0
ifu+v
is not a permutation, and the signature of the permutation otherwise.This is the scalar product of a vertex
u
of the underlying Yang-Baxter graph with the vertexv
in the ‘dual’ Yang-Baxter graph.EXAMPLES:
sage: spc = SymmetricGroupRepresentation([3,2], 'specht') sage: spc.scalar_product((1,0,2,1,0),(0,3,0,3,0)) -1 sage: spc.scalar_product((1,0,2,1,0),(3,0,0,3,0)) 0
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(3),Integer(2)], 'specht') >>> spc.scalar_product((Integer(1),Integer(0),Integer(2),Integer(1),Integer(0)),(Integer(0),Integer(3),Integer(0),Integer(3),Integer(0))) -1 >>> spc.scalar_product((Integer(1),Integer(0),Integer(2),Integer(1),Integer(0)),(Integer(3),Integer(0),Integer(0),Integer(3),Integer(0))) 0
- scalar_product_matrix(permutation=None)[source]¶
Return the scalar product matrix corresponding to
permutation
.The entries are given by the scalar products of
u
andpermutation.action(v)
, whereu
is a vertex in the underlying Yang-Baxter graph andv
is a vertex in the dual graph.EXAMPLES:
sage: spc = SymmetricGroupRepresentation([3,1], 'specht') sage: spc.scalar_product_matrix() [ 1 0 0] [ 0 -1 0] [ 0 0 1]
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(3),Integer(1)], 'specht') >>> spc.scalar_product_matrix() [ 1 0 0] [ 0 -1 0] [ 0 0 1]
- class sage.combinat.symmetric_group_representations.SpechtRepresentations(n, ring=None, cache_matrices=True)[source]¶
Bases:
SymmetricGroupRepresentations_class
- Element[source]¶
alias of
SpechtRepresentation
- sage.combinat.symmetric_group_representations.SymmetricGroupRepresentation(partition, implementation='specht', ring=None, cache_matrices=True)[source]¶
The irreducible representation of the symmetric group corresponding to
partition
.INPUT:
partition
– a partition of a positive integerimplementation
– string (default:'specht'
); one of:'seminormal'
– for Young’s seminormal representation'orthogonal'
– for Young’s orthogonal representation'specht'
– for Specht’s representation'unitary'
– for the unitary representation
ring
– the ring over which the representation is definedcache_matrices
– boolean (default:True
); ifTrue
, then any representation matrices that are computed are cached
EXAMPLES:
Young’s orthogonal representation: the matrices are orthogonal.
sage: orth = SymmetricGroupRepresentation([2,1], "orthogonal"); orth # needs sage.symbolic Orthogonal representation of the symmetric group corresponding to [2, 1] sage: all(a*a.transpose() == a.parent().identity_matrix() for a in orth) # needs sage.symbolic True
>>> from sage.all import * >>> orth = SymmetricGroupRepresentation([Integer(2),Integer(1)], "orthogonal"); orth # needs sage.symbolic Orthogonal representation of the symmetric group corresponding to [2, 1] >>> all(a*a.transpose() == a.parent().identity_matrix() for a in orth) # needs sage.symbolic True
sage: # needs sage.symbolic sage: orth = SymmetricGroupRepresentation([3,2], "orthogonal"); orth Orthogonal representation of the symmetric group corresponding to [3, 2] sage: orth([2,1,3,4,5]) [ 1 0 0 0 0] [ 0 1 0 0 0] [ 0 0 -1 0 0] [ 0 0 0 1 0] [ 0 0 0 0 -1] sage: orth([1,3,2,4,5]) [ 1 0 0 0 0] [ 0 -1/2 1/2*sqrt(3) 0 0] [ 0 1/2*sqrt(3) 1/2 0 0] [ 0 0 0 -1/2 1/2*sqrt(3)] [ 0 0 0 1/2*sqrt(3) 1/2] sage: orth([1,2,4,3,5]) [ -1/3 2/3*sqrt(2) 0 0 0] [2/3*sqrt(2) 1/3 0 0 0] [ 0 0 1 0 0] [ 0 0 0 1 0] [ 0 0 0 0 -1]
>>> from sage.all import * >>> # needs sage.symbolic >>> orth = SymmetricGroupRepresentation([Integer(3),Integer(2)], "orthogonal"); orth Orthogonal representation of the symmetric group corresponding to [3, 2] >>> orth([Integer(2),Integer(1),Integer(3),Integer(4),Integer(5)]) [ 1 0 0 0 0] [ 0 1 0 0 0] [ 0 0 -1 0 0] [ 0 0 0 1 0] [ 0 0 0 0 -1] >>> orth([Integer(1),Integer(3),Integer(2),Integer(4),Integer(5)]) [ 1 0 0 0 0] [ 0 -1/2 1/2*sqrt(3) 0 0] [ 0 1/2*sqrt(3) 1/2 0 0] [ 0 0 0 -1/2 1/2*sqrt(3)] [ 0 0 0 1/2*sqrt(3) 1/2] >>> orth([Integer(1),Integer(2),Integer(4),Integer(3),Integer(5)]) [ -1/3 2/3*sqrt(2) 0 0 0] [2/3*sqrt(2) 1/3 0 0 0] [ 0 0 1 0 0] [ 0 0 0 1 0] [ 0 0 0 0 -1]
The Specht representation:
sage: spc = SymmetricGroupRepresentation([3,2], "specht") sage: spc.scalar_product_matrix(Permutation([1,2,3,4,5])) [ 1 0 0 0 0] [ 0 -1 0 0 0] [ 0 0 1 0 0] [ 0 0 0 1 0] [-1 0 0 0 -1] sage: spc.scalar_product_matrix(Permutation([5,4,3,2,1])) [ 1 -1 0 1 0] [ 0 0 1 0 -1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [-1 0 0 0 -1] sage: spc([5,4,3,2,1]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1] sage: spc.verify_representation() True
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(3),Integer(2)], "specht") >>> spc.scalar_product_matrix(Permutation([Integer(1),Integer(2),Integer(3),Integer(4),Integer(5)])) [ 1 0 0 0 0] [ 0 -1 0 0 0] [ 0 0 1 0 0] [ 0 0 0 1 0] [-1 0 0 0 -1] >>> spc.scalar_product_matrix(Permutation([Integer(5),Integer(4),Integer(3),Integer(2),Integer(1)])) [ 1 -1 0 1 0] [ 0 0 1 0 -1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [-1 0 0 0 -1] >>> spc([Integer(5),Integer(4),Integer(3),Integer(2),Integer(1)]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1] >>> spc.verify_representation() True
The unitary representation:
sage: unitary = SymmetricGroupRepresentation([3,1], "unitary"); unitary Unitary representation of the symmetric group corresponding to [3, 1] sage: unitary_GF49 = SymmetricGroupRepresentation([3,1], "unitary", ring=GF(7**2)); unitary_GF49 Unitary representation of the symmetric group corresponding to [3, 1] sage: unitary_GF49([2,1,3,4]) [6 0 0] [0 1 0] [0 0 1] sage: unitary_GF49([3,2,1,4]) [ 4 2*z2 + 3 0] [5*z2 + 5 3 0] [ 0 0 1] sage: unitary_GF49.verify_representation() True
>>> from sage.all import * >>> unitary = SymmetricGroupRepresentation([Integer(3),Integer(1)], "unitary"); unitary Unitary representation of the symmetric group corresponding to [3, 1] >>> unitary_GF49 = SymmetricGroupRepresentation([Integer(3),Integer(1)], "unitary", ring=GF(Integer(7)**Integer(2))); unitary_GF49 Unitary representation of the symmetric group corresponding to [3, 1] >>> unitary_GF49([Integer(2),Integer(1),Integer(3),Integer(4)]) [6 0 0] [0 1 0] [0 0 1] >>> unitary_GF49([Integer(3),Integer(2),Integer(1),Integer(4)]) [ 4 2*z2 + 3 0] [5*z2 + 5 3 0] [ 0 0 1] >>> unitary_GF49.verify_representation() True
By default, any representation matrices that are computed are cached:
sage: spc = SymmetricGroupRepresentation([3,2], "specht") sage: spc([5,4,3,2,1]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1] sage: spc._cache__representation_matrix {(([5, 4, 3, 2, 1],), ()): [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1]}
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(3),Integer(2)], "specht") >>> spc([Integer(5),Integer(4),Integer(3),Integer(2),Integer(1)]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1] >>> spc._cache__representation_matrix {(([5, 4, 3, 2, 1],), ()): [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1]}
This can be turned off with the keyword
cache_matrices
:sage: spc = SymmetricGroupRepresentation([3,2], "specht", cache_matrices=False) sage: spc([5,4,3,2,1]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1] sage: hasattr(spc, '_cache__representation_matrix') False
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(3),Integer(2)], "specht", cache_matrices=False) >>> spc([Integer(5),Integer(4),Integer(3),Integer(2),Integer(1)]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1] >>> hasattr(spc, '_cache__representation_matrix') False
Note
The implementation is based on the paper [Las].
REFERENCES:
[Las] (1,2)Alain Lascoux, ‘Young representations of the symmetric group.’ http://phalanstere.univ-mlv.fr/~al/ARTICLES/ProcCrac.ps.gz
AUTHORS:
Franco Saliola (2009-04-23)
- class sage.combinat.symmetric_group_representations.SymmetricGroupRepresentation_generic_class(parent, partition)[source]¶
Bases:
Element
Generic methods for a representation of the symmetric group.
- to_character()[source]¶
Return the character of the representation.
EXAMPLES:
The trivial character:
sage: rho = SymmetricGroupRepresentation([3]) sage: chi = rho.to_character(); chi Character of Symmetric group of order 3! as a permutation group sage: chi.values() [1, 1, 1] sage: all(chi(g) == 1 for g in SymmetricGroup(3)) True
>>> from sage.all import * >>> rho = SymmetricGroupRepresentation([Integer(3)]) >>> chi = rho.to_character(); chi Character of Symmetric group of order 3! as a permutation group >>> chi.values() [1, 1, 1] >>> all(chi(g) == Integer(1) for g in SymmetricGroup(Integer(3))) True
The sign character:
sage: rho = SymmetricGroupRepresentation([1,1,1]) sage: chi = rho.to_character(); chi Character of Symmetric group of order 3! as a permutation group sage: chi.values() [1, -1, 1] sage: all(chi(g) == g.sign() for g in SymmetricGroup(3)) True
>>> from sage.all import * >>> rho = SymmetricGroupRepresentation([Integer(1),Integer(1),Integer(1)]) >>> chi = rho.to_character(); chi Character of Symmetric group of order 3! as a permutation group >>> chi.values() [1, -1, 1] >>> all(chi(g) == g.sign() for g in SymmetricGroup(Integer(3))) True
The defining representation:
sage: triv = SymmetricGroupRepresentation([4]) sage: hook = SymmetricGroupRepresentation([3,1]) sage: def_rep = lambda p : triv(p).block_sum(hook(p)).trace() sage: list(map(def_rep, Permutations(4))) [4, 2, 2, 1, 1, 2, 2, 0, 1, 0, 0, 1, 1, 0, 2, 1, 0, 0, 0, 1, 1, 2, 0, 0] sage: [p.to_matrix().trace() for p in Permutations(4)] [4, 2, 2, 1, 1, 2, 2, 0, 1, 0, 0, 1, 1, 0, 2, 1, 0, 0, 0, 1, 1, 2, 0, 0]
>>> from sage.all import * >>> triv = SymmetricGroupRepresentation([Integer(4)]) >>> hook = SymmetricGroupRepresentation([Integer(3),Integer(1)]) >>> def_rep = lambda p : triv(p).block_sum(hook(p)).trace() >>> list(map(def_rep, Permutations(Integer(4)))) [4, 2, 2, 1, 1, 2, 2, 0, 1, 0, 0, 1, 1, 0, 2, 1, 0, 0, 0, 1, 1, 2, 0, 0] >>> [p.to_matrix().trace() for p in Permutations(Integer(4))] [4, 2, 2, 1, 1, 2, 2, 0, 1, 0, 0, 1, 1, 0, 2, 1, 0, 0, 0, 1, 1, 2, 0, 0]
- verify_representation()[source]¶
Verify the representation.
This tests that the images of the simple transpositions are involutions and tests that the braid relations hold.
EXAMPLES:
sage: spc = SymmetricGroupRepresentation([1,1,1]) sage: spc.verify_representation() True sage: spc = SymmetricGroupRepresentation([4,2,1]) sage: spc.verify_representation() True
>>> from sage.all import * >>> spc = SymmetricGroupRepresentation([Integer(1),Integer(1),Integer(1)]) >>> spc.verify_representation() True >>> spc = SymmetricGroupRepresentation([Integer(4),Integer(2),Integer(1)]) >>> spc.verify_representation() True
- sage.combinat.symmetric_group_representations.SymmetricGroupRepresentations(n, implementation='specht', ring=None, cache_matrices=True)[source]¶
Irreducible representations of the symmetric group.
INPUT:
n
– positive integerimplementation
– string (default:'specht'
); one of:'seminormal'
– for Young’s seminormal representation'orthogonal'
– for Young’s orthogonal representation'specht'
– for Specht’s representation'unitary'
– for the unitary representation
ring
– the ring over which the representation is definedcache_matrices
– boolean (default:True
); ifTrue
, then any representation matrices that are computed are cached
EXAMPLES:
Young’s orthogonal representation: the matrices are orthogonal.
sage: orth = SymmetricGroupRepresentations(3, "orthogonal"); orth # needs sage.symbolic Orthogonal representations of the symmetric group of order 3! over Symbolic Ring sage: orth.list() # needs sage.symbolic [Orthogonal representation of the symmetric group corresponding to [3], Orthogonal representation of the symmetric group corresponding to [2, 1], Orthogonal representation of the symmetric group corresponding to [1, 1, 1]] sage: orth([2,1])([1,2,3]) # needs sage.symbolic [1 0] [0 1]
>>> from sage.all import * >>> orth = SymmetricGroupRepresentations(Integer(3), "orthogonal"); orth # needs sage.symbolic Orthogonal representations of the symmetric group of order 3! over Symbolic Ring >>> orth.list() # needs sage.symbolic [Orthogonal representation of the symmetric group corresponding to [3], Orthogonal representation of the symmetric group corresponding to [2, 1], Orthogonal representation of the symmetric group corresponding to [1, 1, 1]] >>> orth([Integer(2),Integer(1)])([Integer(1),Integer(2),Integer(3)]) # needs sage.symbolic [1 0] [0 1]
Young’s seminormal representation.
sage: snorm = SymmetricGroupRepresentations(3, "seminormal"); snorm Seminormal representations of the symmetric group of order 3! over Rational Field sage: sgn = snorm([1,1,1]); sgn Seminormal representation of the symmetric group corresponding to [1, 1, 1] sage: list(map(sgn, Permutations(3))) [[1], [-1], [-1], [1], [1], [-1]]
>>> from sage.all import * >>> snorm = SymmetricGroupRepresentations(Integer(3), "seminormal"); snorm Seminormal representations of the symmetric group of order 3! over Rational Field >>> sgn = snorm([Integer(1),Integer(1),Integer(1)]); sgn Seminormal representation of the symmetric group corresponding to [1, 1, 1] >>> list(map(sgn, Permutations(Integer(3)))) [[1], [-1], [-1], [1], [1], [-1]]
The Specht Representation.
sage: spc = SymmetricGroupRepresentations(5, "specht"); spc Specht representations of the symmetric group of order 5! over Integer Ring sage: spc([3,2])([5,4,3,2,1]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1]
>>> from sage.all import * >>> spc = SymmetricGroupRepresentations(Integer(5), "specht"); spc Specht representations of the symmetric group of order 5! over Integer Ring >>> spc([Integer(3),Integer(2)])([Integer(5),Integer(4),Integer(3),Integer(2),Integer(1)]) [ 1 -1 0 1 0] [ 0 0 -1 0 1] [ 0 0 0 -1 1] [ 0 1 -1 -1 1] [ 0 1 0 -1 1]
The unitary representation.
sage: unitary = SymmetricGroupRepresentations(3, "unitary"); unitary Unitary representations of the symmetric group of order 3! over Symbolic Ring sage: unitary_GF49 = SymmetricGroupRepresentations(4, "unitary", ring=GF(7**2)); unitary_GF49 Unitary representations of the symmetric group of order 4! over Finite Field in z2 of size 7^2 sage: unitary_GF49([3,1])([2,1,3,4]) [6 0 0] [0 1 0] [0 0 1] sage: unitary_GF49([3,1])([3,2,1,4]) [ 4 2*z2 + 3 0] [5*z2 + 5 3 0] [ 0 0 1]
>>> from sage.all import * >>> unitary = SymmetricGroupRepresentations(Integer(3), "unitary"); unitary Unitary representations of the symmetric group of order 3! over Symbolic Ring >>> unitary_GF49 = SymmetricGroupRepresentations(Integer(4), "unitary", ring=GF(Integer(7)**Integer(2))); unitary_GF49 Unitary representations of the symmetric group of order 4! over Finite Field in z2 of size 7^2 >>> unitary_GF49([Integer(3),Integer(1)])([Integer(2),Integer(1),Integer(3),Integer(4)]) [6 0 0] [0 1 0] [0 0 1] >>> unitary_GF49([Integer(3),Integer(1)])([Integer(3),Integer(2),Integer(1),Integer(4)]) [ 4 2*z2 + 3 0] [5*z2 + 5 3 0] [ 0 0 1]
Note
The implementation is based on the paper [Las].
AUTHORS:
Franco Saliola (2009-04-23)
- class sage.combinat.symmetric_group_representations.SymmetricGroupRepresentations_class(n, ring=None, cache_matrices=True)[source]¶
Bases:
UniqueRepresentation
,Parent
Generic methods for the CombinatorialClass of irreducible representations of the symmetric group.
- class sage.combinat.symmetric_group_representations.UnitaryRepresentation(parent, partition)[source]¶
Bases:
SymmetricGroupRepresentation_generic_class
A unitary representation of the symmetric group.
In characteristic zero, this is the same as the orthogonal representation since they are defined over the reals. In positive characteristic, we are able to construct representations when the field is finite, has square order, and the characteristic does not divide the order of the group. In this case, we construct the representation by first constructing the orthogonal representation
and then taking the extended Cholesky decomposition of the unique solution to the equation for all in .- representation_matrix(permutation)[source]¶
Return the matrix representing the
permutation
in this irreducible representation.Note
This method caches the results.
EXAMPLES:
sage: unitary_specht = SymmetricGroupRepresentation([3,1], 'unitary', ring=GF(7**2)) sage: unitary_specht.representation_matrix(Permutation([2,1,3,4])) [6 0 0] [0 1 0] [0 0 1] sage: unitary_specht.representation_matrix(Permutation([3,2,1,4])) [ 4 2*z2 + 3 0] [5*z2 + 5 3 0] [ 0 0 1] sage: all(A * A.H == 1 for A in [unitary_specht.representation_matrix(g) ....: for g in Permutations(4)]) True
>>> from sage.all import * >>> unitary_specht = SymmetricGroupRepresentation([Integer(3),Integer(1)], 'unitary', ring=GF(Integer(7)**Integer(2))) >>> unitary_specht.representation_matrix(Permutation([Integer(2),Integer(1),Integer(3),Integer(4)])) [6 0 0] [0 1 0] [0 0 1] >>> unitary_specht.representation_matrix(Permutation([Integer(3),Integer(2),Integer(1),Integer(4)])) [ 4 2*z2 + 3 0] [5*z2 + 5 3 0] [ 0 0 1] >>> all(A * A.H == Integer(1) for A in [unitary_specht.representation_matrix(g) ... for g in Permutations(Integer(4))]) True
- class sage.combinat.symmetric_group_representations.UnitaryRepresentations(n, ring=None, cache_matrices=True)[source]¶
Bases:
SymmetricGroupRepresentations_class
- Element[source]¶
alias of
UnitaryRepresentation
- class sage.combinat.symmetric_group_representations.YoungRepresentation_Orthogonal(parent, partition)[source]¶
Bases:
YoungRepresentation_generic
- class sage.combinat.symmetric_group_representations.YoungRepresentation_Seminormal(parent, partition)[source]¶
Bases:
YoungRepresentation_generic
- class sage.combinat.symmetric_group_representations.YoungRepresentation_generic(parent, partition)[source]¶
Bases:
SymmetricGroupRepresentation_generic_class
Generic methods for Young’s representations of the symmetric group.
- representation_matrix(permutation)[source]¶
Return the matrix representing
permutation
.EXAMPLES:
sage: orth = SymmetricGroupRepresentation([2,1], "orthogonal") # needs sage.symbolic sage: orth.representation_matrix(Permutation([2,1,3])) # needs sage.symbolic [ 1 0] [ 0 -1] sage: orth.representation_matrix(Permutation([1,3,2])) # needs sage.symbolic [ -1/2 1/2*sqrt(3)] [1/2*sqrt(3) 1/2]
>>> from sage.all import * >>> orth = SymmetricGroupRepresentation([Integer(2),Integer(1)], "orthogonal") # needs sage.symbolic >>> orth.representation_matrix(Permutation([Integer(2),Integer(1),Integer(3)])) # needs sage.symbolic [ 1 0] [ 0 -1] >>> orth.representation_matrix(Permutation([Integer(1),Integer(3),Integer(2)])) # needs sage.symbolic [ -1/2 1/2*sqrt(3)] [1/2*sqrt(3) 1/2]
sage: norm = SymmetricGroupRepresentation([2,1], "seminormal") sage: p = PermutationGroupElement([2,1,3]) sage: norm.representation_matrix(p) [ 1 0] [ 0 -1] sage: p = PermutationGroupElement([1,3,2]) sage: norm.representation_matrix(p) [-1/2 3/2] [ 1/2 1/2]
>>> from sage.all import * >>> norm = SymmetricGroupRepresentation([Integer(2),Integer(1)], "seminormal") >>> p = PermutationGroupElement([Integer(2),Integer(1),Integer(3)]) >>> norm.representation_matrix(p) [ 1 0] [ 0 -1] >>> p = PermutationGroupElement([Integer(1),Integer(3),Integer(2)]) >>> norm.representation_matrix(p) [-1/2 3/2] [ 1/2 1/2]
- representation_matrix_for_simple_transposition(i)[source]¶
Return the matrix representing the transposition that swaps
i
andi+1
.EXAMPLES:
sage: orth = SymmetricGroupRepresentation([2,1], "orthogonal") # needs sage.symbolic sage: orth.representation_matrix_for_simple_transposition(1) # needs sage.symbolic [ 1 0] [ 0 -1] sage: orth.representation_matrix_for_simple_transposition(2) # needs sage.symbolic [ -1/2 1/2*sqrt(3)] [1/2*sqrt(3) 1/2] sage: norm = SymmetricGroupRepresentation([2,1], "seminormal") sage: norm.representation_matrix_for_simple_transposition(1) [ 1 0] [ 0 -1] sage: norm.representation_matrix_for_simple_transposition(2) [-1/2 3/2] [ 1/2 1/2]
>>> from sage.all import * >>> orth = SymmetricGroupRepresentation([Integer(2),Integer(1)], "orthogonal") # needs sage.symbolic >>> orth.representation_matrix_for_simple_transposition(Integer(1)) # needs sage.symbolic [ 1 0] [ 0 -1] >>> orth.representation_matrix_for_simple_transposition(Integer(2)) # needs sage.symbolic [ -1/2 1/2*sqrt(3)] [1/2*sqrt(3) 1/2] >>> norm = SymmetricGroupRepresentation([Integer(2),Integer(1)], "seminormal") >>> norm.representation_matrix_for_simple_transposition(Integer(1)) [ 1 0] [ 0 -1] >>> norm.representation_matrix_for_simple_transposition(Integer(2)) [-1/2 3/2] [ 1/2 1/2]
- class sage.combinat.symmetric_group_representations.YoungRepresentations_Orthogonal(n, ring=None, cache_matrices=True)[source]¶
Bases:
SymmetricGroupRepresentations_class
- Element[source]¶
alias of
YoungRepresentation_Orthogonal
- class sage.combinat.symmetric_group_representations.YoungRepresentations_Seminormal(n, ring=None, cache_matrices=True)[source]¶
Bases:
SymmetricGroupRepresentations_class
- Element[source]¶
alias of
YoungRepresentation_Seminormal
- sage.combinat.symmetric_group_representations.partition_to_vector_of_contents(partition, reverse=False)[source]¶
Return the “vector of contents” associated to
partition
.EXAMPLES:
sage: from sage.combinat.symmetric_group_representations import partition_to_vector_of_contents sage: partition_to_vector_of_contents([3,2]) (0, 1, 2, -1, 0)
>>> from sage.all import * >>> from sage.combinat.symmetric_group_representations import partition_to_vector_of_contents >>> partition_to_vector_of_contents([Integer(3),Integer(2)]) (0, 1, 2, -1, 0)