#Please do not post homework
#Caroline Cote, April 5th, 2026, Assignment 18
HelpHW18:= proc(): print(``): end:

#Problem 1

Moms((x+x^2+x^3+x^4+x^5+x^6)^n, x, 2);
#                    [7    1    (1/2)  (1/2)]
#                    [- n, - 105      n     ]
#                    [2    6                ]

#The mean is 7n/2 and the s.d. is sqrt(sqrt(105*n)/6).

F:=SMoms((x+x^2+x^3+x^4+x^5+x^6)^n, x, 20);

evalf([seq(limit(F[i],n=infinity), i = 3..20)]);
#    [                                                     
#    [0., 3., 0., 15., 0., 105., 0., 945., 0., 10395., 0., 

#                5                 6                  7      
#      1.35135 10 , 0., 2.027025 10 , 0., 3.4459425 10 , 0., 
#
#                   8]
#      6.54729075 10 ]


evalf([seq(subs(x=0, diff(exp(x^2/2), x$k)),k=3..20)]);
#    [                                                     
#    [0., 3., 0., 15., 0., 105., 0., 945., 0., 10395., 0., 

#                5                 6                  7      
#      1.35135 10 , 0., 2.027025 10 , 0., 3.4459425 10 , 0., 

#                   8]
#      6.54729075 10 ]


evalb(%=%%);
#                              true


ave22:= proc(n) local S,s:
S := SYT([n,n,n]):
add(s[2][2], s in S)/nops(S):
end:

seq(avg22(n), n=2..7);
#                   23     405  778  797  6271
#                   --, 5, ---, ---, ---, ----
#                   5      77   143  143  1105

evalf(%);
#    4.600000000, 5., 5.259740260, 5.440559441, 5.573426573, 5.675113122

#this is very off but it's almost like 6 - 2/n for all of them 



#Problem 2
#AppxAve22(n,K): takes K random standard Young tableaux of shape
#[n,n,n], finds its [2,2] entry and takes the average
AppxAve22:= proc(n,K) local k, S, s, sum: 
sum:= 0:
for k from 1 to K do
s:= RandSYT([n,n,n]):
sum:= sum+s[2][2]:
od:
sum/K:
end:

AppxAve22(30,1000);
#                                 781
#                                 ---
#                                 125

AppxAve22(30,10000);
#                                 3901
#                                 ----
#                                 625 

evalf(%); evalf(%%);
#                          6.241600000
#                          6.241600000

#They are very close to each other and they do not match my conjecture


#Problem 3
with(combinat):
#AppxJike(n,x,K): approximates the prob. gen function of the weight 
#pi -> nops(RS(pi)[1]) by taking K random permutations 
#(using randperm(n)
AppxJike := proc(n,x,K) local i,pi, sum:
sum:=0:
for i from 1 to K do
   pi:= randperm(n):
   sum:= sum + x^nops(RS(pi)[1]):
od:
sum/K:
end:

PlotDist(AppxJike(100,x,1000),x);
PlotDist(AppxJike(200,x,1000),x);

#There's not really a way to copy and paste the graphs to here but they look somewhat similar
#For n=100, its a bell curve that peaks around 0.22 and ranges from 12 to 22 ish.
#For n=200, its a bell curve that peaks around 0.2 and ranges from 20 to 30 ish.