#OK to post homework
#Blair Seidler, 2021-03-07 Assignment 13

with(combinat):

Help:=proc():
print(`pnE(n), pnSeqE(N), Franklin(L), FranklinOldMaids(n), SeqFOM(N)`):
end:

# 1.

#pnE(n): Uses the recurrence implied by Euler's Pentagonal Theorem to calculate p(n)
#p(n)=-Sum((-1)^j* p(n-((3*j-1)*j)/2),j=1..infintiy)-Sum((-1)^j* p(n-((3*j+1)*j)/2),j=1..infinity)):
pnE:=proc(n) local p,j:
option remember:
if n<2 then RETURN(1): fi:
if n=2 then RETURN(2): fi:
p:=0:
for j from 1 while n-((3*j-1)*j/2)>=0 do
  p:=p-(-1)^j*pnE(n-((3*j-1)*j)/2):
od:
for j from 1 while n-((3*j+1)*j/2)>=0 do
  p:=p-(-1)^j*pnE(n-((3*j+1)*j)/2):
od:
p:
end:

pnSeqE:=proc(N) local n: [seq(pnE(n),n=0..N)]:end:

#### Running times: restart, read this file, run one of pnSeq, pnSeqGF, psSeqE ####
#		N=1000		N=2000		N=50000		N=100000
# pnSeq		42.796s		429.421s	nope		nope
# pnSeqGF	 3.468s		 25.984s	nope		nope
# pnSeqE	 0.062s		  0.125s	17.812s		58.406s

#Yes, pnSeqE is a lot faster than either of the others.

# 2.

#Franklin(L): inputs a distinct integer partition (i.e. an integer partition with distinct
#parts), and outputs another one with one or less number of parts. It should return FAIL, if
#it not applicable.
#(Implementation via Fabian Franklin's beautiful proof of Euler's Pentagonal Theorem.)
Franklin:=proc(L) local i,diag,small,C:
small:=L[-1]:
diag:=1:
for i from 1 to nops(L)-1 do
  if L[i]=L[i+1]+1 then
    diag:=diag+1:
  else
    break:
  fi:
od:

if diag=nops(L) and ((diag=small) or (diag+1=small)) then
  RETURN(FAIL):
elif diag<small then
  C:=[seq(op(i,L)-1,i=1..diag),seq(op(i,L),i=diag+1..nops(L)),diag]:
else
  C:=[seq(op(i,L)+1,i=1..small),seq(op(i,L),i=small+1..nops(L)-1)]:
fi:

C:

end:

# 3.
#FranklinOldMaids(n)
FranklinOldMaids:=proc(n) local PL,p,S:
PL:=PnD(n,{0}):
S:={}:
for p in PL do
if Franklin(p)=FAIL then
  S:=S union {p}:
fi:
od:
S:
end:

SeqFOM:=proc(N): [seq(FranklinOldMaids(n),n=1..N)]: end:


#### From Homework 12 ####
#PnkD(n,k,DI): The set of partitions of n with largest part k where the 
#difference of two consecutive parts is NOT in DI
PnkD:=proc(n,k,DI) local k1,L,L1:
option remember:

if not (type(n,integer) and type(k,integer) and n>=1 and k>=1 )then
 RETURN([]):
fi:
 
if k>n then
 RETURN([]):
fi:

if k=n then
 RETURN([[n]]):
fi:

L:=[]:

for k1 from min(n-k,k) by -1 to 1 do
 if not(k-k1 in DI) then
   L1:=PnkD(n-k,k1,DI):
   L:=[op(L), seq([k, op(L1[j])],j=1..nops(L1))]:
  fi:
od:

L:

end:

#PnD(n,DI): The set of partitions of n where all parts are congruent to a member of S mod m
PnD:=proc(n,DI) local k:option remember:[seq(op(PnkD(n,n-k+1,DI)),k=1..n)]:end:

#### From C12.txt (my version) ####
Ferr:=proc(L) local i,j:
for i from 1 to nops(L) do
  for j from 1 to L[i] do
    printf("*"):
  od:
  printf("\n"):
od:
end:

#### From C11.txt ####
#pnk(n,k):
pnk:=proc(n,k) local k1:
option remember:
if not (type(n,integer) and type(k,integer) and n>=1 and k>=1 )then
 RETURN(0):
fi:
if n=k then
 RETURN(1):
else
 RETURN(add(pnk(n-k,k1),k1=1..min(k,n-k))):
fi:
end:

pn:=proc(n) local k: option remember: 
if n=0 then
 RETURN(1):
fi:

add(pnk(n,k),k=1..n):end:

pnSeq:=proc(N) local n: [seq(pn(n),n=0..N)]:end:

#pnSeqGF(N):
pnSeqGF:=proc(N) local i,q,f:
f:=taylor(mul(1/(1-q^i),i=1..N),q=0, N+1 ):
[seq(coeff(f,q,i),i=0..N)]:
end: