#####################################################################
## SYT.txt Save this file as   SYT.txt to use it,                    #
# stay in the                                                        #
## same directory, get into Maple (by typing: maple <Enter> )        #
## and then type:  read `SYT.txt` <Enter>                            #
## Then follow the instructions given there                          #
##                                                                   #
## Written by Doron Zeilberger, Rutgers University ,                 #
##  DoronZeil at gmail dot com                                       # 
######################################################################



print(`Written: Dec. 2022- March 2023: tested for Maple 2020 `):
print(`This Version  : March 10, 2023  `):
print():
print(`This is  SYT.txt, a Maple package`):
print(`accompanying Shalosh B. Ekhad and Doron Zeilberger's  article: `):
print(` Experimenting with Standard Young Tableaux `):

print():
print(`The most current version is available on WWW at:`):
print(` http://sites.math.rutgers.edu/~zeilberg/tokhniot/SYT.txt .`):
print(`Please report all bugs to: DoronZeil at gmail dot com .`):

print():
print(`For general help, and a list of the MAIN functions,`):
print(` type "ezra();". For specific help type "ezra(procedure_name);" `):
print(`For a list of the supporting functions type: ezra1();`):
print(`For a list of the Checking functions type: ezraC();`):
print(`For a list of the Dyck functions (i.e. two-rowed Young tableaux) type: ezraD();`):
print(`For a list of the Story functions type: ezraSt();`):
print(`For a list of the SIMULTATION functions type: ezraSi();`):
print():







ezraSi:=proc()
if args=NULL then

print(`The Simulation procedures are: GNW, SiOcGF, SiPr `):

else
ezra(args):

fi:

end:

ezraSt:=proc()
if args=NULL then

print(`The Story procedures are: Paper1, Paper1L, Paper2, Paper2L, Paper3, Paper4, Paper5, Paper6  `):

else
ezra(args):

fi:

end:

ezraC:=proc()
if args=NULL then

print(`The Checking procedures are: CheckPrS, CheckPrSa `):

else
ezra(args):

fi:


end:


ezraD:=proc()
if args=NULL then

print(`The Dyck procedures are: A1nij, Anij, FindZero2, LimMomsD, MinPr2, OcCsD ` ):

else
ezra(args):

fi:


end:

ezra1:=proc()
if args=NULL then

print(`The SUPPORTING procedures are: DrTab, Igor,  Imags, NuSYTc, NuSYT, OcCs, OcGF,   MinPr, Pars2, Pars3, Pr, SSYT, NuSSYT, Shaps, Swee, Ti, YF  `):

else
ezra(args):

fi:


end:

ezra:=proc()
if args=NULL then

print(` SYT.txt: A Maple package for  Experimenting with Standard Young Tableaux`):

print(`The MAIN procedures are:  FindZero,MinPrk, OcGFs, OcGFs1, OcGFs1L,  PrS, PrS1,SYT `):
print(``):


elif nargs=1 and args[1]=Anij then
print(`Anij(n,i,j): Same as PrS([n,n],i,[2,j]). Try:`):
print(`Anij(n,2,1);`):

elif nargs=1 and args[1]=A1nij then
print(`A1nij(n,i,j): Probability of Standard Young tableau of shape [n,n] such that the occupant of [1,i] is larger than the occupant of [2,j]. Try:`):
print(`A1nij(n,3,2);`):

elif nargs=1 and args[1]=CheckPrS then
print(`CheckPrS(S,n,K): Checks procedure PrS(S,i,c2) where S is a symbolic shape depending on the variable n, for i from 2 to K and all c2 possible. Try:`):
print(`CheckPrS([n,n],n,6);`):

elif nargs=1 and args[1]=CheckPrSa then
print(`CheckPrSa(S,n,i,c2,K): Checks procedure PrS(S,i,c2) for n from i to K, where S is a symbolic shape in terms of the symbol n. Try:`):
print(`CheckPrSa([n,n],1,[2,1],9);`):

elif nargs=1 and args[1]=DrTab then
print(`DrTab(T): draws a tableau T. Try:`):
print(`DrTab([[1,2],[3]]);`):



elif nargs=1 and args[1]=FindZero then
print(`FindZero(L,n,K): Given a symbolic shape L pharased in terms of a single variable n, and a positive integer K`):
print(`finds all pairs [i,c] where i<=K, and c is a cell with 2<=c[1]<=nops(L) and c[2]<i such that PrS(L,i,c) goes to zero as n goes to infinity`):
print(`FindZero([n,n],n,K) is the same as FindZero2(K)`):
print(`FindZero([n,n,n],n,3);`):

elif nargs=1 and args[1]=FindZero2 then
print(`FindZero2(K): Finds all pairs [[1,i],[2,j]], 2<=j<i<=K such that the limit of Anij(n,i,j) as n goes to infinity  is 0. Try:`):
print(`FindZero2(10);`):

elif nargs=1 and args[1]=GNW then
print(`GNW(L): A random Standard Young tableau of shape L`):
print(`using the Greene-Nijenhus-Wilf algorithm. For example do:`):
print(`GNW([3,3,3]);`):

elif nargs=1 and args[1]=Igor then
print(`Igor(L,M): NuSSYTc(L,M)/MNC(L). Try:`):
print(`Igor([n,n,n],[2,1]);`):

elif nargs=1 and args[1]=Imags then
print(`Imags(x): All the images of the point x under the Ti's together with their sign. try:`):
print(`Imags([3,2,1]);`):


elif nargs=1 and args[1]=LimMomsD then
print(`LimMomsD(i,K1,K2): inputs a symbol i and  positive integers K1, K2, outputs a list of length K1 whose first entry is the limiting expectation (as n goes to infinity) of`):
print(`the occupant of cell [1,i] in a standard Young tableau of shape [n,n], followed by the variance, and scaled third-through-K th moment. It also returns the asympotics expressions to order K2`):
print(`and the scaling limits as i goes to infinity, and d the limiting scaled moments up to the K1, and the floating-point version, and finally the asymptotic to order K2 of the scaled moments. Try:`):
print(`LimMomsD(i,4,6);`):

elif nargs=1 and args[1]=MinPr then
print(`MinPr(L): The minimum of the absolute value of Pr(L,c1,c2) over all cells c1 and c2, followed by the set of champions. `):
print(`WARNING. Since it uses Pr it is VERY SLOW. Try:`):
print(`MinPr([4,4]);`):

elif nargs=1 and args[1]=MinPr2 then
print(`MinPr2(n): For numeric (pos. integer) n, Same as MinPr([n,n]) but much faster, using Anij(n,i,j). Try: MinPr2(30);`):


elif nargs=1 and args[1]=MinPrk then
print(`MinPrk(k,n1): For numeric (pos. integer) n1. Similar to  MinPr([n$k])]) but c1 must be in the first row. Hopefully faster (but not as fast as MinPr2(n1) when k=2). Try: `):
print(`MinPrk(3,10);`):

elif nargs=1 and args[1]=NuSYT then
print(`NuSYT(L): The NUMBER of all the standard Young Tableaux of shape L. Using dynamical programming. For checking purposes only`):
print(`For example, try:`):
print(` NuSYT([2,2]);`):


elif nargs=1 and args[1]=NuSSYT then
print(`NuSSYT(L,M): the NUMBER of all  Skew standard Young Tableaux of shape L/M, using Dynamical programming. For checking only.`):
print(`For example, try: NuSSYT([2,2],[1]);`):

elif nargs=1 and args[1]=NuSSYTc then
print(`NuSSYTc(L,M): the NUMBER of all  Skew standard Young Tableaux of shape L/M`):
print(`using The formula. Try:`):
print(`For example, try: NuSSYTc([2,2],[1,1]);`):


elif nargs=1 and args[1]=OcCs then
print(`OcCs(L,k,c): Given a symbolic shape L, a numeric integer k and a numeric cell c, finds the rational function in the symbols used in L`):
print(`that describes the exact probability that a random Young tableau of shape L, T, will have  T[c[1]][c[2]]=k. Try:`):
print(`OcCs([n,n],3,[1,2]);`):

elif nargs=1 and args[1]=OcCsD then
print(`OcCsD(n,k,i): Same as OcCs([n,n],k,[1,i]) but faster. Try:`):
print(`OcCsD(n,3,2);`):

elif nargs=1 and args[1]=OcGF then
print(`OcGF(L,c,x): Given a shape L and a cell in it,and a variable x, outputs the prob. generating function, in the variable x, for the occupant of cell c in a random Young tableau of shape L.`):
print(`USES BRUTE FORCE. FOR CHECKING PURPOSES ONLY. Try: `):
print(`OcGF([4,4,4],[2,2],x);`):


elif nargs=1 and args[1]=OcGFs then
print(`OcGFs(L,c,x,K): Given a symbolic shape L and a numeric, cell ,and a variable x, outputs the prob. generating function, in the variable x, for the occupant of cell c in a random Young tableau of shape L`):
print(`conditioned that it is <=K. It also gives you the probability of it being <=K.  Try:`):
print(`OcGFs([n,n],[1,5],x,9);`):

elif nargs=1 and args[1]=OcGFs1 then
print(`OcGFs1(L,i,x): Given a symbolic shape L and a numeric positive integer i, and a variable x, outputs the prob. generating function, in the variable x, for the occupant of cell [1,i] in a random Young tableau of shape L`):
print(`It is the same as  OcGFs(L,[1,i],x,nops(L)*(i-1)+1)[1]: Try:`):
print(`OcGFs1([n,n,n,n],2,x);`):


elif nargs=1 and args[1]=OcGFs1L then
print(`OcGFs1L(L,n,i,x,K): Given a symbolic shape L phrased in terms of the  variable n, outputs the limiting distribution as n goes to infinity `):
print(`of the occupants of the cell [1,i] in a random Young tableau of shape L`):
print(`followed by a list consisting of the average, s.d, and the scaled moments through the K-th.`):
print(`Try:`):
print(`OcGFs1L([n,n],n,10,x,4);`):



elif nargs=1 and args[1]=Paper1 then
print(`Paper1(R,n,K1,K2): inputs a positive integer R>=2, a symbol, n, and integers K1 and K2, outputs a paper with explicit expressions in n`):
print(`for the probability distribution of "occupant of cell [1,i]" for a random Young Tableau of shape [n$R], for all i from 2 to K1. It also gives`):
print(`the limiting distribution as n goes to infinity, and average, s.d. and the (scaled) moments up to the K2-th. Try:`):
print(`Paper1(2,n,10,4);`):


elif nargs=1 and args[1]=Paper1L then
print(`Paper1L(R,K1,K2): An abbreviated version of Paper1(R,n,K1,K2), only giving the limits as n goes to infinity`):
print(`Paper1L(2,10,4);`):


elif nargs=1 and args[1]=Paper2 then
print(`Paper2(R,n,K1): inputs a positive integer R>=2, a symbol, n, and integer K1  outputs a paper with explicit expressions in n`):
print(`for the Pr(T[1,i]>T[k,j])-Pr(T[k,j]>T[1,i]) (alias 2*Pr(T[1,i]>T[k,j])-1) for a random Young tableau of shape [n$R]`):
print(`for all i from 2 to K1, k from 2 to R, and j from 1 to i-1.`):
print(` It also gives the limiting quantity and the "cut-off" place where the limit switches from being pos. to negative. Try:`):
print(`Paper2(2,n,10);`):


elif nargs=1 and args[1]=Paper2L then
print(`Paper2L(R,K1): An abbreviated version of Paper2(n,R,K1), only stating limiting values`):
print(`for the limit of Pr(T[1,i]>T[k,j])-Pr(T[k,j]>T[1,i]) (alias 2*Pr(T[1,i]>T[k,j])-1) for a random Young tableau of shape [n$R]`):
print(`as n goes to infinity. Try:`):
print(`Paper2L(2,10);`):


elif nargs=1 and args[1]=Paper3 then
print(`Paper3(S,n,n0,i,K): Given a symbolic shape S, phrased in terms of the variable n, a positive integer n0, a positive integer i, and a large positive integer K`):
print(`illustates both the efficiently of symbolic computation, and that of simulation using the beautiful Greene-Nijenhuis-Wilf algorithm to generate`):
print(`a random Standard Young Tableau of, by comparing the exact probability distribution of the occupants of the cell [1,i], and what you get`):
print(`by sampling K random standard Young tableaux using the GNW algorithm. The former uses procedure OcGFs1 (q.v) and the later procedure SiOcGF. Try:`):
print(`Paper3([3*n,2*n],n,10,5,100);`):


elif nargs=1 and args[1]=Paper4 then
print(`Paper4(S,n,n0,i,c2,K): Given a symbolic shape S, phrased in terms of the variable n, a positive integer n0, a positive integer i, and a numeric cell c2`):
print(`and a large integer K`):
print(`illustates both the efficiently of symbolic computation, and that of simulation using the beautiful Greene-Nijenhuis-Wilf algorithm to `):
print(`a random Standard Young Tableau of, by comparing the sorting probability of cell [1,i] vs. cell c2, and what you get`):
print(`by sampling K random standard Young tableaux using the GNW algorithm. The former uses procedure PrS (q.v) and the later proceaure SiPr (q.v.). Try:`):
print(`Paper4([3*n,2*n],n,10,2,[2,1],100);`):

elif nargs=1 and args[1]=Paper5 then
print(`Paper5(k,N): inputs a positive integer k>1 and a positive integer N outputs the minimal sorting probabilities of the shape [n$k]  for n from 1 to N.`):
print(`If k>2 it is only an upper bound, since it only takes into account comparisons with the cells of the first row. Try:`):
print(`Paper5(3,10):`):

elif nargs=1 and args[1]=Paper6 then
print(`Paper6(i,K1,K2): A verbose form of LimMomsD(i,K1,K2) (q.v.). Try:`):
print(`Paper6(i,10,10):`):

elif nargs=1 and args[1]=Pars2 then
print(`Pars2(n,s): all the partitions of n with at most s parts. Try:`):
print(`Pars2(6,2);`):

elif nargs=1 and args[1]=Pars3 then
print(`Pars3(n,s,c): all the partitions of n with at most s parts where cell c2 is a corner cell. Try:`):
print(`Pars3(7,3,[2,2]);`):

elif nargs=1 and args[1]=Pr then
print(`Pr(L,c1,c2):Given a shape L and two cells c1 and c2 finds Pr(T[c1]>T[c2])-Pr(T[c1]<T[c2]) (alias 2*Pr(T[c1]>T[c2])-1). `):
print(` BY PURE BRUTE FORCE, for CHECKING PURPOSES ONLY. Try:`):
print(`Pr([3,3,3],[1,2],[2,1]);`):



elif nargs=1 and args[1]=PrS1 then
print(`PrS1(S,i,k,c2): Inputs a symbolic shape,  numeric pos. integers i and k, and a pair numeric c2 indicating a location such that c2[1]<=R and c2[2]<i, outputs`):
print(`the rational function, in the variables featuring in S,  describing the probability that a random standard Young tableau of the shape S would have cell [1,i] occupied by k`):
print(`and the occupant of cell c2 be smaller than k. Try:`):
print(`PrS1([n,n,n],2,3,[2,1]);`):


elif nargs=1 and args[1]=PrS then
print(`PrS(S,i,c2): Inputs symbolic shape S, a  numeric pos. integer i, and a numeric pair c2 indicating a location such that c2[1]<=nops(S) and c2[2]<i, outputs`):
print(`a closed form expression as a rational function in the symbols featuring in S of Pr(S,[1,i],c2); for any specialization of S. Try:`):
print(`PrS([n,n,n],2,[2,1]);`):

elif nargs=1 and args[1]=Shaps then
print(`Shaps(n,r,c2,s): all the shapes of size n with largest part  r and at most s parts, that contain the cell c2. Try:`):
print(`Shaps(10,5,[3,2],s);`):



elif nargs=1 and args[1]=SiPr then
print(`SiPr(L,c1,c2,K):Same as Pr(L,c1,c2) but approximately,simulating with K tries, using the Greene-Nijenhuis-Wilf algoirthm. Try:`):
print(`SiPr([3,3,3],[1,2],[2,1],1000);`):

elif nargs=1 and args[1]=SiOcGF then
print(`SiOcGF(L,c,x,K): Given a shape L and a cell in it,and a variable x, outputs the APPROXIMATE prob. generating function, in the variable x, for the occupant of cell c in a random Young tableau of shape L. `):
print(`using simulation (via the Greene-Nijenhuis-Wilf algorithm), with K tries. Just for fun! Try:`):
print(`SiOcGF([4,4,4],[2,2],x,1000);`):



elif nargs=1 and args[1]=Swee then
print(`Swee(L,M): Igor(L,M)/Igor(L). Try:`):
print(`Swee([n,n,n],[2,1]);`):

elif nargs=1 and args[1]=SYT then
print(`SYT(L): The set of all the standard Young Tableaux of shape L`):
print(`For example, try:`):
print(` SYT([2,2]);`):


elif nargs=1 and args[1]=SSYT then
print(`SSYT(L,M): the set of all  Skew standard Young Tableaux of shape L/M`):
print(`For example, try: SSYT([2,2],[1,1]);`):

elif nargs=1 and args[1]=Ti then
print(`Ti(x,i): Given a point x :makes x[i], x[i+1]-1 and x[i+1] x[i]+1. Try:`):
print(`Ti([a,b,c],2);`):

elif nargs=1 and args[1]=YF then
print(`YF(L): The Young-Frobenius formula. Try:`):
print(`YF([5,3,2]);`):

else
 print(`There is no such thing as`, args):

fi:


end:


with(combinat):


###START FROM GreeneNiejenhuisWilf

#OneStepGNW1(L,cell): Given a shape L and a cell=[i,j]
#decides where in the hook to go
OneStepGNW1:=proc(L,cell) local H,i,j,i1,j1:

i:=cell[1]:
j:=cell[2]:

if not (1<=i and i<=nops(L) and 1<=j and j<=L[i]) then
 ERROR(`Bad input`):
fi:

H:=[seq([i,j1],j1=j+1..L[i])]:

for i1 from i+1 to nops(L) while j<=L[i1] do
 H:=[op(H),[i1,j]]:
od:

if H=[] then
 RETURN(FAIL):
fi:

H[rand(1..nops(H))()]:
end:




#OneStepGNW(L): Given a shape L applies one step
#of the Greene-Nijenhius-Wilf algorithm to decide
#where to put n. Try:
#OneStepGNW([2,2]);
OneStepGNW:=proc(L) local cell,n,cell1,cell2,T1,i1,j1:
n:=convert(L,`+`):
T1:=[seq(seq([i1,j1],j1=1..L[i1]),i1=1..nops(L))]:
if n=1 then
 RETURN([1,1]):
fi:


cell:=T1[rand(1..n)()]:

cell1:=OneStepGNW1(L,cell):

if cell1=FAIL then
 RETURN(cell):
fi:

while cell1<>FAIL do
cell2:=OneStepGNW1(L,cell1):
cell:=cell1:
cell1:=cell2:
od:
cell:

end:


#GNW(L): A random Standard Young tableau of shape L
#using the Greene-Nijenhus-Wilf algorithm. For example do:
#GNW([3,3,3]);
GNW:=proc(L) local L1,cell,T1,n,i:

n:=convert(L,`+`):
if n=1 then
RETURN([[1]]):
fi:

cell:=OneStepGNW(L):
i:=cell[1]:

if L[i]>1 then
L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:
else
L1:=[op(1..i-1,L),op(i+1..nops(L),L)]:
fi:


T1:=GNW(L1):

if L[i]>1 then
[op(1..i-1,T1),[op(T1[i]),n],op(i+1..nops(T1),T1)]:
else
[op(1..i-1,T1),[n]]:
fi:

end:

#SYT(L): the set of all  standard Young Tableaux of shape L
#For example, try: SYT([2,2]);
SYT:=proc(L) local gu,n,L1,mu,mu1,i:
option remember:

if L=[] then
 RETURN({[]}):
fi:

n:=convert(L,`+`):

gu:={}:

for i from 1 to nops(L) do

if i=nops(L) or L[i]>L[i+1] then
 L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:
 if L1[nops(L1)]=0 then
  L1:=[op(1..nops(L1)-1,L1)]:
 fi:
mu:=SYT(L1):

if nops(L1)<nops(L) then
 gu:=gu union {seq([op(mu1),[n]],mu1 in mu)}:
else
 gu:=gu union 
{
seq(
[op(1..i-1,mu1),
[op(mu1[i]),n],
op(i+1..nops(mu1),mu1)],
mu1 in mu)
}:
fi:

fi:
od:

gu:

end:

##end from Greene-Nijenhuis-Wilf





#SYT(L): the set of all  standard Young Tableaux of shape L
#For example, try: SYT([2,2]);
SYT:=proc(L) local gu,n,L1,mu,mu1,i:
option remember:

if L=[] then
 RETURN({[]}):
fi:

n:=convert(L,`+`):

gu:={}:

for i from 1 to nops(L) do

if i=nops(L) or L[i]>L[i+1] then
 L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:
 if L1[nops(L1)]=0 then
  L1:=[op(1..nops(L1)-1,L1)]:
 fi:
mu:=SYT(L1):

if nops(L1)<nops(L) then
 gu:=gu union {seq([op(mu1),[n]],mu1 in mu)}:
else
 gu:=gu union 
{
seq(
[op(1..i-1,mu1),
[op(mu1[i]),n],
op(i+1..nops(mu1),mu1)],
mu1 in mu)
}:
fi:

fi:
od:

gu:

end:




#NuSYT(L): the Number of all  standard Young Tableaux of shape L, by Dynamical Programming. Try
#For example, try: NuSYT([2,2]);
NuSYT:=proc(L) local gu,n,L1,i:
option remember:

if L=[] then
 RETURN(1):
fi:

n:=convert(L,`+`):

gu:=0:

for i from 1 to nops(L) do

if i=nops(L) or L[i]>L[i+1] then
 L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:

 if L1[nops(L1)]=0 then
  L1:=[op(1..nops(L1)-1,L1)]:

 fi:
gu:=gu+NuSYT(L1):
fi:




od:
gu:
end:







#YF(L): The Young-Frobenius formula
YF:=proc(L) local i,j,n,k:
n:=convert(L,`+`):
k:=nops(L):
n!/mul((L[i]+k-i)!,i=1..k)*mul(mul((L[j]+k-j)-(L[i]+k-i),j=1..i-1),i=1..k):
end:


#Pr(L,c1,c2):Given a shape L and two cells c1 and c2 finds Pr(T[c1]>T[c2])-Pr(T[c1]<T[c2]) (alias 2*Pr(T[c1]>T[c2])-1). Try:
#Pr([3,3,3],[1,2],[2,1]);
Pr:=proc(L,c1,c2) local gu,co,gu1:

if not (c1[1]<=nops(L) and c2[1]<=nops(L) and c1[2]<=L[c1[1]] and c2[2]<=L[c2[1]]) then
 print(`Bad input`):
 RETURN(FAIL):
fi:


gu:=SYT(L):
co:=0:

for gu1 in gu do

if gu1[c1[1]][c1[2]]>gu1[c2[1]][c2[2]] then
 co:=co+1:
fi:

od:

2*co/YF(L)-1:
end:



#MinPr(L): The minimum of the absolute value of Pr(L,c1,c2) over all cells c1 and c2, followed by the set of champions. Try
#MinPr([4,4]);
MinPr:=proc(L) local c1,c2,i,j,i1,j1,rec,cha,hal:

rec:=1:

cha:={}:

for i from 1 to nops(L) do
for j from 1 to L[i] do
 c1:=[i,j]:
 for i1 from i+1 to nops(L) do
  for j1 from 1 to min(j-1,L[i1]) do
    c2:=[i1,j1]:

     hal:=abs(Pr(L,c1,c2)):
   
   if   hal<rec then
      rec:=hal:
      cha:={[c1,c2]}:
  elif   hal=rec then
        cha:= cha union {[c1,c2]}:
  fi:
 od:
od:
od:
od:
rec,cha:
end:


#DrTab(T): draws a tableau T. Try:
#DrTab([[1,2],[3]]);
DrTab:=proc(T) local i:
for i from 1 to nops(T) do
lprint(op(T[i])):
od:
end:

#SSYT(L,M): the set of all  Skew standard Young Tableaux of shape L/M
#For example, try: SSYT([2,2],[1,1]);
SSYT:=proc(L,M) local gu,n,mu,mu1,i,L1:
option remember:

if not (nops(M)<=nops(L) and   min(seq(L[i]-M[i],i=1..nops(M)))>=0) then
 RETURN({}):
fi:

if L=[] and M=[] then
 RETURN({[]}):
fi:


if L=M then
 RETURN({[seq([0$M[i]],i=1..nops(M))]}):
fi:


n:=convert(L,`+`)-convert(M,`+`):

gu:={}:

for i from 1 to nops(L) do

if i=nops(L) or L[i]>L[i+1] then
 L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:
 if L1[nops(L1)]=0 then
  L1:=[op(1..nops(L1)-1,L1)]:
 fi:
mu:=SSYT(L1,M):



if nops(L1)<nops(L) then
 gu:=gu union {seq([op(mu1),[n]],mu1 in mu)}:

else
 gu:=gu union 
{
seq(
[op(1..i-1,mu1),
[op(mu1[i]),n],
op(i+1..nops(mu1),mu1)],
mu1 in mu)
}:

fi:

fi:
od:

gu:

end:




#NuSSYT(L,M): the NUMBER of all  Skew standard Young Tableaux of shape L/M
#using Dynamical programming.
#For example, try: NuSSYT([2,2],[1,1]);
NuSSYT:=proc(L,M) local gu,n,i,L1:
option remember:

if not (nops(M)<=nops(L) and   min(seq(L[i]-M[i],i=1..nops(M)))>=0) then
 RETURN(0):
fi:

if L=[] and M=[] then
 RETURN(1):
fi:


if L=M then
 RETURN(1):
fi:


n:=convert(L,`+`)-convert(M,`+`):

gu:=0:

for i from 1 to nops(L) do

if i=nops(L) or L[i]>L[i+1] then
 L1:=[op(1..i-1,L),L[i]-1,op(i+1..nops(L),L)]:
 if L1[nops(L1)]=0 then
  L1:=[op(1..nops(L1)-1,L1)]:
 fi:
gu:=gu+NuSSYT(L1,M):
fi:

od:

gu:

end:


#Ti(x,i): Given a point x :makes x[i], x[i+1]-1 and x[i+1] x[i]+1. Try:
#Ti([a,b,c],2);
Ti:=proc(x,i) :
[op(1..i-1,x),x[i+1]-1,x[i]+1,op(i+2..nops(x),x)]:
end:



#Imags(x): All the images of the point x under the Ti's together with their sign. try:
#Imags([3,2,1]);
Imags:=proc(x) local gu,i,gu1,n,mu,lu:
n:=nops(x):
gu:={[x,1]}:
mu:=gu union {seq(seq([Ti(gu1[1],i),-gu1[2]],gu1 in gu),i=1..n-1)}:

while gu<>mu do
lu:=mu union {seq(seq([Ti(gu1[1],i),-gu1[2]],gu1 in mu),i=1..n-1)}:
gu:=mu:
mu:=lu:


od:
mu:
end:


MNC:=proc(L) local i: 
if {seq(type(L[i],integer),i=1..nops(L))}={true} and  min(op(L))<0 then
 0:
else
convert(L,`+`)!/mul(L[i]!,i=1..nops(L))
fi:
end:

#NuSSYTc(L,M): the NUMBER of all  Skew standard Young Tableaux of shape L/M
#using The formula. Try:
#For example, try: NuSSYTc([2,2],[1,1]);

NuSSYTc:=proc(L,M) local i,gu,x,gu1:

if not (nops(M)<=nops(L) and   min(seq(L[i]-M[i],i=1..nops(M)))>=0) then
 RETURN(0):
fi:

x:=[op(M),0$(nops(L)-nops(M))]:
gu:=Imags(x):

add(gu1[2]*MNC(L-gu1[1]),gu1 in gu):
end:




#Igor(L,M): NuSSYTc(L,M)/MNC(L). Try:
#Igor([n,n,n],[2,1]);

Igor:=proc(L,M) local gu,x,gu1:

x:=[op(M),0$(nops(L)-nops(M))]:
gu:=Imags(x):

normal(add(gu1[2]*simplify(MNC(L-gu1[1])/MNC(L)),gu1 in gu)):
end:




#Swee(L,M): Igor(L,M)/Igor(L,[]). Try:
#Swee([n,n,n],[2,1]);
Swee:=proc(L,M):
normal(Igor(L,M)/Igor(L,[])):
end:


#Pars1(n,r): all the partitions of n with largest part exactly r
Pars1:=proc(n,r) local gu,r1,mu,mu1:
option remember:

if r>n then
 RETURN({}):
fi:

if r=n then
 RETURN({[r]}):
fi:

gu:={}:

for r1 from 1 to r do
 mu:=Pars1(n-r,r1):
 gu:=gu union {seq([r,op(mu1)],mu1 in mu)}:
od:

gu:

end:



#Shaps(n,r,c2,s): all the shapes of size n with largest part  r and at most s parts, that contain the cell c2. Try:
#Shaps(10,5,[3,2],3);
Shaps:=proc(n,r,c2,s) local mu,mu1,gu:
mu:=Pars1(n,r):
gu:={}:

for mu1 in mu do
 if c2[1]<=nops(mu1) and c2[2]<=mu1[c2[1]] and nops(mu1)<=s then
   gu:=gu union {mu1}:
 fi:
od:

gu:
end:



#PrS1(S,i,k,c2): Inputs symbolic shape S, a pos. integer k and a pair c2 indicating a location such that c2[1]<=nops(S) and c2[2]<i, outputs
#the rational function describing the probability that a random standard Young tableau of shape S would have cell [1,i] occupied by k
#and the occupant of cell c2 be smaller than k. Try:
#PrS1(3,n,2,3,[2,1]);
PrS1:=proc(S,i,k,c2) local mu,mu1,gu,mu1a:

if not (c2[1]<=nops(S) and c2[2]<i) then
 RETURN(FAIL):
fi:

mu:=Shaps(k-1,i-1,c2,nops(S)):

gu:=0:

for mu1 in mu do
 mu1a:=[mu1[1]+1,op(2..nops(mu1),mu1)]:
 gu:=normal(gu+YF(mu1)*Swee(S,mu1a)):
end:
normal(simplify(gu)):
end:


#PrS(S,i,c2): Inputs a symbolic shape S, a Numeric pos. integer i, and a pair c2 indicating a location such that c2[1]<=nops(S) and c2[2]<i, outputs
#a closed form expression as a rational function in n of Pr(S,[1,i],c2); Try:
#PrS[n,n],2,[2,1]);
PrS:=proc(S,i,c2) local k,gu,R:
option remember:
R:=nops(S):
if (nops(c2)>nops(S) or c2[2]>=i) then
 RETURN(FAIL):
fi:
 
gu:=add(PrS1(S,i,k,c2),k=i..(i-1)*R+1):
factor(normal(2*gu-1)):

end:



#CheckPrSa(R,n,i,c2,K): Checks procedure PrS(R,n,i,c2) for n from i to K. Try:
#CheckPrSa(2,n,1,[2,1],9);
CheckPrSa:=proc(S,n,i,c2,K) local gu,j:
gu:=PrS(S,i,c2):
evalb({seq(subs(n=j,gu)-Pr(subs(n=j,S),[1,i],c2),j=i..K)}={0}):
end:


#CheckPr(S,n,K): Checks procedure PrS(S,i,c2) for i from 2 to K and all c2 possible. Try:
#CheckPr([n,n],n,6);
CheckPrS:=proc(S,n,K) local i,c2,j1,j2:

for i from 2 to K do
 for j1 from 2 to nops(S) do
  for j2 from 1 to i-1 do
   c2:=[j1,j2]:
   if not CheckPrSa(S,n,i,c2,K+1) then
     print([S,i,c2]):
     RETURN(false):
   fi:
 od:
od:
od:
true:
end:



#AnijOld(n,i,j): Same as PrS([n,n],i,j) (or Pr([n,n],[1,i],[2,j]), for numeric n). Try:
#AnijOld(n,2,1);
AnijOld:=proc(n,i,j)  local k,gu:
if not (i>1 and j>=1 and j<i) then
 RETURN(FAIL):
fi:
gu:=(n+1)/binomial(2*n,n)*add((k-1)!*(2*n-k)!*(2*i-k+1)*(2*i-k)/i!/(k-i)!/(n-i)!/(n-k+i+1)!,k=i+j..2*i-1):
factor(simplify(2*gu-1)):
end:


#FindZero2(K): Finds all pairs [[1,i],[2,j]], 2<=j<i<=K such that the limit of Anij(n,i,j) as n goes to infinity  is 0. Try:
#FindZero2(10);
FindZero2:=proc(K) local i,j,gu,mu,n:
gu:={}:
for i from 2 to K do
 for j from 1 to i-1 do
  mu:=Anij(n,i,j):
  if limit(mu,n=infinity)=0 then
    gu:=gu union {[[1,i],[2,j]]}:
  fi:
 od:
od:
gu:
end:



#A1nij(n,i,j): Number of Standard Young tableau of shape [n,n] such that the occupant of [1,i] is larger than the occupant of [2,j]. Try:
#A1nij(n,2,1);
A1nij:=proc(n,i,j)  local k,gu:
#if not (i>1 and j>=1 and j<i) then
# RETURN(FAIL):
#fi:
#gu:=add(simplify((k-1)!*(2*n-k)!*(2*i-k+1)*(2*i-k)/i!/(k-i)!/(n-i)!/(n-k+i+1)!/(2*n)!*n!*(n+1)!),k=i+j..2*i-1):
add(OcCsD(n,k,i),k=i+j..2*i-1):
end:


Anij:=proc(n,i,j)
factor(simplify(factor(simplify(2*A1nij(n,i,j)-1)))):
end:

#MinPr2(n): For numeric (pos. integer) n, Same as MinPr([n,n]) but much faster, using Anij(n,i,j). Try: MinPr2(30);
MinPr2:=proc(n) local i,j,si, aluf,hal:
aluf:={[[1,2],[2,1]]}:
si:=abs(Anij(n,2,1)):

for i from 3 to n do
 for j from 1 to i-1 do
  hal:=abs(Anij(n,i,j)):
   if  hal<si then
   aluf:={[[1,i],[2,j]]}:
    si:=hal:
   elif  hal=si then
     aluf:=aluf union {[[1,i],[2,j]]}:
  fi:
 od:
od:
si,aluf:
end:




#MinPr2Old(n): For numeric (pos. integer) n, Same as MinPr([n,n]) but much faster, using Anij(n,i,j). Try: MinPr2(30);
MinPr2Old:=proc(n) local i,j,si, aluf,hal:
aluf:={[[1,2],[2,1]]}:
si:=abs(Anij(n,2,1)):

for i from 3 to n do
 for j from 1 to i-1 do
  hal:=abs(AnijOld(n,i,j)):
   if  hal<si then
   aluf:={[[1,i],[2,j]]}:
    si:=hal:
   elif  hal=si then
     aluf:=aluf union {[[1,i],[2,j]]}:
  fi:
 od:
od:
si,aluf:
end:





#OcGF(L,c,x): Given a shape L and a cell in it,and a variable x, outputs the prob. generating function, in the variable x, for the occupant of cell c in a random Young tableau of shape L. Try:
#OcGF([4,4,4],[2,2],x);
OcGF:=proc(L,c,x) local gu,gu1:
if not (c[1]>=1 and c[1]<=nops(L) and c[2]<=L[c[1]]) then
 print(`cell `, c, `does not belong to shape`, L):
 RETURN(FAIL):
fi:


gu:=SYT(L):

add(x^gu1[c[1]][c[2]],gu1 in gu)/nops(gu):

end:





#FindZero(L,n,K): Given a symbolic shape L pharased in terms of a single variable n, and a positive integer K
#finds all pairs [i,c] where i<=K, and c is a cell with 2<=c[1]<=nops(L) and c[2]<i such that PrS(L,i,c) goes to zero to zero as n goes to infinity
#FindZero([n,n],n,K) is the same as FindZero2(K)
#FindZero([n,n,n],n,3);
FindZero:=proc(L,n,K) local i,j1,j2,gu,mu,c:
gu:={}:
for i from 2 to K do
 for j1 from 2 to nops(L) do
  for j2 from 1 to i-1 do
    c:=[j1,j2]:
  mu:=PrS(L,i,c):
  if limit(mu,n=infinity)=0 then
    gu:=gu union {[[1,i],c]}:
  fi:
 od:
od:
od:
gu:
end:




#Pars2(n,s): all the partitions of n with at most s parts. Try:
#Pars2(6,2);
Pars2:=proc(n,s) local gu,gu1,i1:
option remember:

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

if n<s then
 RETURN(Pars2(n,n)):
fi:



gu:=Pars2(n-s,s):

{seq([seq(gu1[i1]+1,i1=1..nops(gu1)),1$(s-nops(gu1))],gu1 in gu)} union Pars2(n,s-1):

end:

#Pars3(n,s,c): all the partitions of n with at most s parts where cell c2 is a corner cell. Try:
#Pars3(7,3,[2,2]);

Pars3:=proc(n,s,c) local gu,mu,mu1:

if not (c[1]>=1 and c[1]<=s) then
 RETURN(FAIL):
fi:

mu:=Pars2(n,s):

gu:={}:
for mu1 in mu do
if nops(mu1)>=c[1] and mu1[c[1]]=c[2] and (nops(mu1)=c[1] or mu1[c[1]+1]<mu1[c[1]]) then
 gu:=gu union {mu1}:
fi:
od:
gu:
end:



#OcCs(L,k,c): Given a symbolic shape L, a numeric integer k and a numeric cell c, finds the rational function in the symbols used in L
#that describes the exact probability that a random Young tableau of shape L, T, will have  T[c[1]][c[2]]=k. Try:
#OcCs([n,n],3,[1,2]);
OcCs:=proc(L,k,c) local gu,mu,mu1,mu1a:

mu:=Pars3(k,nops(L),c):


gu:=0:

for mu1 in mu do
mu1a:=[op(1..c[1]-1,mu1),mu1[c[1]]-1,op(c[1]+1..nops(mu1),mu1)]:
if mu1a[nops(mu1a)]=0 then
  mu1a:=[op(1..nops(mu1a)-1,mu1a)]:
fi:
gu:=normal(gu+YF(mu1a)*Swee(L,mu1)):
od:

gu:=normal(simplify(gu,factorial)):
factor(gu):

end:



#OcGFs(L,c,x,K): Given a symbolic shape L and a numeric, cell ,and a variable x, outputs the prob. generating function, in the variable x, for the occupant of cell c in a random Young tableau of shape L
#conditioned that it is <=K. It also gives you the probability of it being <=K.  Try:
#OcGFs([n,n],[1,5],x,9);
OcGFs:=proc(L,c,x,K) local k,gu,tot,i:

gu:=add(OcCs(L,k,c)*x^k,k=c[1]*c[2]..K):

tot:=normal(subs(x=1,gu)):
if tot=0 then
RETURN(0,0):
fi:

gu:=gu/tot:
gu:=add(factor(coeff(gu,x,i))*x^i,i=0..degree(gu,x)):
add(factor(coeff(gu,x,i))*x^i,i=ldegree(gu,x)..degree(gu,x)),factor(simplify(tot)):
end:





#OcGFs1(L,i,x): Given a symbolic shape L and a numeric positive integer i, and a variable x, outputs the prob. generating function, in the variable x, for the occupant of cell [1,i] in a random Young tableau of shape L
#It is the same as  OcGFs(L,[1,i],x,nops(L)*(i-1)+1)[1]: Try:
#OcGFs1([n,n,n,n],2,x);
OcGFs1:=proc(L,i,x) option remember:
OcGFs(L,[1,i],x,(i-1)*nops(L)+1)[1]:
end:



#Stat(F,t,K): Given a generating function F of t, for a r.v. finds the statistical information. Try:
#Stat((1_t)^10000,4);
Stat:=proc(F,t,K) local gu,mu,sd,k,f:
f:=F/subs(t=1,F):
mu:=subs(t=1,diff(f,t)):
gu:=[mu]:
f:=f/t^mu:
f:=t*diff(t*diff(f,t),t):
sd:=sqrt(subs(t=1,f)):
gu:=[op(gu),sd]:
for k from 3 to K do
 f:=t*diff(f,t):
gu:=[op(gu),subs(t=1,f)/sd^k]:
od:

gu:


end:


#OcGFs1L(L,n,i,x,K): Given a symbolic shape L phrased in terms of the  variable n, outputs the limiting distribution as n goes to infinity 
#of the occupants of the cell [1,i] in a random Young tableau of shape L
#followed by a list consisting of the average, s.d, and the scaled moments through the K-th.
#Try:
#OcGFs1L([n,n],n,10,x,4);
OcGFs1L:=proc(L,n,i,x,K)  local gu,lu,i1:
gu:=OcGFs1(L,i,x):
lu:=expand(limit(gu,n=infinity)):
[seq([i1,coeff(lu,x,i1)],i1=ldegree(lu,x)..degree(lu,x))],Stat(lu,x,K):
end:



#Paper1(R,n,K1,K2): inputs a positive integer R>=2, a symbol, n, and integers K1 and K2, outputs a paper with explicit expressions in n
#for the probability distribution of "occupant of cell [1,i]" for a random Young Tableau of shape [n$R], for all i from 2 to K1. It also gives
#the limiting distribution as n goes to infinity, and average, s.d. and the (scaled) moments up to the K2-th. Try:
#Paper1(2,n,10,4);
Paper1:=proc(R,n,K1,K2) local gu,lu,i,x,j,mu:
print(``):
print(`The Distribution of the Occupant of Cell [1,i] in a random Standard Young tableau of shape`, [n$R], `and its Limiting behavior as n goes to infinity for i from 2 to`, K1):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

for i from 2 to K1 do
print(``):
print(`---------------------------------------------`):
print(``):
print(`The occupants of cell`, [1,i], `in a standard Young tableau of shape`, [n$R], `are all the integers from`, i, ` to `, R*(i-1)+1):
gu:=OcGFs1([n$R],i,x):
print(` The probability distribution is`):
print(``):

print([seq(factor(coeff(gu,x,j)),j=ldegree(gu,x)..degree(gu,x))]):
print(``):
print(`and in Maple notation`):
print(``):
lprint([seq(factor(coeff(gu,x,j)),j=ldegree(gu,x)..degree(gu,x))]):
print(``):
mu:=Stat(gu,x,2):

print(``):
lprint(`The average is`):
print(``):
lprint(factor(mu[1]) ):
print(``):
print(`and the variance is`):
print(``):
lprint( factor(mu[2]^2)):
print(``):
lu:=OcGFs1L([n$R],n,i,x,K2):
print(``):
print(`as n goes to infinity, the distribution is`):
print(``):
lprint(evalf(lu[1])):
print(``):
print(`The limiting average, standard deviation up to the`, K2, ` -th scaled-moment are`):
print(``):
print(lu[2]):
print(``):
print(`and in floating-point`):
print(``):
lprint(evalf(lu[2])):
print(``):
print(`Here is a plot`):
print(``):
print(plot(lu[1]));
print(``):

od:

end:





#Paper2(R,n,K1): inputs a positive integer R>=2, a symbol, n, and integer K1  outputs a paper with explicit expressions in n
#for the Pr(T[1,i]>T[k,j])-Pr(T[k,j]>T[1,i]) (alias 2*Pr(T[1,i]>T[k,j])-1) for a random Young tableau of shape [n$R]
#for all i from 2 to K1, k from 2 to R, and j from 1 to i-1.
# It also gives the limiting quantity and the "cut-off" place where the limit switches from being pos. to negative. Try:
#Paper2(2,n,10);
Paper2:=proc(R,n,K1) local gu,lu,i,j,i1,k:
print(``):
print(`The Sorting Probabilities  of The entries in the first row vs. those not related to it in lower rows`):
print(`in a random Standard Young tableau of shape`, [n$R], `and its Limiting behavior as n goes to infinity for i from 2 to`, K1):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

for i from 2 to K1 do
print(``):
print(`---------------------------------------------`):
print(``):
for k from 2 to R do

if R=2 and k=2 then
gu:=[seq(Anij(n,i,j),j=1..i-1)]:
else
gu:=[seq(PrS([n$R],i,[k,j]),j=1..i-1)]:
fi:

lu:=[seq(limit(gu[j],n=infinity),j=1..nops(gu))]:

print(`The rational functions describing the sorting probabilities of the cell`, [1,i], `vs. those in the`, k, ` -th row  from j=1 to j=`, i-1, `are as follws `):
print(``):
print(gu):
print(``):
print(`and in Maple notation`):
print(``):
lprint(gu):
print(``):

print(`The limits, as n goes to infinity are`):
print(``):
print(lu):
print(``):
print(``):
print(`and in Maple notation`):
print(``):
lprint(lu):
print(``):
print(`and in floating point`):
print(``):
lprint(evalf(lu)):
print(``):
for i1 from 1 while lu[i1]>0 do od:
print(``):
print(`The cut off is at j=`,i1):
print(``):
od:

od:

print(`-------------------------`):
print(``):
print(`This ends this article that took`, time(), `seconds to produce `):
end:





#SiOcGF(L,c,x,K): Given a shape L and a cell in it,and a variable x, outputs the APPROXIMATE prob. generating function, in the variable x, for the occupant of cell c in a random Young tableau of shape L. 
#using simulation (via the Greene-Nijenhuis-Wilf algorithm), with K tries. Just for fun! Try:
#SiOcGF([4,4,4],[2,2],x,1000);
SiOcGF:=proc(L,c,x,K) local i:
if not (c[1]>=1 and c[1]<=nops(L) and c[2]<=L[c[1]]) then
 print(`cell `, c, `does not belong to shape`, L):
 RETURN(FAIL):
fi:

evalf(add(x^GNW(L)[c[1]][c[2]],i=1..K)/K):

end:




#SiPr(L,c1,c2,K):Same as Pr(L,c1,c2) but approximately,simulating with K tries, using the Greene-Nijenhuis-Wilf algoirthm. Try:
#SiPr([3,3,3],[1,2],[2,1],1000);
SiPr:=proc(L,c1,c2,K) local co,gu1,i:

if not (c1[1]<=nops(L) and c2[1]<=nops(L) and c1[2]<=L[c1[1]] and c2[2]<=L[c2[1]]) then
 print(`Bad input`):
 RETURN(FAIL):
fi:

co:=0:

for i from 1 to K do

gu1:=GNW(L):

if gu1[c1[1]][c1[2]]>gu1[c2[1]][c2[2]] then
 co:=co+1:
fi:

od:

co:=evalf(co/K):
2*co-1:
end:





#Paper1L(R,K1,K2): An abbreviated version of Paper1(R,n,K1,K2), only giving the limits as n goes to infinity
#Paper1L(2,10,4);
Paper1L:=proc(R,K1,K2) local lu,i,x,n:
print(``):
print(`The Limiting Distribution, as n goes to infinity, of the Occupant of Cell [1,i] in a random Standard Young tableau of shape`, [n$R], ` from i=2 to i=`, K1):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

for i from 2 to K1 do
print(``):
print(`---------------------------------------------`):
print(``):
print(`The occupants of cell`, [1,i], `in a standard Young tableau of shape`, [n$R], `are all the integers from`, i, ` to `, R*(i-1)+1):

lu:=OcGFs1L([n$R],n,i,x,K2):
print(``):
print(`as n goes to infinity, the distribution is`):
print(``):
lprint(evalf(lu[1])):
print(``):
print(`The limiting average, standard deviation up to the`, K2, ` -th scaled-moment are`):
print(``):
lprint(evalf(lu[2])):
print(``):
print(`Here is a plot`):
print(``):
print(plot(lu[1]));
print(``):

od:

end:




#Paper2L(R,K1): An abbreviated version of Paper2(n,R,K1), only stating limiting values
#for the limit of Pr(T[1,i]>T[k,j])-Pr(T[k,j]>T[1,i]) (alias 2*Pr(T[1,i]>T[k,j])-1) for a random Young tableau of shape [n$R]
#as n goes to infinity. Try:
#Paper2L(2,10);
Paper2L:=proc(R,K1) local n,gu,lu,i,j,i1,k:
print(``):
print(`The Limiting Sorting Probabilities  of The entries in the first row vs. those not related to it in lower rows`):
print(`in a random Standard Young tableau of shape`, [n$R], `as n goes to infinity for i from 2 to`, K1):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

for i from 2 to K1 do
print(``):
print(`---------------------------------------------`):
print(``):
for k from 2 to R do

if R=2 and k=2 then
gu:=[seq(Anij(n,i,j),j=1..i-1)]:
else
gu:=[seq(PrS([n$R],i,[k,j]),j=1..i-1)]:
fi:

lu:=[seq(limit(gu[j],n=infinity),j=1..nops(gu))]:

print(`The limits of the sorting probabilities of the cell`, [1,i], `vs. those in the`, k, ` -th row  from j=1 to j=`, i-1, ` as n goes to infinity, are as follows `):
print(``):
lprint(evalf(lu)):
print(``):
for i1 from 1 while lu[i1]>0 do od:
print(``):
print(`The cut off is at j=`,i1):
print(``):
od:

od:

print(`-------------------------`):
print(``):
print(`This ends this article that took`, time(), `seconds to produce `):
end:


#Paper3(S,n,n0,i,K): Given a symbolic shape S, phrased in terms of the variable n, a positive integer n0, a positive integer i, and a large positive integer K
#illustates both the efficiently of symbolic computation, and that of simulation using the beautiful Greene-Nijenhuis-Wilf algorithm to generate
#a random Standard Young Tableau of, by comparing the exact probability distribution of the occupants of the cell [1,i], and what you get
#by sampling K random standard Young tableaux using the GNW algorithm. The former uses procedure OcGFs1 (q.v) and the later procedure SiOcGF. Try:
#Paper3([3*n,2*n],n,10,5,100);
Paper3:=proc(S,n,n0,i,K) local gu,S0,t0,x,gu0,lu:
print(``):

if i>n0 then
 print(i, `should be smaller than`, n0):
  RETURN(FAIL):
fi:

S0:=subs(n=n0,S):
print(`The Exact Probability Distirbution of the occupants of cell`, [1,i], ` in a random Standard Young Tableau of shape`, S0, `using the exact symbolic expression`):
print(``):
print(`Versus the approximation obtained by simulation via the Greene-Nijenhuis-Wilf algorithm`, K, `times `):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):
gu:=OcGFs1(S,i,x):
print(`The occupants of cell`, [1,i], `in a standard Young tableau of symbolic shape`, S, `are all the integers from`, i, ` to `, nops(S)*(i-1)+1):
t0:=time():
gu:=OcGFs1(S,i,x):
print(` The probability generating function, using the variable x  (i.e. the polynomial whose coeff. of x^j is the probability that Y[1,i]=j) is`):
print(``):

lprint(gu):
print(``):
print(`BTW this took`, t0, `seconds to compute `):
print(``):
t0:=time():
print(`We are interested in what happens when n=`, n0, `i.e. for the specific shape`, S0):
print(``):
gu0:=sort(expand(evalf(subs(n=n0,gu)))):
print(``):
print(`Then the prob. gen. function is`):
print(``):
print(gu0):
print(``):
print(`Note that there are`, YF(S0), `Standard Young Tableaux of shape`, S0, `so it would be impractical to compute this directly. Let's do it by simulation`):
print(``):
print(`Now let's compare it with the approximation that you get by generating`, K, `standard Young tableau of that shape, namely`, S0, `using the Green-Nijenhuis-Algorithm`):
print(``):
lu:=sort(SiOcGF(S0,[1,i],x,K)) :
print(``):
print(lu):
print(``):
print(`BTW this took`, time()-t0, ` seconds `):
print(``):
print(`The difference is`):
print(``):
print(gu0-lu):
print(``):

end:





#Paper4(S,n,n0,i,c2,K): Given a symbolic shape S, phrased in terms of the variable n, a positive integer n0, a positive integer i, and a numeric cell c2
#and a large integer K
#illustates both the efficiently of symbolic computation, and that of simulation using the beautiful Greene-Nijenhuis-Wilf algorithm to 
#a random Standard Young Tableau of, by comparing the sorting probability of cell [1,i] vs. cell c2, and what you get
#by sampling K random standard Young tableaux using the GNW algorithm. The former uses procedure PrS (q.v) and the later procedure SiPr (q.v.). Try:
#Paper4([3*n,2*n],n,10,2,[2,1],100);
Paper4:=proc(S,n,n0,i,c2,K) local gu,S0,t0,gu0,lu:
print(``):

if i>n0 then
 print(i, `should be smaller than`, n0):
  RETURN(FAIL):
fi:

if not (c2[1]>1 and c2[2]<i) then
 print(c2, `is not unrelated to`, [1,i]):
 RETURN(FAIL):
fi:

S0:=subs(n=n0,S):
print(`The Exact Sorting Probability Distirbution of cell`, [1,i], ` vs. the cell `, c2, `in a random Standard Young Tableau of shape`, S0, `using the exact symbolic expression`):
print(``):
print(`Compared to  the approximation obtained by simulation via the Greene-Nijenhuis-Wilf algorithm`, K, `times `):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):
gu:=PrS(S,i,c2):

print(`The exact value of the prob. that the occupant of cell`, [1,i], `is larger than the occupant of cell `, c2, `minus the prob. that the opposite is true for symbolic shape`, S, ` is `):
t0:=time():
gu:=PrS(S,i,c2):
print(``):
lprint(gu):
print(``):
print(`BTW this took`, t0, `seconds to compute `):
print(``):
t0:=time():
print(`We are interested in what happens when n=`, n0, `i.e. for the specific shape`, S0):
print(``):
gu0:=evalf(subs(n=n0,gu)):
print(``):
print(`the value is`):
print(``):
print(gu0):
print(``):
print(`Note that there are`, YF(S0), `Standard Young Tableaux of shape`, S0, `so it would be impractical to compute this directly. Let's do it by simulation`):
print(``):
print(`Now let's compare it with the approximation that you get by generating`, K, `standard Young tableau of that shape, namely`, S0, `using the Green-Nijenhuis-Algorithm`):
print(``):
lu:=SiPr(S0,[1,i],c2,K):
print(``):
print(lu):
print(``):
print(`BTW this took`, time()-t0, ` seconds `):
print(``):
print(`The ratio is`):
print(``):
print(gu0/lu):
print(``):

end:









#MinPrk(k,N): For numeric (pos. integer) N, Similar to MinPr([n$k])]) but c1 is in the first row
# MinPrk(3,30);
MinPrk:=proc(k,n1) local i,j,si, aluf,hal,n,k1:


aluf:={[[1,2],[2,1]]}:
si:=abs(subs(n=n1,PrS([n$k],2,[2,1]))):

for i from 2 to n1 do
 for k1 from 2 to k do
 for j from 1 to i-1 do
  hal:=abs(subs(n=n1,PrS([n$k],i,[k1,j]))):
   if  hal<si then
   aluf:={[[1,i],[k1,j]]}:
    si:=hal:
   elif  hal=si then
     aluf:=aluf union {[[1,i],[k1,j]]}:
  fi:
 od:
od:
od:
si,aluf:
end:


#Paper5(k,N): inputs a positive integer k>1 and a positive integer N outputs the minimal sorting probabilities of the shape [n$k]  for n from 1 to N.
#If k>2 it is only an upper bound, since it only takes into account comparisons with the cells of the first row. Try:
#Paper5(3,10):
Paper5:=proc(k,N) local gu,n:

print(`The Minimal sorting probabilities, and the champion pairs for Standard Young Tableau of shape`, [n$k] , `for n from 1 to`, N):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

for n from 2 to N do

if k>2 then
gu:=MinPrk(k,n)
else
gu:=MinPr2(n)
fi:

if k=2 then
print(`The minimal sorting probability of Standard Young tableaux of shape`, [n$k], ` is `, gu[1], `and in floating point`, evalf(gu[1]), `and this happens with the pairs of cells`, gu[2]):
else
print(`The minimal sorting probability of Standard Young tableaux of shape`, [n$k], ` but where one of the cells is in the first row is `, gu[1], `and in floating point`, evalf(gu[1]), `and this happens with the pairs of cells`, gu[2]):
fi:
od:
print(``):
print(`-------------------------------------`):
print(``):
print(`This took`, time(), `seconds. `):

end:


###new stuff






#OcCsD(n,k,i): Same as OcCs([n,n],k,[1,i]) but faster. Try:
#OcCsD(n,3,2);
OcCsD:=proc(n,k,i):

#if k<i or k>2*i-1 then
# 0:
#else
factor(expand((k-1)!*(2*n-k)!*(2*i-k+1)*(2*i-k)/i!/(k-i)!/(n-i)!/(n-k+i+1)!/(2*n)!*n!*(n+1)!)):
#fi:
end:

#LimMomsD(i,K1,K2): inputs a symbol i and a positive integers K1, K2, outputs a list of length K1 whose first entry is the limiting expectation (as n goes to infinity) of
#the occupant of cell [1,i] in a standard Young tableau of shape [n,n], followed by the variance, and scaled third-through-K th moment. 
#It also returns the asympotics to order K2, and the limiting scaled moments up to the K1, and the floating-point version, and finally the asymptotic to order K2 of the scaled moments
#Try:
#LimMomsD(i,4);
LimMomsD:=proc(i,K1,K2) local gu,k,n,j,mu,che,mu1,ku,ku1,r,T,kuA,i1,ru,ru1,vu,vu1:

gu:=OcCsD(n,k,i):

gu:=limit(gu,n=infinity):

che:=simplify(sum(gu,k=i..2*i-1)):
if che<>1 then
 print(`Something bad happened`):
 RETURN(FAIL):
fi:

mu:=convert(simplify(sum(k*gu,k=i..2*i-1)),factorial):

T[0]:=1:
T[1]:=mu:

for j from 2 to K1 do
 mu1:=convert(simplify(sum(k^j*gu,k=i..2*i-1)),factorial):
T[j]:=mu1:
od:

ku:=[mu]:

for j from 2 to K1 do
ku1:=simplify(expand(add(binomial(j,r)*T[r]*(-mu)^(j-r),r=0..j))):
ku:=[op(ku),ku1]:
od:

kuA:=[]:

for j from 1 to K1 do

 ku1:=asympt(ku[j],i,K2+j):
 ku1:=add(simplify(op(i1,ku1)),i1=1..nops(ku1)):
 kuA:=[op(kuA),ku1]:
od:

ku,kuA:

ru:=[0,1]:
vu:=[0,1];
for j from 3 to K1 do
 ru1:=factor(simplify(limit(ku[j]/ku[2]^(j/2),i=infinity))):
 vu1:=asympt(ku[j]/ku[2]^(j/2),i,K2+j):
 vu1:=add(simplify(op(i1,vu1)),i1=1..nops(vu1)):

ru:=[op(ru),ru1]:
vu:=[op(vu),vu1]:

od:
[ku,kuA,ru,evalf(ru),vu]:
end:






#Paper6(i,K1,K2): A verbose form of LimMomsD(i,K1,K2) (q.v.). Try:
#Paper6(i,10,10):
Paper6:=proc(i,K1,K2) local gu,j,t0:

t0:=time():
print(`The Statistics of the Limiting Occupancy of the cell [1,i] in a 2-rowed Standard Young tableau of shape [n,n], as n goes to infinity, and its Asympotic behavior as i goes to infinity`):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):
gu:=LimMomsD(i,K1,K2):

print(`In this article we will discuss the limiting statistics of the occupant of the cell [1,i] in a 2-rowed standard Young  tableau of shape [n,n] as n goes to infinity, and later as i goes to infinity`):
print(``):
print(`Of course the occupant can be anywhere between i and 2*i-1, regardless of n. As n goes to infinity`):

print(` the expected size of the occupant of cell [1,i] is`):
print(``):
print(gu[1][1]):
print(``):
print(`and in Maple notation`):
print(``):
lprint(gu[1][1]):
print(``):
print(`This is asymptotically, as i goes to infinity`):
print(``):
print(gu[2][1]):
print(`and in Maple notation`):
print(``):
lprint(gu[2][1]):
print(``):


print(` the variance of the occupant of cell [1,i] is`):
print(``):
print(gu[1][2]):
print(``):
print(`and in Maple notation`):
print(``):
lprint(gu[1][2]):
print(``):
print(`This is asymptotically, as i goes to infinity`):
print(``):
print(gu[2][2]):
print(`and in Maple notation`):
print(``):
lprint(gu[2][2]):
print(``):

for j from 3 to K1 do
print(`The `, j, `-th moment about the mean is`):
lprint(gu[1][j]):
print(``):
print(`This is asymptotically, as i goes to infinity`):
print(``):
print(gu[2][j]):
print(`and in Maple notation`):
print(``):
lprint(gu[2][j]):
print(``):
od:

print(`The limiting scaled moments about the mean, as i goes to infinity are `):
print(``):
print(gu[3]):
print(``):
print(`and in Maple notation`):
print(``):
lprint(gu[3]):
print(`and in floating point`):
lprint(evalf(gu[3])):
print(``):
print(`In more detail, the asymptotics is`):
print(``):
print(gu[5]):
print(``):
print(`and in Maple notation`):
print(``):
lprint(gu[5]):
print(``):
print(`and in floating point`):
lprint(evalf(gu[5])):
print(``):
print(`---------------------------------------------------------------`):
print(``):
print(`This ends this article that took`, time()-t0, `seconds to generate. `):
print(``):
end: