# OK to post homework
# Lucy Martinez, 03-28-2026, Assignment 17

with(combinat):


# Question 1:
# Write a procedure RandSYT(L)
#  that inputs a partition L and outputs a uniformly-at-random
#  standard Young tableau of shape L, by implementing
#  The Greene-Nijenhuis-Wilf proof of the Hook-Lenght formula
#  (described at the bottom of p.105 and top of p. 106 (pp. 2 and 3 of the .pdf),
#  the proof is nice but irrelevant to the algorithm) 
# Hints: 
#  1. If you have a set S, to get a uniformly-at-random member of S do S[rand(1..nops(S))()].
#  2. Write a procedure OneTrial, that inputs a partition L and
#     picks a random corner-cell, following the "kick the buck" algorithm of
#     the paper. Let L1 be the shape obtained by removing that corner,
#     then recursively find RandSYT(L1) and put "n" in the removed corner.
#     You may have to add a new row with only n in it 
RandSYT:=proc(L1) local box,c,n,L,T:
n:=add(box, box in L1):

if n=0 then
 return([]):
fi:

c:=OneTrial(L1):
L:=DeleteCell(L1,c):
T:=RandSYT(L):
T:=AddEntry(T,c,n):
T:
end:


OneTrial:=proc(L) local S,ra,H,i,j:
S:=Cells(L):
ra:=S[rand(1..nops(S))()]:

while nops(Hook(L, ra)) > 1 do
 H:=Hook(L, ra) minus {ra}:
 ra:=H[rand(1..nops(H))()]:
od:
ra:
end:


#DeleteCell(L,c): Given a partition L and a cell c,
# removes the cell c -- it can be any cell and if the cell
# is not a corner, then it will no longer return a SYT
DeleteCell:=proc(L,c) local i, L1:
i:=c[1]:
if not member(c, Cells(L)) then
 return(FAIL):
fi:
L1:=[op(1..i-1, L), L[i]-1, op(i+1..nops(L), L)]:

if L1[i]=0 then
 if i<nops(L1) then
  L1:=[op(1..i-1, L1), op(i+1..nops(L1), L1)]:
 else
  L1:=[op(1..i-1, L1)]:
 fi:
fi:
L1:
end:


#AddEntry(T,c,n): Given a Tableaux T, a cell c and an integer n,
# inserts the element "n" in the cell c of the Tableaux T
AddEntry:=proc(T,c,n) local i,j,T1:
i:=c[1]:
j:=c[2]:
T1:=[op(T)]:

if i>nops(T1) then 
 T1:=[op(T1),[n]]:
else
 T1[i]:=[op(T1[i]),n]:
fi:
T1:
end:

#Question 2:
# By running RandSYT([5$5]) 10000 times, estimate the probability that "2"
#  is in location [1,2] (rather than [2,1]). 
# ANSWER: Should be about 1/2 of the time.
#  ProbRSYT([5,5,5,5,5],10000) returned 4959/10000
ProbRSYT:=proc(L,K) local co,i,T:
co:=0:
for i from 1 to K do
 T:=RandSYT(L):
 if T[1,2]=2 then
  co:=co+1:
 fi:
od:
co/K:
end:


#########################From previous class:
# C17.txt, March 23, 2026

Help17:=proc():
 print(`Cells(L), Conj(L), Hook(L,c), HL(L,c), NuSYTc(L)`):
 print(`HLc(L,c), NuSYTcc(L)`):
end:

#S:=SYT([3,3,3]):
# nops(S)=42
#One way that we could get a random Young Tableaux is to call in Maple: 
# S[rand(1..nops(S))()];
#However, SYT(L) is constructing all Young Tableaux - not efficient to do this
#We will set up for next week so we can construct a way to generate random Young Tableaux
#
#Robinson-Frame-Thrall
#Suppose we have placed the following boxes:
# ***
# ***
# ***
#NuSYT(L)=n!/mul(nops(Hook(c)), c a cell in the shape L)
# Let L be of shape:
# ***
# ***
# ***
# Here, n=9 (9 cells)
# To count the hook of any cell: count the the number of 
#  stars in the same column and to the right of it (row)
# So now, the hook length for cell=(1,1) is 5 
#  the hook length for cell (1,2) is 4
#  the hook length for cell (2,2) is 3
# If we do it all, we get
#  9!/(5*4*3*4*3*2*3*2*1)=42


#Cells(L): the set of n (=sum(L)) cells [i,j] in the shape L
Cells:=proc(L) local k,i,j:
k:=nops(L):
{seq(seq([i,j], j=1..L[i]),i=1..k)}:
end:

#Def: Conjugate of a partition: each column (left to right) becomes a row
# from top to bottom
# Example: 
# Say we have the following diagram
# 331
# ***
# ***
# *
# Then the conjugate is 322
# ***
# **
# **

#Conj(L): Given a partition L outputs the conjugate of L
Conj:=proc(L) local k,L1,i1,C1,i:
option remember:

if L=[] then
 return([]):
fi:

k:=nops(L):
L1:=[seq(L[i]-1,i=1..k)]:
for i1 from 1 to nops(L1) while L1[i1]>0 do od:
i1:=i1-1:
L1:=[op(1..i1,L1)]:
C1:=Conj(L1):
[k,op(C1)]:
end:


#Hook(L,c): The set of the cells in the hook corresponding
# to the cell c=[i,j]
# i.e. the set of cells to the right and to the bottom of call
Hook:=proc(L,c) local k,i,j,C,i1,j1,i2:
k:=nops(L):
i:=c[1]:
j:=c[2]:
if not (i>=1 and i<=k) then
 return(FAIL):
fi:

if not (j>=1 and j<=L[i]) then
 return(FAIL):
fi:

C:={seq([i,j1], j1=j..L[i])}:

for i2 from i to k while j<=L[i2] do od:
i2:=i2-1:
C:=C union {seq([i1,j],i1=i..i2)}:
C:
end:


#HL(L,c): the hook-length of the cell c in the shape L
HL:=proc(L,c) local n,C,i:
nops(Hook(L,c)):
end:

#HLc(L,c): the hook-length of the cell c in the shape L
# done cleverly
HLc:=proc(L,c) local i,j,L1:
L1:=Conj(L):
i:=c[1]:
j:=c[2]:
L[i]-j+L1[j]-i+1:
end:


#NuSYTc(L): implementing the Frame-Robinson-Thrall hook lenght formula
NuSYTc:=proc(L) local n,C,i,c:
n:=add(L[i],i=1..nops(L)):
C:=Cells(L):
n!/mul(HL(L,c),c in C):
end:


#NuSYTcc(L): implementing the Frame-Robinson-Thrall hook lenght formula
# done cleverly
NuSYTcc:=proc(L) local n,C,i,c:
n:=add(L[i],i=1..nops(L)):
C:=Cells(L):
n!/mul(HLc(L,c),c in C):
end:


#######################################
Help16:=proc():
 print(`RS1(Y,i), RS1left(pi), RS(pi)`):
end:

#We will finish the RSK algorithm
# Example:
# Given the following partial STY 
# 4 5 7
# 6 8
# Want to insert 2:
# 2 5 7  and bumpee=4 so now we have to insert 4 next
# 6 8    
# Want to insert 4 (after second row):
# 2 5 7
# 4 8    and bumpee=6 so now we have to insert 6 next
# Want to insert 6 (after third row):
# 2 5 7
# 4 8
# 6
#The previous SYT is the output  

#Now given a permutation pi=31452
#Begin STY:
# 3
# insert 1 and bump 3
# 1
# 3
# insert 4
# 1
# 3 4
# insert 5
# 1
# 3 4 5
# insert 2
# 1 2
# 3 4 5 


#RS1(Y,i): Inputs a partial Young tableaux and another integer i NOT yet in Y
# and places it in the right place, by a bumping process
# it returns a tableaux with one more box followed by the row 
RS1:=proc(Y,i) local k, NewY,lucy,bumpee,i1,j:
if Y=[] then
 return([[i]],1):
fi:
k:=nops(Y):
lucy:=RS11(Y[1],i):
NewY[1]:=lucy[1]:
bumpee:=lucy[2]:

if bumpee=0 then 
 return([NewY[1], op(2..k,Y)],1):
fi:

for i1 from 2 to k while bumpee<>0 do
 lucy:=RS11(Y[i1],bumpee):
 NewY[i1]:=lucy[1]:
 bumpee:=lucy[2]:
 if bumpee=0 then
  return( [seq(NewY[j],j=1..i1), op(i1+1..k,Y)],i1):
 fi:
od:
[seq(NewY[j],j=1..k),[bumpee]],k+1:
end:



#RS1left(pi): inputs a permutation pi (of size nops(pi)) and outputs ONE
# STY of whatever shape (with n boxes)
RS1left:=proc(pi) local Y,i:
Y:=[]:
for i from 1 to nops(pi) do
 Y:=RS1(Y,pi[i])[1]:
od:
Y:
end:



#RS(pi): inputs a permutation pi of ({1,2,...,n:=nops(pi)) and outputs
# a pair of SYT of the SAME shape (with n boxes)
#The Robinson–Schensted algorithm
RS:=proc(pi) local Y1, i,Yr,eaea,p:
if pi=[] then
 return([[],[]]):
fi:
Y1:=[[pi[1]]]:
Yr:=[[1]]:

for i from 2 to nops(pi) do
 eaea:=RS1(Y1,pi[i]):
 Y1:=eaea[1]:
 p:=eaea[2]:
 if p<=nops(Yr) then
  Yr:=[op(1..p-1,Yr),[op(Yr[p]),i],op(p+1..nops(Yr),Yr)]:
 else
  Yr:=[op(Yr),[i]]:
 fi:
od:
[Y1,Yr]:
#Yr is usually what they call Q
end:



############################
# C14.txt, March 09, 2026

Help14:=proc():
 print(`NuSYT(L), SYTpairs(n), NuSYTpairs(n)`):
 print(`RS11(a,i)`):
end:


#NuSYT(L): The NUMBER of Standard Young tableaux of shape L
NuSYT:=proc(L) local i, k,S,L1:
option remember:
k:=nops(L):

if k=0 then
  RETURN(1):
fi:
S:=0:

#now we look for all the legal rows where the element n can be placed
for i from 1 to k-1 do
if L[i]>L[i+1] then
 L1:=[op(1..i-1,L),L[i]-1,op(i+1..k,L)]:
 S:=S+NuSYT(L1): 
fi:
od:

if L[k]>1 then
L1:=[op(1..k-1,L),L[k]-1]:

S:=S+NuSYT(L1):
else
L1:=[op(1..k-1,L)]:
S:=S+NuSYT(L1):
fi:
S:
end:


#Pairs of SYT of the same shape
# We are interested in pairs (Y1,Y2) such that
# shape of Y1 equals shape of Y2
# For n=2:
#  SYT([1,1])={[ [1], [2] ] } so there is only 1 of these STY of shape [1,1]
#
#  SYTpairs(2)= { [[[1, 2]], [[1, 2]]], [ [[1], [2]], [[1], [2]] ]}

#SYTpairs(n):
SYTpairs:=proc(n) local S,P,p,S1,s1,s2:
option remember:
P:=Par(n):
S:={}:
for p in P do
 S1:=SYT(p): #the set of Young Tableaux of the partition p (shape p)
 S:=S union {seq(seq([s1,s2], s1 in S1),s2 in S1)}:
od:
S:
end:

#NuSYTpairs(n):
NuSYTpairs:=proc(n) local S,P, p:
option remember:
P:=Par(n):
S:=0:
for p in P do
S:=S+NuSYT(p)^2:
od:
S:
end:



#RSK: Find a mapping from the set of permutations to the set of STYpairs(n)
#Input: permutation of {1,2,...,n}
#Output: a pair of SYTs of the SAME shape (with n boxes)
#Example: Given [1,4,5,9] then where can we place 10 so that
# it is still increasing? --> put it at the end
# What if I want to place the number 3 using RSK such that if I place it
#  somewhere other than at the end, we must bump the number that was in it
#  So, [1,4,5,9] becomes [1,3,5,9]
# 

#RS11(a,i): inputs an INCREASING list of positive integers a and
# another positive integer i, outputs a pair a1,j where 
# (usually a1 is of the same length as a) and i is placed where it belongs
# and j is the entry that it bumped, unless i is larger than all the members of a
# (i.e. larger than a[-1]) then a1 is [op(a), i], and j is 0
# For example: RS1([11,12,13],4)= [4,12,13], 11
#  where the second output is the element that got bumped out
RS11:=proc(a,i) local k,j:
k:=nops(a):
for j from 1 to k while a[j]<i do
od:
if j=k+1 then
 return([op(a),i],0):
fi:
[op(1..j-1,a),i,op(j+1..k,a)], a[j]:
end:



#Now that we have the helper function RS11(a,i), we do the following:
#
#Input a Young Tableaux (YT) with distinct entries and an integer i
# NOT yet in the YT, outputs a larger YT with ONE more box, obtained as follows:
#  Try to insert i in the first row if it fits w/o bumping it
#  otherwise put it in its rightful place and bump somebody
#  now the one that got bumped has to bump someone else
#  and place the bumpee in the second row


#Example of the algorithm above
#Y= [ [1,4,8], [2,7], [3] ]
#i=5

#148 bumpee is 8
#27
#3

#becomes 

#145
#278
#3

#######################################
# C13.txt, March 05, 2026

Help13:=proc():
 print(`ParN(n,k), PFG(L), SYT(L), PSYT(Y)`):
end:


#ParN(n,k): outputs the partitions of n for which the first part equal to k
ParN:=proc(n,k) local s,S,T:
S:=Par(n):
T:={}:
for s in S do
 if s[1]=k then
  T:=T union {s}:
 fi:
od:
T:
end:


#PFG(L): The Ferrers graph of the partition L
PFG:=proc(L) local i:
for i from 1 to nops(L) do
 lprint(1$L[i]):
od:
end:



#SYT(L): the set of standard young tableaux of shape L
SYT:=proc(L) local i,s1,k,n,S,L1,S1:
option remember:
k:=nops(L):
n:=add(L[i],i=1..k):

if k=0 then
 return({[]}):
fi:

S:={}:

#now we look for all the legal rows where the element n can be placed
for i from 1 to k-1 do
 if L[i]>L[i+1] then
  L1:=[op(1..i-1,L),L[i]-1,op(i+1..k,L)]:
  S1:=SYT(L1):
  S:=S union {seq([op(1..i-1,s1) , [op(s1[i]),n], op(i+1..k,s1)] ,s1 in S1)}:
 fi:
od:

if L[k]>1 then 
 L1:=[op(1..k-1,L),L[k]-1]:
 S1:=SYT(L1):
 S:=S union {seq([op(1..k-1,s1) , [op(s1[k]),n]],s1 in S1)}:
else
 L1:=[op(1..k-1,L)]:
 S1:=SYT(L1):
 S:=S union {seq([op(1..k-1,s1) , [n] ],s1 in S1)}:
fi:
S:
end:


#PSYT(Y): prints the SYT Y
PSYT:=proc(Y) local i:
for i from 1 to nops(Y) do
 lprint(op(Y[i])):
od:
end:



#Def: We are going to construct the Standard Young tableaux
# of shape L=[3,3,2], which is a partition of 8

# The following is the Ferrers graph for [3,3,2]
# 111
# 111
# 11

# Now, imagine each number "1" is an empty box. You fill each box
#  with the numbers {1,...,8} such that every row and every column
#  are increasing from left to right and top to bottom, respectively
# One way is the following:
# 123
# 456
# 78

# Example 2:
# Let L=[2,2] so n=4
# This is the shape:
# 11
# 11
#
# So the Young Tableaux are the following 2:
#
# 12
# 34
#
# AND
#
# 13
# 24

# Where can n (the tallest person) sit w/o disturbing shorter people?
#

# Below is a recap from last week:
#
# Recall that the number of partitions of n, denoted by p(n) satisfies
# sum(p(n)*q^n,n=0..infinity)=Prod(1/(1-q^k),k=1..infinity)
#
# G.F. of distinct partitions Prod((1+q^k),k=1..infinity)
# sum(d(n)*q^n,n=0..infinity)=(1+q)*(1+q^2)*(1+q^3)*...
#  where d(n) is the number of partitions with distinct parts
#
# sum(o(n)*q^n,n=0..infinity)=Prod(1/(1-q^(2*k-1)),k=1..infinity)=
#  1/((1-q)*(1-q^3)*(1-q^5)*...)
#  where o(n) is the number of partitions where all the parts are odd 
#
# Recall that (1+z)=(1-z^2)/(1-z)
# SO,
# 1/((1-q)*(1-q^3)*(1-q^5)*...) = (1-q^2)*(1-q^4)*(1-q^6)*(1-q^8)*.../((1-q)*(1-q^2)*(1-q^3)*(1-q^4)...)
#                           = ( (1-q^2)/(1-q) )*( (1-q^4)/(1-q^2) )*( (1-q^6)/(1-q^3) )*...
#                           =  (1+q)*(1+q^2)*(1+q^3)*...
# On the left hand side we have the G.F. of the partitions with odd parts
# On the right hand side we have the G.F. of the partitions with distinct parts
# For more detail, see: https://sites.math.rutgers.edu/~zeilberg/numtheory/L22.pdf




###################################
# C12.txt, March 02, 2026

Help12:=proc():
 print(`Park(n,k), Par(n)`):
end:

#Park(n,k): The set of partitions of n into exactly k parts
Park:=proc(n,k) local S,k1,S1,s1:
option remember:

if n<k then
   RETURN({}):
fi:

if k=0 then
 if n=0 then
   RETURN({[]}):
   else
    RETURN({}):
fi:
fi:

S:={}:
for k1 from 0 to k do
 S1:=Park(n-k,k1):
 S:= S union {seq([op(s1+[1$k1]),1$(k-k1)], s1 in S1)}:
od:
S:
end:


#Par(n): The set of partitions of n, the clever way.
Par:=proc(n) local k:
{seq(op(Park(n,k)), k=1..n)}:
end: