# NOT OK to post homework
# Michael Saunders, 10/4/2020, Homework 7
#
# load packages
#
with(combinat):


# 1.
# Examine the Maple code of M7... especially the procedures not discussed in
# lecture, GFseq(f, x, N):
#
# GFseq(f,x,N): Given a rational function f in the variable x, (whose denominator does not vanish at 0) [or a more general function/expression)
# returns the first N+1 terms, starting at n=0 of the coefficients of the sequence that has it has a GENERATING FUNCTION
# Equivalently, returns the first n+1 terms. stating at n=0 of its Taylor coefficients. Try:
# GFseq(x/(1-x-x^2),x,100);
#
GFseq:=proc(f,x,N) local  f1,i:
		f1:=taylor(f,x=0,N+1):
		[seq(coeff(f1,x,i),i=0..N)]:
end:
#
# GFseq(f,x,N) first uses taylor(f,x=0,N+1) to compute the Taylor series
# expansion of the expression f--of the order N+1--with respect to the variable x.
# The expansion of the Taylor series is assigned to variable f1.
# Then, [seq(coeff(f1,x,i), i=0..N)] builds a list of the sequence created by
# extracting the coefficients of x^i in the polynomial f1, for i=0..N.
#
# Compute the following...
#
GFseq(1/(1-x-x^2-x^3),x,30);
# [1, 1, 2, 4, 7, 13, 24, 44, 81, 149, 274, 504, 927, 1705, 3136, 5768, 10609, 
# 19513, 35890, 66012, 121415, 223317, 410744, 755476, 1389537, 2555757, 4700770,
# 8646064, 15902591, 29249425, 53798080]
#
GFseq(1/(1-x+x^2-x^3),x,30);
# [1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1,
# 0, 0, 1, 1, 0]
#
GFseq(exp(-x)/(1-x),x,30);
# [1, 0, 1/2, 1/3, 3/8, 11/30, 53/144, 103/280, 2119/5760, 16687/45360, 16481/
# 44800, 1468457/3991680, 16019531/43545600, 63633137/172972800, 2467007773/
# 6706022400, 34361893981/93405312000, 15549624751/42268262400, 8178130767479/
# 22230464256000, 138547156531409/376610217984000, 92079694567171/250298560512000
# , 4282366656425369/11640679464960000, 72289643288657479/196503623737344000, 
# 6563440628747948887/17841281393295360000, 39299278806015611311/
# 106826515449937920000, 9923922230666898717143/26976017466662584320000, 
# 79253545592131482810517/215433472824041472000000, 5934505493938805432851513/
# 16131658445064225423360000, 14006262966463963871240459/
# 38072970106357874688000000, 461572649528573755888451011/
# 1254684545727217532928000000, 116167945043852116348068366947/
# 315777214062132212662272000000, 3364864615063302680426807870189/
# 9146650338351415815045120000000]


# 2.
# Are GFseq(x/1-x-x^2),x,30) and SeqRex([1,1],[0,1],30) identical? 
# Why or why not?
#
# First, include procedures Rec(P,n) and SeqRec(P,N) from
# sites.math.rutgers.edu/~zeilberger/Combo20/M7.txt
#
# Rec(P,n): inputs a pair P=[L,INI], where   L=[a1,a2,..,ak] (where k=nops(L)), of numbers and another list INI=[e_0,...,e_(k-1)] of numbers of the SAME size
# and a positive ineteger n, outputs the n-th term x(n) of the sequence defined by the recurrence
# x(n)=a1*x(n-1)+ ...+ ak*x(n-k), for n>=k with the INITIAL CONDITIONS x(i)=e_{i} for i=0..k-1. For example to get
# the 100-th Fibonacci number, try:
# Rec([1,1],[0,1],100);
#
Rec:=proc(P,n)  local L, INI,k,i:
		option remember:
		# We first unpack P
		L := P[1]:
		INI := P[2]:
		k := nops(L):
		if nops(INI)<>k then
				print(INI, `should have`, k, `entries `):
				RETURN(FAIL):
		fi:
		if n<k then
				# For this case we use the initial condition, since list in Maple start at index 1, we need to access INI[n+1]
				RETURN(INI[n+1]):
		else
		#We use the definig recurrence
				add(L[i]*Rec(P,n-i),i=1..k):
		fi:
end:
#
# SeqRec(P,N): The first N+1 terms (starting at n=0) of the sequence Rec(P,n)
#
SeqRec:=proc(P,N) local n:
		[seq(Rec(P,n),n=0..N)]:
end:
#
evalb(GFseq(x/(1-x-x^2),x,30) = SeqRec([[1,1],[0,1]],30));
# true
#
# I could write more... but essentially, SeqRec([[1,1],[0,1]],30) builds the
# first 30 terms of the Fibonnaci sequence using the Reccurence Relatioin defined
# in the procedure, and GFseq(x/(1-x-x^2),x,30) does the same thing using 
# the Ordinary Generating Function for Fibonnaci sequence.


# 3.
# 
# In the lecture we started with a sequence that is known to satisfy a linear
# recurrence equation with constant coefficients and showed how to derive an
# explicit expression for the generating function
#
# Sum(a(n)*x^n, n=0..infinity) = 1/(1-3*x+5*x^2)
#
# as a certain rational function whose denominator ONLY depends on the
# coefficients of the reccurende (c1, c2, ..., ck) and whose numerator also
# depends on the initical conditions. Show how to go the other way.
#
# I have no idea what's going on... I'm having a hard time following the
# lectures. I am more of a "learn-by-reading" type... I've never learned much from other people
# telling/showing me something in a lecture, unless they are following the
# structure of a text book that I can follow along with.
#
# In the Homework description, you say "In other words a(n) is the coefficient
# of x^n in the Taylor expansion (around x=0) of 1/(1-3*x+5*x^2)
#
# [In Maple, you can do coeff(taylor(1/(1-3*x+5*x^2),x=0,n+1),x,n) to compute
# it.] "
#
# However, there seems to be an error in the Maple code you gave us, and I can't
# be bothered to continue to fix it, while also trying to figure out what is
# happening. 
#
coeff(taylor(1/(1-3*x+5*x^2),x=0,n+1),x,n);
# Error, invalid input: taylor expects its 3rd argument, n, to be of type
# {infinity, nonnegint}, but received n+1
#
# Nope. I'm over it.