###################################################################### ## BounceArea.txt: Save this file as BounceArea.txt to use it, # # stay in the # ## same directory, get into Maple (by typing: maple ) # #then type: read `AreaBounce.txt`: # ## Then follow the instructions given there # ## # ## Written by AJ Bu and Doron Zeilberger, Rutgers University , # ## DoronZeil@gmail.com . # ###################################################################### with(combinat): Digits:=100: print(`First Written: July 2022: tested for Maple 20 `): print(): print(`This is BounceArea.txt, A Maple package`): print(`accompanying AJ Bu and Doron Zeilberger's article: `): print(`Using Symbolic Computation to Explore Generalized Dyck Paths and Their Areas`): print(`-------------------------------`): print(`-------------------------------`): print(): print(`The most current version is available on WWW at:`): print(` http://www.math.rutgers.edu/~zeilberg/tokhniot/BounceArea.txt .`): print(`Please report all bugs to: DoronZeil at gmail dot com .`): print(): print(`-------------------------------`): print(`For general help, and a list of the MAIN functions,`): print(` type "ezra();". For specific help type "ezra(procedure_name);" `): print(`-------------------------------`): print(`For a list of the supporting functions type: ezra1();`): print(`For help with a specific procedure, type "ezra(procedure_name);"`): print(): print(`-------------------------------`): ezra1:=proc() if args=NULL then print(` The supporting procedures are: `): print(` Area, AreaGF, Bo, BoPath, BoSeq, BoPSeq, Combs, Pl,PlBo, AreaGf, NextPt, PeakSeq, WtCo, WtoP `): else ezra(args): fi: end: ezra:=proc() if args=NULL then print(` BounceArea.txt: A Maple package for experimenting with the Area and Bounce Statistics on Dyck Paths`): print(`The main procedures are : ABgf, Alpha2, DW, Fn, Hn, RandW `): elif nargs=1 and args[1]=ABgf then print(`ABgf(n,q,t): The weight-enumerator of Dyck paths of length n according to the weight q^Area(W)*t^Bo(W) where Bo(W) is the bounce of W. DONE DIRECTLY, rather than using Fn(n,q,t) or Hn(n,q,t).Try:`): print(`ABgf(5,q,t);`): elif nargs=1 and args[1]=Alpha2 then print(`Alpha2(f,x1,x2,K): The list of lists L such that if f(x1,x2) is some enumerating generating function L[i+1][j+1] is the (i,j) scaled moment except`): print(`that L[1][2] and L[2][1] are the averages and L[1][3], L[3][1] are strandard deviations`): elif nargs=1 and args[1]=Area then print(`Area(W): The area of the Dyck path W. Try: `): print(`Area([1,1,1,-1,-1,1,1,-1,-1,-1]);`): elif nargs=1 and args[1]=AreaGF then print(`AreaGF(n,q): The area generating function for Dyck path of semi-length n. Try: `): print(`AreaGF(5,q);`): elif nargs=1 and args[1]=Bo then print(`Bo(W): The Bounce of the Dyck path W. Try: `): print(`Bo([1,1,1,-1,-1,1,1,-1,-1,-1]);`): elif nargs=1 and args[1]=BoPath then print(`BoPath(W), the bounce path of W. Try:`): print(`BoPath([1,1,-1,-1]);`): elif nargs=1 and args[1]=BoSeq then print(`BoSeq(W): The bounce Sequence of the Dyck work W`): elif nargs=1 and args[1]=BoPSeq then print(`BoPSeq(W): The co-bounce Sequence of the Dyck work W`): elif nargs=1 and args[1]=Comps then print(`Comps(n,k): The set of compositions of n into k parts. Try:`): print(`Comps(10,3);`): elif nargs=1 and args[1]=DW then print(`DW(n): The set of Dyck words of semi-length n. Try:`): print(` DW(5); `): elif nargs=1 and args[1]=Fnk then print(`Fnk(n,q,t): Implementing Theorem 13 Jim Haglund's formula Theorem 13 (p. 25) in: "Catalan Paths and q,t-Enumeration. Try:`): print(`Fnk(4,2,q,t):`): elif nargs=1 and args[1]=Fn then print(`Fn(n,q,t): ABgf(n,q,t) via Jim Haglund's formula Theorem 13 (p. 25) in: "Catalan Paths and q,t-Enumeration. Try:`): print(`Fn(4,q,t):`): elif nargs=1 and args[1]=Hnk then print(`Hnk(n,q,t): Implementing Corollary 6,Jim Haglund's formula Theorem 13 (p. 25) in: "Catalan Paths and q,t-Enumeration. the contribution from compsoitions into k parts. Try:`): print(`Hnk(4,2,q,t):`): elif nargs=1 and args[1]=Hn then print(`Hn(n,q,t): ABgf(n,q,t) via Jim Haglund's Corollary 6 (p. 25) in: "Catalan Paths and q,t-Enumeration. Try:`): print(`Hn(4,q,t):`): elif nargs=1 and args[1]=NextPt then print(`NextPt(P,a): Given a path and an integer a finds the b such that if you take a bliiard ball starting at (a,0) and hit the path you would come to (b,0). Try: `): print(`NextPt(WtoP([1,-1,1,-1,1,-1]),0);`): elif nargs=1 and args[1]=PeakSeq then print(`PeakSeq(W): The sequence of projections of the peaks of W on the x-axis divided by 2, except for the last peak. Note that BoPseq(W) is a subsequence. Try:`): print(`PeakSeq([1,1,-1,-1,1,1,-1,-1]);`): elif nargs=1 and args[1]=Pl then print(`Pl(W): plots the path of the Dyck word W. Try:`): print(`Pl([1,1,1,-1,-1,1,1,-1,-1,-1]);`): elif nargs=1 and args[1]=PlBo then print(`PlBo(W,eps): plots the path of the Dyck word W along with the Bounce path with epsilon difference. Try:`): print(`PlBo([1,1,1,-1,-1,1,1,-1,-1,-1],0.1);`): elif nargs=1 and args[1]=RandW then print(`RandW(n): A random Dyck path of semi-length b. Try:`): print(`RandW(50);`): elif nargs=1 and args[1]=WtCo then print(`WtCo(A,q,t): The weight of the composition A , i.e. the summand in Eq. (91) of Haglund's paper. Try:`): print(`WtCo([2,3],q,t);`): elif nargs=1 and args[1]=WtoP then print(`WtoP(W): The path corresponding to the Dyck work W. Try: `): print(`WtoP([1,1,-1,-1]);`): else print(`There is no such thing as`, args): fi: end: Pl:=proc(W): plot(WtoP(W),scaling=constrained):end: DW:=proc(n) local k,gu,gu1,gu2,g1,g2: option remember: if n=0 then RETURN({[]}): fi: gu:={}: for k from 1 to n do gu1:=DW(k-1): gu2:=DW(n-k): gu:=gu union {seq(seq([1,op(g1),-1,op(g2)],g1 in gu1),g2 in gu2)}: od: gu: end: #WtoP(W): The path corresponding to the Dyck work W. Try: #WtoP([1,1,-1,-1]); WtoP:=proc(W) local P,i: P:=[[0,0]]: for i from 1 to nops(W) do P:=[op(P),P[-1]+[1,W[i]]]: od: P: end: #Area(W) Area:=proc(W) local i,P: P:=WtoP(W): (add(P[i][2],i=2..nops(P)-1)-nops(W)/2)/2: end: AreaGF:=proc(n,q) local gu,gu1: gu:=DW(n): add(q^Area(gu1),gu1 in gu): end: #NextPt(P,a): Given a path and an integer a finds the b such that if you take a bliiard ball starting at (a,0) and hit the path you would come to (b,0). Try #NextP(WtoP([1,-1,1,-1,1,-1]), NextPt:=proc(P,a) local i,j: for i from 0 while not member([a+i,i],P) do od: for j from i while member([a+j,j],P) do od: a+2*(j-1): end: #Bo(W): The bounce of the Dyck work W Bo:=proc(W) local n,P,a,co: n:=nops(W)/2: P:=WtoP(W): co:=0: a:=NextPt(P,0) : while a<2*n do co:=co+n-a/2: a:=NextPt(P,a): od: co: end: #ABgf(n,q,t): The weight-enumerator of Dyck paths of length n according to the weight q^Area(W)*t^Bo(W) where Bo(W) is the bounce of W. DONE DIRECTLY, rather than using Gn(n,q,t) or Hn(n,q,t).Try: #ABgf(5,q,t); ABgf:=proc(n,q,t) local gu,gu1: gu:=DW(n): add(q^Area(gu1)*t^Bo(gu1),gu1 in gu): end: #BoSeq(W): The bounce Sequence of the Dyck work W BoSeq:=proc(W) local n,P,a,PI: n:=nops(W)/2: P:=WtoP(W): a:=NextPt(P,0) : PI:=[]: while a<2*n do PI:=[op(PI),n-a/2]: a:=NextPt(P,a): od: PI: end: #BoPSeq(W): The bounce Pre-Sequence of the Dyck work W BoPSeq:=proc(W) local n,P,a,PI: n:=nops(W)/2: P:=WtoP(W): a:=NextPt(P,0) : PI:=[]: while a<2*n do PI:=[op(PI),a/2]: a:=NextPt(P,a): od: PI: end: #PeakSeq(W): The sequence of projections of W on the x-axis. Try: #PeakSeq([1,1,-1,-1,1,1,-1,-1]); PeakSeq:=proc(W) local P,i,gu: P:=WtoP(W): gu:=[]: for i from 2 to nops(P)-1 do if P[i-1][2]Cn(n) then RETURN(FAIL): fi: r:=rand(1..Cn(n))(): for i from 1 to nops(gu) while add(gu[j],j=1..i)[2*n,0] do lu:=NextPts(P,P1[1]): gu:=[op(gu),op(lu)]: P1:=lu[2]: od: gu: end: PlBo:=proc(W,eps) local Peps,i: Peps:=BoPath(W): Peps:=[seq(Peps[i]+[eps,eps],i=1..nops(Peps))]: plot([WtoP(W),Peps],color=[blue,red],scaling=constrained): end: ez:=proc(D): print(`Drs(f,x1,x2,r,s); Alpha2(f,x1,x2,K) `): end: #Drs(f,x1,x2,r,s) (x1d/dx1)^r(x2d/dx2)^s f(x1,x2) at x1=1 x2=1 Drs:=proc(f,x1,x2,r,s) local i,f1: f1:=f: for i from 1 to r do f1:=x1*diff(f1,x1): od: for i from 1 to s do f1:=x2*diff(f1,x2): od: subs({x1=1,x2=1},f1): end: #Alpha2(f,x1,x2,K): The list of lists L such that if f(x1,x2) is some enumerating generating function L[i+1][j+1] is the (i,j) scaled moment except #that L[1][2] and L[2][1] are the averages and L[1][3], L[3][1] are strandard deviations Alpha2:=proc(f,x1,x2,K) local f1,T,i,j,av1,av2,sd1,sd2: f1:=f/subs({x1=1,x2=1},f): av1:=subs({x1=1,x2=1},diff(f1,x1)): av2:=subs({x1=1,x2=1},diff(f1,x2)): T[0,0]:=1: T[1,0]:=av1: T[0,1]:=av2: f1:=f1/x1^av1/x2^av2: sd1:=sqrt(subs({x1=1,x2=1},x1*diff(x1*diff(f1,x1),x1))): sd2:=sqrt(subs({x1=1,x2=1},x2*diff(x2*diff(f1,x2),x2))): T[2,0]:=sd1: T[0,2]:=sd2: for i from 3 to K do T[i,0]:=Drs(f1,x1,x2,i,0)/sd1^i: T[0,i]:=Drs(f1,x1,x2,0,i)/sd2^i: od: for i from 1 to K do for j from 1 to K do T[i,j]:=Drs(f1,x1,x2,i,j)/sd1^i/sd2^j: od: od: evalf([seq([seq(T[i,j],j=0..K)],i=0..K)]): end: