#OK to post homework
#William Wang, 9/23/2020, Assignment 6

#1.
#PIE is only of theoretical interest because as the number of subsets increases, the time it takes to calculate PIE increases tremendously. `^()Therefore at larger values of n, it is more efficient to calculate the union of a set of subsets by looking at each one by one.

with(combinat):

FP := proc(pi) local i:
coeff(add(evalb(pi[i] = i), i = 1 .. nops(pi)), true, 1):
end:

Pnx := proc(n, x) local S, pi:
S := permute(n):
add(x^FP(pi), pi in S):
end:

PnxF := proc(n, x) local k:
expand(add(binomial(n, k)*(x - 1)^k*(n - k)!, k = 0 .. n)):
end:

[seq(Pnx(n, x), n = 1 .. 7)];
    [    2       3             4      2            
    [x, x  + 1, x  + 3 x + 2, x  + 6 x  + 8 x + 9, 

       5       3       2              
      x  + 10 x  + 20 x  + 45 x + 44, 

       6       4       3        2                
      x  + 15 x  + 40 x  + 135 x  + 264 x + 265, 

       7       5       4        3        2                ]
      x  + 21 x  + 70 x  + 315 x  + 924 x  + 1855 x + 1854]


[seq(PnxF(n, x), n = 1 .. 7)];
    [    2       3             4      2            
    [x, x  + 1, x  + 3 x + 2, x  + 6 x  + 8 x + 9, 

       5       3       2              
      x  + 10 x  + 20 x  + 45 x + 44, 

       6       4       3        2                
      x  + 15 x  + 40 x  + 135 x  + 264 x + 265, 

       7       5       4        3        2                ]
      x  + 21 x  + 70 x  + 315 x  + 924 x  + 1855 x + 1854]


time([seq(Pnx(n, x), n = 1 .. 8)]);
                             0.687

time([seq(PnxF(n, x), n = 1 .. 8)]);
                               0.

time([seq(PnxF(n, x), n = 1 .. 30)]);
                             0.015

#time([seq(Pnx(n,x),n=1..30)]) will take a much longer time*(>0.015) since Pnx has to figure out all of the permutations, whereas PnxF uses PIE to efficiently find the coefficients of each term x^(i)





#2.
#In PnxF(n,x), the coefficient of x`^(0) for any polynomial is P(0). However, any value of x raised to the 0^(th) power just equals 1, and multipled by its coefficient is just its coefficient. Therefore, the coefficient of x^(0) for any polynomial in x is simply P(0).





#3.
[seq(coeff(PnxF(n, x), x, 0), n = 0 .. 12)];
[1, 0, 1, 2, 9, 44, 265, 1854, 14833, 133496, 1334961, 14684570, 

  176214841]

#The above sequence is A000166 in the OEIS


[seq(coeff(PnxF(n, x), x, 1), n = 0 .. 12)];
[0, 1, 0, 3, 8, 45, 264, 1855, 14832, 133497, 1334960, 14684571, 

  176214840]



[seq(coeff(PnxF(n, x), x, 2), n = 0 .. 12)];
  [0, 0, 1, 0, 6, 20, 135, 924, 7420, 66744, 667485, 7342280, 

    88107426]

#The above sequence is A000387 in the OEIS
#The largest i for which there is a sequence in the OEIS: i = 2





#4.
SumTools[Hypergeometric][ZeilbergerRecurrence](binomial(n, k)*(-1)^k*(n - k)!, n, k, d, 0 .. n);
                                               n
               (-n - 1) d(n) + d(n + 1) = -(-1) 

#Rearrange as d(n+1)-(n+1)d(n) = -(-1)^n
#Plug in n-1 for n to get an equivalent equation --> d(n) - n*d(n-1)= -(-1)^(n - 1)
#Plug in for n-1 again --> d(n-1)-(n-1)*d(n-2)= -(-1)^(n - 2)
#Add the two above equations --> d(n)+d(n-1)-n*d(n-1)-(n-1)*d(n-2)=0
#Simplify and get the other recurrence

#d(n)-(n-1)*d(n-1)-(n-1)*d(n-2)=0