# OK to post homework
# Lucy Martinez, 01-31-2025, Assignment 3

with(combinat):

# Problem 4:Convince yourself that the following one-line procedure,
#  for p rational between 0 and 1 and n positive integer is an estimate,
#  with large enough K (say at least 200), for the probability of
#  connectivity of a random graph with n edges where every edge
#  is picked with prob. p (the Edros-Renyi model)
PC:=proc(n,p,K) local i:
 evalf(coeff(add(IsCo(RG(n,p)),i=1..K),true,1)/K):
end:

# With n=10,20,30, and p=1/10, 2/10, ... and K=200 experiment
#  and try to estimate the threshhold for connectivity,
#  i.e. when it changes being very unlikely to very likely.
# Answer: Threshold should be around p=3/10
# PC(10,1/10,200)=0.01000000000
# PC(20,1/10,200)=0.04000000000
# PC(30,1/10,200)=0.2450000000

# PC(10,2/10,200)=0.2150000000
# PC(20,2/10,200)=0.8000000000
# PC(30,2/10,200)=0.9650000000

# PC(10,3/10,200)=0.6550000000
# PC(20,3/10,200)=0.9800000000
# PC(30,3/10,200)=0.9950000000

# PC(10,4/10,200)=0.9150000000
# PC(20,4/10,200)=1
# PC(30,4/10,200)=1

# PC(10,5/10,200)=0.9800000000
# PC(20,5/10,200)=1
# PC(30,5/10,200)=1

# PC(10,6/10,200)=0.9950000000
# PC(20,6/10,200)=1
# PC(30,6/10,200)=1

# PC(10,7/10,200)=1




######################


# Problem 5: Procedure NuCCc(n) outputs the first n terms (starting at n=1)
#  of OEIS sequence A1187 that counts the sequence let's call it C[n],
#  of the number of connected labeled graphs with n vertices.
#  It is based on the formula
#  Sum(C[n]*z^n/n!,n=0..infinity)=log(Sum(2^(n*(n-1)/2)*z^n/n!,n=0..infinity)
#  Let P[n](x) be the polynomial in x whose coefficient of x^r
#  is the number of labeled connected graphs on n vertices
#  and EXACTLY r edges. Convince yourself (or believe) that
#  using generating functions one has the formula
#  Sum(P[n](x)*z^n/n!,n=0..infinity)=log(Sum((1+x)^(n*(n-1)/2)*z^n/n!,n=0..infinity)
#  Write a procedure NuCGx(n,x) that inputs a pos. integer n
#  and a variable x, and outputs a list of length n
#  whose i-th entry is P[i](x)
NuCGx:=proc(n,x) local f,L,i,P:
f:=log(sum((1+x)^binomial(i,2)*z^i/i!,i=0..n)):
f:=taylor(f,z=0,n+1):
P:=n!*coeff(f,z,n):
L:=[seq(coeff(P,x,i),i=1..degree(P))]:
end:


######################


# Problem 6: Convince yourself that the lowest power (lowdegree in Maple)
#  of P[n](x) is n-1. Using NuCGx(n,x)
#  (make sure to put "option remember" there),
#  Write a procedure NuCGk(n,k) that inputs a positive integer n
#  and a non-negative integer k and outputs the sequence of coefficients
#  of x^(n-1+k) of P[n](x)

#  What OEIS number is NuCGk(n,0)?
#  What OEIS number is NuCGk(n,1)? (if it exists)
#  Keep upping k and see the largest k for which NuCGk(n,k) is in the OEIS
NuCGk:=proc(n,k) local L,x:
L:=NuCGx(n,x):

end:


######################Procedures from C3.txt
# Neis(G): inputs a graph G=[n,E] and outputs a list N of lenght n such that
#  N[i] is the set of neighbors of vertex i
Neis:=proc(G) local n, E,N,i,e:
n:=G[1]:
E:=G[2]:

for i from 1 to n do
 N[i]:={}:
od:

for e in E do
 N[e[1]]:=N[e[1]] union {e[2]}:
 N[e[2]]:=N[e[2]] union {e[1]}:
od:
[seq(N[i],i=1..n)]:
end:


#CC(G,i): The connected component of vertex i
CC:=proc(G,i) local n,N,C1,C2,C3,i1:

if not (type(G,list) and nops(G)=2 and type(G[1],posint) and type(G[2],set)) then
 return(FAIL):
fi:

n:=G[1]:
N:=Neis(G):
#C1 is the current friend
C1:={i}:
#C2 is the current friends
C2:=C1 union {seq(op(N[i1]),i1 in C1)}:

while C1<>C2 do
 #C3 is the friends of friends
 C3:=C2 union {seq(op(N[i1]),i1 in C2)}:
 C1:=C2:
 C2:=C3:
od:
C2:
end:



#IsCo(G): is the graph G=[n,E] connected?
IsCo:=proc(G) local n,i,V:
n:=G[1]:
V:={seq(i,i=1..n)}:
evalb(CC(G,1)=V):
end:


#LC(p): inputs a rational number between 0 and 1 and outputs true with prob. p
LC:=proc(p) local a,b,ra:
if not (type(p,fraction) and p>= 0 and p<=1) then
 return fail:
fi:

a:=numer(p):
b:=denom(p):

ra:=rand(1..b)():

if ra<=a then
 true:
else
 false:
fi:
end:


RG:=proc(n,p) local E,i,j:
E:={}:
for i from 1 to n do
 for j from i+1 to n do
  if LC(p) then
   E:=E union {{i,j}}:
  fi:
 od:
od:
[n,E]:
end: