# OK to post

# Assignment 10
# February 22, 2022
# Kayla Gibson

# 1. Write a procedure CountInstabilities(M,W,pi)
# CountInstabilities takes as input M and W lists of preferences, and counts how many of the 
# n*(n-1) ordered pairs of husbands 1<=i1,i2<=n are such that i1 is married to j1 and i2 is 
# married to j2, but i1 would rather be married to j2 (instead of his current wife j1) and j2
# would rather be married to i1 (rather than her current husband i2)

CountInstabilities:=proc(M,W,pi) local i,j,L,Instabilities,n:
Instabilities:=[]:
n:=nops(M[1]):
for i from 1 to n do
   for j from 1 to n do
      L:=IsStable1(M,W,pi,i,j):
      if L=false and i<>j then
         Instabilities:=[op(Instabilities),[i,j]]:
      fi:
   od:
od:
RETURN(Instabilities):
end:


# 2. Write a procedure GFstab(n,K,x) 
# GFstab(n,K,x)  takes as input a positive integer n, a large positive integer K, and a variable
# name, say x, it generates random ranking tables RN(n) K times and a random matching, pi, and 
# outputs the polynomial in x (of degree <=n(n-1)), whose coefficient of x^i is the number of 
# those random choices that lead to exactly i instabilities. By taking the first and second
# derivatives w.r.t. x, estimate the expectation and variance of the "number of instabilities".

GFstab:=proc(n,K,x) local i,j,k,R,pn,coeffs:
coeffs:=[seq(0,j=1..(n*(n-1))+1)]:
for i from 1 to K do
   R:=RT(n):
   pn:=randperm(n):
   j:=nops(CountInstabilities(R[1],R[2],pn))+1:
   coeffs[j]+=1:
od:
diff(sum(coeffs[k]*(x^(k-1)),k=1..n*(n-1)+1),x):
end: