#Ok to post
#Michael Yen, Assignment 18, 11/15/20

#Lecture 18: Due Nov. 15, 9:00pm. Email ShaloshBEKhad@gmail.com an attachment

#hw18FirstLast.txt

#Indicate whether it is OK to post

#1) [Better Version of Nov. 11, 2020 ]
#Generalize procedure NuConG(N) of M18.txt to write a procedure

#NuKcomponentsGraphs(N,K)

#that inputs a positive integer N and another positive integer K and outputs the first N terms of the 
#sequence enumerating the number of labeled graphs on n vertices with Exactly K componets. Make sure 
#that NuKComponentsGraphs(N,1) equals NuConG(N) for N ≤ 6
#Hint: Of course, you first need to write a procedure AllKcomponentsGraphs(n,K) by tweaking procedure #AllConGraphs(n) and change the line nops(CCs(s))=1 to nops(CCs(s))=k

AllKComponentsGraphs:=proc(n,K) local S,s,C:
S:=AllGraphs(n):
C:={}:

for s in S do
 if nops(CCs(s))=K then
  C:=C union {s}:
 fi:
od:

C:
end:

NuKComponentsGraphs:=proc(N,K) local n1:
[seq(nops(AllKComponentsGraphs(n1,K)),n1=1..N)]:
end:

#Are the following sequences in the OEIS? If yes, what are the A-numbers

#NuConGgen(5,2), NuConGgen(5,3)

[0, 1, 3, 19, 230] Yes, A323875
[0, 0, 1, 6, 55] Yes, A323876

#2) 
#Using procedures RandGr(n,K) and CCs(G), in M18.txt

#write a procedure

#AveNuCC(n,K,M)

#that generates M random graphs on n vertices with K edges, and for each of them finds the number of 
#components and take the average over these M random graphs. Run it several times with n=30 and M=1000 
#to check that you get close answers. If not what M would give such agreement?

AveNuCC:=proc(n,K,M) local i, G, numComponents, average:
average:=0:
for i from 1 to M do
 G:=RandGr(n,K):
 numComponents:=nops(CCs(G)):
 average:=average+numComponents:
od:
average/M:
end:

#Running with AveNuCC(30,100,1000):
#509/500=1.018
#203/200=1.015
#126/125=1.008
#Yes, these all look really close.

#What about lower K?
#Running with AveNuCC(30,50,1000):
#891/500=1.782
#1803/1000=1.803
#181/100=1.81
#They are still close.

#3) Using the above,say, write a procedure
#EstimateCutOff(n,M)

#that inputs a positive integer n and outputs an estimate for the "threshold" for K, when most graphs 
#with K edges start to be connected (say when the average number of connected components is less than 
#1.05)

EstimateCutOff:=proc(n,M) local k, ave, beforeAve, start, fin:
start:=0:
fin:=binomial(30,2):
k:=ceil((start+fin)/2):
ave:=AveNuCC(n,k,M):
beforeAve:=AveNuCC(n,k-1,M):
while ave>=1.05 or beforeAve<1.05 and start < fin do:
 if ave>=1.05 then
  start:=k:
  k:=ceil((start+fin)/2):
 else 
  fin:=k:
  k:=ceil((start+fin)/2):
 fi:
 ave:=AveNuCC(n,k,M):
 beforeAve:=AveNuCC(n,k-1,M):
od:
k:
end:

#EstimateCutOff(30,1000)
#	82