###################################################################### ##WALKSab: Save this file as WALKSab # ## To use it, stay in the # ##same directory, get into Maple (by typing: maple ) # ##and then type: read WALKSab # ##Then follow the instructions given there # ## # ##Written by Doron Zeilberger, Rutgers University , # #zeilberg at math dot rutgers dot edu # ###################################################################### #Created: June 3, 2009 print(`Created: June 3, 2009`): print(` This is WALKSab `): print(`That automatically generates, and solves`): print(`the algebraic equations for the generating functions`): print(`let's call them F[0,i](x[a],x[-b]), `): print(`for i=0,1, ..., b-1 `): print(`for walks in the upper plane starting at the origin`): print(`with steps (1,a) and (1,-b), such that their graphs`): print(`are always weakly above the origin and end-up`): print(`at a point where y=i (and x anything)`): print(``): print(`It is one of the packages that accompany the article `): print(` "An Experimental Mathematics Perspective of the Old, and Still open Question, of when to Stop" `): print(`by Luis Medina and and Doron Zeilberger,`): print(`available from both authors' websites.`): print(``): print(`Please report bugs to zeilberg at math dot rutgers dot edu`): print(``): print(`The most current version of this package and paper`): print(` are available from`): print(`http://www.math.rutgers.edu/~zeilberg/ .`): print(`For a list of the procedures type ezra();, for help with`): print(`a specific procedure, type ezra(procedure_name); .`): print(``): with(combinat): ezra1:=proc() if args=NULL then print(` The supporting procedures are: `): print(` CheckGFs, Fab, F1ab,OneStep, Tay `): else ezra(args): fi: end: ezra:=proc() if args=NULL then print(`The main procedures are: Eqs, GFs, GFsAll, GFspq, Tay1,Tay1d `): elif nops([args])=1 and op(1,[args])=CheckGFs then print(`CheckGFs(a,b,N): checks GFs against Tay `): print(`For example, try:`): print(`CheckGFs(1,1,20);`): print(`CheckGFs(3,2,20);`): elif nops([args])=1 and op(1,[args])=Eqs then print(`Eqs(x,a,b,F): The set of equations expressed in terms`): print(`of pairs [LHS, RHS] for walks in the first quadrant`): print(`with steps [1,a], [1,-b] always staying weakly above the x-axis`): print(`F[i,j](resp. G[i,j]) denote walks that start at [0,i] `): print(`and end at the line y=i that are weakly (resp. strictly) `): print(`above the x-axis.`): print(`For example, try:`): print(`Eqs(x,1,1,F);`): print(`Eqs(x,2,3,F);`): elif nops([args])=1 and op(1,[args])=Fab then print(`Fab(a,b,n,i): the number of walks with steps (1,a) (1,-b) from`): print(`the origin to (n,i) staying weakly above the origin.`): print(`For example, try:`): print(`Fab(1,1,10,0);`): elif nops([args])=1 and op(1,[args])=F1ab then print(`F1ab(a,b,n): the number of walks with steps (1,a) (1,-b) from`): print(`the origin with n steps staying weakly above the origin.`): print(`and at the last step, going down`): print(`For example, try:`): print(`F1ab(1,1,10);`): elif nops([args])=1 and op(1,[args])=GFs then print(`GFs(x,a,b): The generating functions `): print(`F[0,j](j=0..b-1) for walks `): print(`with steps [1,a], [1,-b] always staying weakly above the x-axis`): print(`and ending at a point where y=j.`): print(`It also returns, as the second component, the full generating`): print(`function in x[a], x[-b] for all `): print(` walks that are always weakly above the x-axis`): print(` but at the very last step go below the x-axis and stop.`): print(`For example, try:`): print(`GFs(x,1,1);`): print(`GFs(x,2,3);`): elif nops([args])=1 and op(1,[args])=GFsAll then print(`GFsAll(x,a,b): All The generating functions `): print(`F[i,j] for walks `): print(`with steps [1,a], [1,-b] always staying weakly above the x-axis`): print(`and ending at a point where y=j.`): print(`It also returns, as the second component, the full generating`): print(`function in x[a], x[-b] for all `): print(` walks that are always weakly above the x-axis`): print(` but at the very last step go below the x-axis and stop.`): print(`For example, try:`): print(`GFsAll(x,1,1);`): print(`GFsAll(x,2,3);`): elif nops([args])=1 and op(1,[args])=GFspq then print(`GFspq(a,b,p,q): The probability of eventually exiting the`): print(`positive real axis if you start at the origin`): print(`with steps a (prob. p) and -b (prob. q). For example, try: `): print(`GFspq(1,1,1/2,1/2);`): print(`GFspq(2,3,1/2,1/2);`): elif nops([args])=1 and op(1,[args])=OneStep then print(`OneStep(S,t,n): given a triple S=[Trans,Var,Current]`): print(`a variable t, and a current subsitution with Taylor`): print(`series up to t^n, upgrades it to t^(n+1) . `): print(`For example, try:`): print(`OneStep(Eqst(x,1,1,F,t),t,1);`): elif nops([args])=1 and op(1,[args])=Tay then print(`Tay(x,a,b,t,N): The Taylor expansion up to order N`): print(`of the GFs(x,a,b)`): print(`For example, try:`): print(`Tay(x,1,1,t,10);`): elif nops([args])=1 and op(1,[args])=Tay1 then print(`Tay1(a,b,t,N): The Taylor expansion up to order N`): print(`of the generating function for walks`): print(`that are above the x-axis until the last step, and then `): print(`go below and quit, with steps (1,a) and (1,-b). `): print(`This is done via the algebraic system`): print(`For example, try:`): print(`Tay1(1,1,t,10);`): elif nops([args])=1 and op(1,[args])=Tay1d then print(`Tay1d(a,b,t,N): The Taylor expansion up to order N`): print(`of the generating function for walks`): print(`that are above the x-axis until the last step, and then `): print(`go below and quit, with steps (1,a) and (1,-b). `): print(`This is done directly`): print(`For example, try:`): print(`Tay1d(1,1,t,10);`): else print(`There is no ezra for`,args): fi: end: #EqsOld(x,a,b,F,G): The set of equations expressed in terms #of pairs [LHS, RHS] for walks in the first quadrant #with steps [1,a], [1,-b] always staying weakly above the x-axis #F[i,j](resp. G[i,j]) denote walks that start at [0,i] #and end at the line y=i that are weakly (resp. strictly) above the x-axis. #For example, try: #Eqs(x,1,1,F,G); #Eqs(x,2,3,F,G); EqsOld:=proc(x,a,b,F,G) local eq,AlreadyDone,StillToDo,pt,i,j: eq:={}: StillToDo:={F[0,0]}: AlreadyDone:={}: while StillToDo<>{} do pt:=StillToDo[1]: i:=op(1,pt): j:=op(2,pt): if op(0,pt)=F then if i=0 and j=0 then AlreadyDone:= AlreadyDone union {F[0,0]}: eq:=eq union {[F[0,0],1+G[0,0]*F[0,0]]}: StillToDo:=StillToDo union {G[0,0],F[0,0]} minus AlreadyDone: elif i>0 and j=0 then AlreadyDone:= AlreadyDone union {F[i,0]}: eq:=eq union {[F[i,0],G[i,0]*F[0,0]]}: StillToDo:=StillToDo union {G[i,0],F[0,0]} minus AlreadyDone: elif i=0 and j>0 then AlreadyDone:= AlreadyDone union {F[0,j]}: eq:=eq union {[F[0,j],F[0,0]*G[0,j]]}: StillToDo:=StillToDo union {F[0,0],G[0,j]} minus AlreadyDone: elif i>0 and j>0 then AlreadyDone:= AlreadyDone union {F[i,j]}: eq:=eq union {[F[i,j],G[i,j]+G[i,0]*F[0,0]*G[0,j]]}: StillToDo:=StillToDo union {G[i,j],G[i,0],F[0,0],G[0,j]} minus AlreadyDone: else ERROR(`Something is wrong`): fi: elif op(0,pt)=G then if i=0 and j=0 then AlreadyDone:= AlreadyDone union {G[0,0]}: eq:=eq union {[G[0,0],x[a]*G[a,b]*x[-b]]}: StillToDo:=StillToDo union {G[a,b]} minus AlreadyDone: elif i>0 and j=0 then AlreadyDone:= AlreadyDone union {G[i,0]}: eq:=eq union {[G[i,0],G[i,b]*x[-b] ]}: StillToDo:=StillToDo union {G[i,b]} minus AlreadyDone: elif i=0 and j>0 then AlreadyDone:= AlreadyDone union {G[0,j]}: eq:=eq union {[G[0,j],x[a]*G[a,j]]}: StillToDo:=StillToDo union {G[a,j]} minus AlreadyDone: elif i>0 and j>0 then AlreadyDone:= AlreadyDone union {G[i,j]}: eq:=eq union {[G[i,j],F[i-1,j-1]]}: StillToDo:=StillToDo union {F[i-1,j-1]} minus AlreadyDone: else ERROR(`Something is wrong`): fi: else ERROR(`Something is wrong`): fi: od: eq,AlreadyDone: end: #GFs(x,a,b): The generatinf functions #F[0,j](j=0..b-1) for walks #with steps [1,a], [1,-b] always staying weakly above the x-axis #and ending at a point where y=j #For example, try: #GFs(x,1,1); #GFsOld(x,2,3); GFsOld:=proc(x,a,b) local F,G,gu,eq,var,i,j: gu:=EqsOld(x,a,b,F,G): var:=gu[2]: eq:=gu[1]: eq:={seq(eq[i][1]-eq[i][2],i=1..nops(eq))}: var:=solve(eq,var): [seq(subs(var,F[0,j]),j=0..b-1)]: end: #Eqs(x,a,b,F): The set of equations expressed in terms #of pairs [LHS, RHS] for walks in the first quadrant #with steps [1,a], [1,-b] always staying weakly above the x-axis #F[i,j] denote walks that start at [0,i] #and end at the line y=i that are weakly (resp. strictly) above the x-axis. #it also returns the initial conditions #For example, try: #Eqs(x,1,1,F); #Eqs(x,2,3,F); Eqs:=proc(x,a,b,F) local eq,AlreadyDone,StillToDo,pt,i,j,Ini,a1,i1,j1: eq:={}: StillToDo:={F[0,0]}: AlreadyDone:={}: while StillToDo<>{} do pt:=StillToDo[1]: i:=op(1,pt): j:=op(2,pt): if i=0 and j=0 then AlreadyDone:= AlreadyDone union {F[0,0]}: eq:=eq union {[F[0,0],1+x[a]*F[a-1,b-1]*x[-b]*F[0,0]]}: StillToDo:=StillToDo union {F[a-1,b-1],F[0,0]} minus AlreadyDone: elif i>0 and j=0 then AlreadyDone:= AlreadyDone union {F[i,0]}: eq:=eq union {[F[i,0],F[i-1,b-1]*x[-b]*F[0,0]]}: StillToDo:=StillToDo union {F[i-1,b-1],F[0,0]} minus AlreadyDone: elif i=0 and j>0 then AlreadyDone:= AlreadyDone union {F[0,j]}: eq:=eq union {[F[0,j],F[0,0]*x[a]*F[a-1,j-1]]}: StillToDo:=StillToDo union {F[0,0],F[a-1,j-1]} minus AlreadyDone: elif i>0 and j>0 then AlreadyDone:= AlreadyDone union {F[i,j]}: eq:=eq union {[F[i,j],F[i-1,j-1]+F[i-1,b-1]*x[-b]*F[0,0]*x[a]*F[a-1,j-1]]}: StillToDo:=StillToDo union {F[i-1,j-1],F[i-1,b-1],F[0,0],F[a-1,j-1]} minus AlreadyDone: else ERROR(`Something is wrong`): fi: od: AlreadyDone:=convert(AlreadyDone,list): Ini:={}: for a1 in AlreadyDone do i1:=op(1,a1): j1:=op(2,a1): if i1=j1 then Ini:=Ini union {a1=1}: else Ini:=Ini union {a1=0}: fi: od: eq,AlreadyDone,Ini: end: #GFs(x,a,b): The generating functions #F[0,j](j=0..b-1) for walks #with steps [1,a], [1,-b] always staying weakly above the x-axis #and ending at a point where y=j #For example, try: #GFs(x,1,1); #GFs(x,2,3); GFs:=proc(x,a,b) local F,gu,eq,var,i,j: gu:=Eqs(x,a,b,F): var:=convert(gu[2],set): eq:=gu[1]: eq:={seq(eq[i][1]-eq[i][2],i=1..nops(eq))}: var:=solve(eq,var): [seq(subs(var,F[0,j]),j=0..b-1)], normal(x[-b]*add(subs(var,F[0,j]),j=0..b-1)): end: #Fab(a,b,n,i): the number of walks with steps (1,a) (1,-b) from #the origin to (n,i) staying weakly above the origin #For example, try: #Fab(1,1,10,0); Fab:=proc(a,b,n,i) option remember: if n<0 or i<0 then RETURN(0): fi: if n=0 then if i=0 then RETURN(1): else RETURN(0): fi: fi: Fab(a,b,n-1,i-a)+Fab(a,b,n-1,i+b): end: #GFspq(a,b,p,q): The probability of eventually exiting the #positive real axis if you start at the origin #with steps a (prob. p) and -b (prob. q) #GFspq(1,1,1/2,1/2); #GFspq(2,3,1/2,1/2); GFspq:=proc(a,b,p,q) local F,gu,eq,var,i,j,x,v: gu:=Eqs(x,a,b,F): var:=gu[2]: eq:=gu[1]: eq:=subs({x[a]=p,x[-b]=q},eq): eq:={seq(eq[i][1]-eq[i][2],i=1..nops(eq))}: var:={solve(eq,var)}: {seq(q*add(subs(v,F[0,j]),j=0..b-1),v in var)}: end: #GFsAll(x,a,b): The generating functions #F[0,j](j=0..b-1) for walks #with steps [1,a], [1,-b] always staying weakly above the x-axis #and ending at a point where y=j #For example, try: #GFsAll(x,1,1); #GFsAll(x,2,3); GFsAll:=proc(x,a,b) local F,gu,eq,var,i,v,var1: gu:=Eqs(x,a,b,F): var:=gu[2]: eq:=gu[1]: eq:={seq(eq[i][1]-eq[i][2],i=1..nops(eq))}: var1:=solve(eq,var): {seq([v,subs(var1,v)],v in var)}: end: #Eqst(x,a,b,F,t): The set of equations expressed in terms #of pairs [LHS, RHS] for walks in the first quadrant #with steps [1,a], [1,-b] always staying weakly above the x-axis #F[i,j] denote walks that start at [0,i] #and end at the line y=i that are weakly (resp. strictly) above the x-axis. #it also returns the initial conditions #For example, try: #Eqst(x,1,1,F,t); #Eqst(x,2,3,F,t); Eqst:=proc(x,a,b,F,t) local gu: gu:=Eqs(x,a,b,F): [subs({x[a]=t*x[a],x[-b]=t*x[-b]},gu[1]),gu[2],gu[3]]: end: #OneStep(S,t,n): given a triple S=[Trans,Var,Current] #a variable t, and a current subsitution with Taylor #series up to t^n, upgrades it to t^(n+1) #For example, try: #OneStep(Eqst(x,1,1,F,t),t,1); OneStep:=proc(S,t,n) local gu1,gu2,gu3,g1,kha,i,kha1: gu1:=S[1]: gu2:=S[2]: gu3:=S[3]: kha:={}: for g1 in gu1 do kha1:=expand(subs(gu3,g1[2])): kha1:=add(expand(coeff(kha1,t,i))*t^i,i=0..n): kha:=kha union {g1[1]=kha1}: od: [gu1,gu2,kha]: end: #Tay(x,a,b,t,N): The Taylor expansion up to order N #of the GFs(x,a,b) #For example, try: #Tay(x,1,1,t,10); Tay:=proc(x,a,b,t,N) local S,F,n,j: S:=Eqst(x,a,b,F,t): for n from 1 to N do S:=OneStep(S,t,n) : od: [[seq(expand(subs(S[3],F[0,j])),j=0..b-1)], expand(t*x[-b]*add(subs(S[3],F[0,j]),j=0..b-1))]: end: #CheckGFs(a,b,N): cheks GFs CheckGFs:=proc(a,b,N) local gu,i1,j,x,t,mu: gu:=GFs(x,a,b) : gu:=subs({x[a]=t*x[a],x[-b]=t*x[-b]},gu[1]), subs({x[a]=t*x[a],x[-b]=t*x[-b]},gu[2]): gu:=[[seq(taylor(gu[1][j],t=0,N+2),j=1..nops(gu[1]))],taylor(gu[2],t=0,N+2)]: gu:=[[seq(add(coeff(gu[1][j],t,i1)*t^i1,i1=0..N),j=1..nops(gu[1]) )], add(coeff(gu[2],t,i1)*t^i1,i1=0..N)]: mu:=Tay(x,a,b,t,N+2): mu:=[[seq(add(coeff(mu[1][j],t,i1)*t^i1,i1=0..N),j=1..nops(mu[1]) )], add(coeff(mu[2],t,i1)*t^i1,i1=0..N)]: expand([seq(gu[1][j]-mu[1][j],j=1..nops(mu[1]))]),expand(gu[2]-mu[2]): end: #Tay1(a,b,t,N): The Taylor expansion up to order N #of the GFs(x,a,b)[2] #For example, try: #Tay1(1,1,t,10); Tay1:=proc(a,b,t,N) local S,F,n,j,gu,x,i1: S:=subs({x[a]=1,x[-b]=1},Eqst(x,a,b,F,t)): for n from 1 to N do S:=OneStep(S,t,n) : od: gu:=expand(t*add(subs(S[3],F[0,j]),j=0..b-1)): add(coeff(gu,t,i1)*t^i1,i1=1..N): end: #F1ab(a,b,n): the number of walks with steps (1,a) (1,-b) from #the origin to (n,i) staying weakly above the origin and then #going down #For example, try: #F1ab(1,1,10); F1ab:=proc(a,b,n) local i: add(Fab(a,b,n-1,i),i=0..b-1) : end: #Tay1d(a,b,t,N): The Taylor expansion up to order N #of the GFs(x,a,b)[2], done directly #For example, try: #Tay1d(1,1,t,10); Tay1d:=proc(a,b,t,N) local i: add(F1ab(a,b,i)*t^i,i=1..N): end: