> #Weiji Zheng, 10/11/2020, Assignment9
;
> #OK TO POST HOMEWORK
;
> 
;
> #Q1
;
> 
;
> #DiagSeq2(f,x,y,N): Given a rational function f of the variables x and y such that the denominator does not vanish at (0,0) outputs the list whose i-th entry is the coefficient of x^(i-1)*y^(i-1) in the 2-variable power-series expansion of f
;
> DiagSeq2:=proc(f,x,y,N) local i:
> 
> if subs({x=0,y=0},denom(f))=0 then
>  print(`The denominator of`, f, `should have a non-zero constant term `):
>  RETURN(FAIL):
> fi:
> 
> [seq(coeff(taylor(coeff(taylor(f,x=0,N+1),x,i),y=0,N+1),y,i),i=0..N)]:
> 
> end:
> #DiagWalks2D(S,N): Inputs  a set of "atomic steps" in the 2-dimenional square lattice, with positve steps, of the form [a,b]
> #For example for a forward moving Knight the set if {[2,1],[1,2]} for a forward moving King it  is {[1,0],[0,1],[1,1]} Outputs the first N+1 terms of the sequence
> #Number of walks from [0,0] to [n,n] using the atomic steps. 
;
> DiagWalks2D:=proc(S,N) local s,x,y:
> if not (type(S, set) and type(N,integer) and N>=0) then
>  print(`Bad input`):
>  RETURN(FAIL):
> fi:
> 
> DiagSeq2(1/(1-add(x^s[1]*y^s[2], s in S)),x,y,N  ) :
> end:
> #DiagSeq3(f,x,y,z,N): Given a rational function f of the variables x, y, and z, such that the denominator does not vanish at (0,0,0),outputs the list whose i-th entry is the coefficient of x^(i-1)*y^(i-1)*z^(i-1) in the 3-variable power-series expansion of f
;
> DiagSeq3:=proc(f,x,y,z,N) local i:
> 
> if subs({x=0,y=0,z=0},denom(f))=0 then
>  print(`The denominator of`, f, `should have a non-zero constant term `):
>  RETURN(FAIL):
> fi:
> 
> [seq(coeff(taylor(coeff(taylor(coeff(taylor(f,x=0,N+1),x,i),y=0,N+1),y,i),z=0,N+1),z,i),i=0..N)]:
> 
> end:
> #DiagWalks3D(S,N): Inputs  a set of "atomic steps" in the 2-dimenional square lattice, with positve steps, of the form [a,b,c]
> #Number of walks from [0,0,0] to [n,n,n] using the atomic steps.
;
> DiagWalks3D:=proc(S,N) local s,x,y,z:
> if not (type(S, set) and type(N,integer) and N>=0) then
>  print(`Bad input`):
>  RETURN(FAIL):
> fi:
> 
> DiagSeq3(1/(1-add(x^s[1]*y^s[2]*z^s[3], s in S)),x,y,z,N  ) :
> end:
> #1)
;
> f := 1/(1-x-7*y-7*x^3+5*y^3)
f := 1/(-7*x^3+5*y^3-x-7*y+1)
;
> DiagSeq2(f,x,y,50)
[1, 14, 294, 16374, 665350, 23650536, 941988474, 39113859456, 1588244047866, 64722442252520, 2675703077583984, 111082334912162964, 4618596160389759226, 192735096918227269196, 8069571430745268285240, 338607042626572315139988, 14236385679386983181084526, 599702721280980695042944356, 25304348920523585816931359976, 1069259822330280233475611896380, 45242280882256424868700738138080, 1916584443066697037915361241288800, 81279989102054920078053191882552220, 3450402613508075932982797286421529320, 146605818427540305329861221255058109150, 6234423888553177529339205049510415036748, 265324228336476876507185996544143602984284, 11299730491842235369282382982751938461870964, 481556991097713286466976871734386406598224544, 20535051636638344078332720035585666662971121560, 876182286565610418328312734287751589253587554684, 37404815671871856525583656127751307259551719301672, 1597644994688587355147626770506885227816600628188686, 68271497022277564797795782155825479015550505179200636, 2918722397548720963551209642884982408226712976381945100, 124832959048336808030543731932438532987035941771153182068, 5341197820756931889618125366069824319880625158547225286224, 228618633865923382510794839736180030598999874467981168648424, 9789027125357076074524405200997972936008909064347024251943204, 419289925364883353985801230106383468733392299299572510941673160, 17965071622312307833625649781325340441618274666784962981132122772, 769974405511114929934016537556468400064571564844912610703743620944, 33010348092987337031729997075564619288762783243353916843072240922304, 1415612887297408184212262991579818185037154253986678383521153187662584, 60723110857624565490888600853038736769353396856530934577405565776920820, 2605394528643749649306039515659542652782413512182708637934439637217145000, 111814512190402422385596718261204903741947590782375773701784694286339068200, 4799804638698727301406272506160479790128536230867368502772519390666170803560, 206084561627067581616074338345751559894385252750787860312192046630243265539610, 8850339343476969379450304976598339089935807403870426331301682494991354047340300, 380157223292136448996063169461959991237321785947444528460503383239410999897982900]
;
>  #The coefficient is 380157223292136448996063169461959991237321785947444528460503383239410999897982900
;
> #2)
;
> DiagWalks3D({[1,9,1], [7, 0, 0], [0,9,9]},31)
[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
;
> 
;
> #Q2
;
> 
;
> #Write a Maple procedure WordsMod(a,S,x) that inputs a positive integer a, subsest S of {1,2,..., a} and a variable x and outputs the rational function that is the generating function of finite sequence of positive integers that leave remainder that belongs to S when divided by a. Use the Maple "factor" to make it as nice as possible.
;
> 
;
> WordsMod := proc(a,S,x) local s, f1, f2:
> option remember:
> 
> f1 := add(x^s, s in S)/(1-x^a):
> f2 := normal(1/(1-f1)):
> 
> f2:
> end:
> 
;
> f := WordsMod(5,{1,2,3,4,5},x)
f := (x-1)/(2*x-1)
;
> seq(coeff(taylor(f,x=0,101),x,i),i=0..100)[101]
633825300114114700748351602688
;
> #There are 633825300114114700748351602688 ways
;
> 
;
> #Q3
;
> 
;
> ResComps:= proc(S,x) local s, f1, f2, f3:
> option remember:
> 
> f1 := add(x^s,s in S):
> f2 := x/(1-x):
> f3 := 1/(1-f2+f1):
> 
> f3:
> end:
> 
;
> f :=ResComps({2,3},x)
f := 1/(1-x/(1-x)+x^3+x^2)
;
> seq(coeff(taylor(f,x=0,301),x,i),i=0..300)[301]
179789662603291780480882832586094549526183607154633616127794901
;
> #number of walks := 179789662603291780480882832586094549526183607154633616127794901
;
> 
;
> #Q4
;
> 
;
> #1£©
;
> DiagWalks2D({[0,1],[1,0]},40)
[1, 2, 6, 20, 70, 252, 924, 3432, 12870, 48620, 184756, 705432, 2704156, 10400600, 40116600, 155117520, 601080390, 2333606220, 9075135300, 35345263800, 137846528820, 538257874440, 2104098963720, 8233430727600, 32247603683100, 126410606437752, 495918532948104, 1946939425648112, 7648690600760440, 30067266499541040, 118264581564861424, 465428353255261088, 1832624140942590534, 7219428434016265740, 28453041475240576740, 112186277816662845432, 442512540276836779204, 1746130564335626209832, 6892620648693261354600, 27217014869199032015600, 107507208733336176461620]
;
> gfun[listtorec]([1, 2, 6, 20, 70, 252, 924, 3432, 12870, 48620, 184756, 705432, 2704156, 10400600, 40116600, 155117520, 601080390, 2333606220, 9075135300, 35345263800, 137846528820, 538257874440, 2104098963720, 8233430727600, 32247603683100, 126410606437752, 495918532948104, 1946939425648112, 7648690600760440, 30067266499541040, 118264581564861424, 465428353255261088, 1832624140942590534, 7219428434016265740, 28453041475240576740, 112186277816662845432, 442512540276836779204, 1746130564335626209832, 6892620648693261354600, 27217014869199032015600, 107507208733336176461620],ff(n))
[{(-n-2)*ff(n+2)+(6+4*n)*ff(n+1), ff(0) = 1, ff(1) = 2}, ogf]
;
> #A000984
;
> 
;
> #2£©
;
> DiagWalks2D({[0,1],[1,0],[1,1]},40)
[1, 3, 13, 63, 321, 1683, 8989, 48639, 265729, 1462563, 8097453, 45046719, 251595969, 1409933619, 7923848253, 44642381823, 252055236609, 1425834724419, 8079317057869, 45849429914943, 260543813797441, 1482376214227923, 8443414161166173, 48141245001931263, 274738209148561921, 1569245074591690083, 8970232353223635949, 51313576749006450879, 293733710358893793729, 1682471873186160624243, 9642641465118083682429, 55294491352291112747007, 317241780630136241094657, 1820991621200098527441027, 10457362855894862001750669, 60078868458555230983696959, 345299757825442889707393857, 1985346154034162284608274707, 11419126147845924397833900957, 65701922725618214591910684159, 378150244155138145169182750209]
;
> gfun[listtorec]([1, 3, 13, 63, 321, 1683, 8989, 48639, 265729, 1462563, 8097453, 45046719, 251595969, 1409933619, 7923848253, 44642381823, 252055236609, 1425834724419, 8079317057869, 45849429914943, 260543813797441, 1482376214227923],ff(n))
[{(-n-3)*ff(n+3)+(15+6*n)*ff(n+2)+(-n-2)*ff(n+1), ff(0) = 1, ff(1) = 3, ff(2) = 13}, ogf]
;
> #A001850
;
> 
;
> #3)
;
> DiagWalks2D({[0,1],[1,0],[0,2],[2,0]},40)
[1, 2, 14, 84, 556, 3736, 25612, 177688, 1244398, 8777612, 62271384, 443847648, 3175924636, 22799963576, 164142004184, 1184574592592, 8567000931404, 62073936511496, 450518481039956, 3274628801768744, 23833760489660324, 173679413875623368, 1267013689048017584, 9252299435205985664, 67626504432024377756, 494710324956646794296, 3621791112234327295592, 26534383313499716907504, 194529413506239838951024, 1427026630616364232856416, 10474450985957996054150576, 76924819005707350396260192, 565226994478396276874021358, 4155153128570549776415060268, 30559526512591226909700271992, 224848729555165636463415722784, 1655038891063501813947452278468, 12186822164518997043970184661512, 89769447797439589097829820560392, 661476500151656115107632518838128, 4875741284021393250978952096137064]
;
> gfun[listtorec]([1, 2, 14, 84, 556, 3736, 25612, 177688, 1244398, 8777612, 62271384, 443847648, 3175924636, 22799963576, 164142004184, 1184574592592, 8567000931404, 62073936511496, 450518481039956, 3274628801768744, 23833760489660324],ff(n))
[{(352*n^3+2768*n^2+7088*n+5920)*ff(n+1)+(-484*n^3-4048*n^2-11184*n-10204)*ff(n+2)+(-352*n^3-3120*n^2-9114*n-8766)*ff(n+3)+(55*n^3+515*n^2+1570*n+1560)*ff(n+4), ff(0) = 1, ff(1) = 2, ff(2) = 14, ff(3) = 84}, ogf]
;
> #A036692
;
> 
;
> #4)
;
> DiagWalks2D({[0,1],[0,2],[1,0],[2,0],[1,1],[2,2]},40)
[1, 3, 22, 165, 1327, 10950, 92045, 783579, 6733966, 58294401, 507579829, 4440544722, 39000863629, 343677908223, 3037104558574, 26904952725061, 238854984979423, 2124492829796598, 18927927904130617, 168888613467092895, 1508973226894216106, 13498652154574126523, 120886709687492946083, 1083687170568092836350, 9723660300694365146989, 87322023363899721467343, 784797155029928933462614, 7058384549327774062817985, 63525027490219905080597647, 572078727948766828614679830, 5154896021503193637126695629, 46475144470335949808346622587, 419221380265186873217492981134, 3783331574173309950816296187025, 34158702327707968075053323650645, 308541107632605100591602320772210, 2788040608584603973862476702488265, 25202873888846971990310549726981859, 227906767043212151853904258820045110, 2061638450059846741700344736532466833, 18655567955631232569396589476504942877]
;
> gfun[listtorec]([1, 3, 22, 165, 1327, 10950, 92045, 783579, 6733966, 58294401, 507579829],ff(n))
> 

[{(33095335382797003565*n^2+244905287147564379923*n+361612491948565735644)*ff(n+1)+(-19561680998939410557*n^2-81380566882732553103*n-85543083305737218636)*ff(n+2)+(-28132660851759867329*n^2-185298876422503244507*n-301973593491614385996)*ff(n+3)+(3260032226768207496*n^2+22577196217358605284*n+38148269241143101200)*ff(n+4), ff(0) = 1, ff(1) = 3, ff(2) = 22, ff(3) = 165}, ogf]
;
> #A192365
;
> 
;
> #Q5
;
> 
;
> #1)
;
> DiagWalks3D({[0,0,1],[0,1,0],[0,0,1]},50)
[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
;
> gfun[listtorec]([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],ff(n))
[{-ff(n+1), ff(0) = 1}, ogf]
;
> #Sorry Dr.Z, my computer cannot calculate the below  DiagWalks3D(S,N)
;
> 
;
> #Q6
;
> #Just add one additional process to Diag3
;
> DiagSeq3:=proc(f,x,y,z,zz,N) local i:
> 
> if subs({x=0,y=0,z=0,zz=0},denom(f))=0 then
>  print(`The denominator of`, f, `should have a non-zero constant term `):
>  RETURN(FAIL):
> fi:
> 
> [seq(coeff(taylor(coeff(taylor(coeff(taylor(coeff(taylor(f,x=0,N+1),x,i),y=0,N+1),y,i),z=0,N+1),z,i),zz=0,N+1),zz,i),i=0..N)]:
> 
> end:
> DiagWalks4D:=proc(S,N) local s,x,y,z,zz:
> if not (type(S, set) and type(N,integer) and N>=0) then
>  print(`Bad input`):
>  RETURN(FAIL):
> fi:
> 
> DiagSeq4(1/(1-add(x^s[1]*y^s[2]*z^s[3]*zz^s[4], s in S)),x,y,z,zz,N) :
> end:
> 
;