###################################################################### ## ResPermsRM.txt: Save this file as ResPermsRM.txt to use it, # # stay in the # ## same directory, get into Maple (by typing: maple ) # #then type: read `ResPerms.txt`: # ## Then follow the instructions given there # ## # ## Written by George Spahn and Doron Zeilberger, Rutgers University ,# ## DoronZeil@gmail.com . # ###################################################################### #read Pdata: with(combinat): Digits:=100: print(`First Written: Nov. 2022: tested for Maple 20 `): print(): print(`This is ResPermsRM.txt, A second Maple package`): print(`accompanying George Spahn and Doron Zeilberger's article: `): print(` Counting Permutations such that avoid restrictions of the kind "pi[i+s1]-pi[i] is never s2. and |pi[i+s1]-pi[i]| is never s2" `): print(`It contatins everything in the Maple package ResPerms.txt but has new stuff, that can be consulted by typing ezraRM();`): print(`incorporating Rintaro Matsuo's brilliant algorithm and mapping`): print(`-------------------------------`): print(`-------------------------------`): print(): print(`The most current version is available on WWW at:`): print(` http://www.math.rutgers.edu/~zeilberg/tokhniot/ResPerms.txt .`): print(`Please report all bugs to: zeilberg at math dot rutgers dot edu .`): 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(`-------------------------------`): print(`-------------------------------`): print(`For a list of the Dave Robbins style functions type: ezraR();`): print(`For help with a specific procedure, type "ezra(procedure_name);"`): print(): print(`-------------------------------`): print(`-------------------------------`): print(`For a list of the GENERALIZED functions type: ezraG();`): print(`For help with a specific procedure, type "ezra(procedure_name);"`): print(): print(`-------------------------------`): print(`-------------------------------`): print(`For a list of the functions coming from the Multi-variable approach, type: ezraM();`): print(`For help with a specific procedure, type "ezra(procedure_name);"`): print(): print(`-------------------------------`): print(`-------------------------------`): print(`For a list of the functions coming from generalized inclusion-exclusion (keeping track of not only those with no sins), type ezraT();`): print(`For help with a specific procedure, type "ezra(procedure_name);"`): print(): print(`-------------------------------`): print(`-------------------------------`): print(`For a list of the functions relatted to Rintaro Matsuo's approach, type ezraRM`): print(`For help with a specific procedure, type "ezra(procedure_name);"`): print(): print(`-------------------------------`): ezraRM:=proc() if args=NULL then print(` The Rintaro Matsuo procedures are: `): print(` A22rm, A22rmBF, AnrsRM, B22rm, B22rmBF, FindVi1, FindVr, IsRM2, matsuo, RIN, RIN1, RINv, RINv1, RM2, RMg, RMperm, U22rm, U22seqRM, U22seqRMie, U22seqSR, V22seqF, V22seqRM, V22seqRMie `): else ezra(args): fi: end: ezraR:=proc() if args=NULL then print(` The Dave Robbins style procedures are: `): print(` DE1s, f2nR,fsnR, Ker1s, Paper1s, REC1s, V12R, V1sR, WtTilR `): else ezra(args): fi: end: ezraT:=proc() if args=NULL then print(` The procedures keeping track of the number of sins are `): print(` AnRbfT, SinsR, UrsT, VrsT`): else ezra(args): fi: end: ezraG:=proc() if args=NULL then print(` The Generalized procedures are: AnrsBF, AnRbf `): print(` `): else ezra(args): fi: end: ezraM:=proc() if args=NULL then print(` The Multi-variable procedures are: a11R, a22, a22R, Data22, Data22s `): print(` `): else ezra(args): fi: end: ezra1:=proc() if args=NULL then print(` The supporting procedures are: Ank, Bnk, Anrs, AnrsBF, ASrs, BnrsBF, Coe,f1n, f2n, findrec, fsn, GO1s, GG, IsBad1, IsBad, Ope21KK, Ope21SBE, Ope22KK, Ope22KKv, Ope22SBE, SeqFromRec, Til, U11i, U22i, U22nr, U22nri, Ursi, V11f, WtTil `): print(` `): else ezra(args): fi: end: ezra:=proc() if args=NULL then print(` ResPermsRM.txt: A Maple package for counting permutations with the type of restrictions of the form Pi[i+r]-pi[i]<>s and |Pi[i+r]-pi[i]|<>s `): print(`The main procedures are: a11G, An, NuAnrs, U11, U12, U1sEN, U1sF, U22, U22seq, U22seqF, Urs, Urr, V21seq, V21seqF, Vrs, Vrr `): elif nargs=1 and args[1]=a11G then print(`a11G(n): Same as V11f(n) or Vrr(1,n) but using George Spahn's scheme. Try:`): print(`[seq(a11G(n),n=1..30)];`): elif nargs=1 and args[1]=a11R then print(`a11R(n,R): The number of permutations of {1^R,...n^R} such that pi[i+1]-pi[i]<>1 . Try:`): print(`[seq(a11R(i,2),i=1..10)];`): elif nargs=1 and args[1]=a22 then print(`a22(n): The number of permutations of {1,2,...,n} such that pi[i+2]-pi[i]<>2 . So far n must be <=8. It uses the reduced data-base. Try:`): print(`[seq(a22(i,2),i=1..8)];`): elif nargs=1 and args[1]=a22R then print(`a22R(n,R): The number of permutations of {1^R,...n^R} such that pi[i+2]-pi[i]<>2 . So far n must be <=8. It uses the FULL data-base. Try:`): print(`[seq(a22R(i,2),i=1..8)];`): elif nargs=1 and args[1]=A22rm then print(`A22rm(n): The image of Anrs(n,2,2) under the Rintaro Matsuo transform. Try:`): print(`A22rm(7);`): elif nargs=1 and args[1]=A22rmBF then print(`A22rmBF(n): Same as A22rm(n) but by brute force, for checking only. Try:`): print(`evalb(A22rm(6)=A22rmBF(6));`): elif nargs=1 and args[1]=An then print(`An(n):The number of permutations of {1,..,n} such that pi[i+1]-pi[i]<>1 for i=2..n, using Inclusion-Exclusion. Try:`): print(`[seq(An(i),i=1..10)];`): elif nargs=1 and args[1]=Ank then print(`Ank(n,k): The weight-enumeartor of objects with clusters, using the recurrence. Needed for An(n)`): elif nargs=1 and args[1]=Anrs then print(`Anrs(n,r,s): The set of permutations of n such that pi[i+r]-pi[i]<>s for i=1..n-r, using dynamical programming`): elif nargs=1 and args[1]=AnrsBF then print(`AnrsBF(n,r,s): The set of permutations of n such that pi[i+r]-pi[i]<>s for i=1..n-r, using brute force`): elif nargs=1 and args[1]=AnrsRM then print(`AnrsRM(n,r,s): The image of Anrs(n,r,s) under the Rintaro Matsuo transform pi->RMg(pi,r,s). Try:`): print(`AnrsRM(7,2,2); It should be the same as A22rm(7);`): elif nargs=1 and args[1]=AnRbf then print(`AnRbf(n,R): Given sets of pairs {[r,s]}, R, and a positive integer n, outputs the set of permutations such that for all pairs=[r,s] in R `): print(`pi[i+r]-pi[i]<>s. Try:`): print(`AnRbf(7,{[1,1],[2,2]});`): elif nargs=1 and args[1]=AnRbfT then print(`AnRbfT(n,R,t): Given a set of pairs and a positive integer n, outputs the weight-enumerator according to the number of violations`): print(`Try:`): print(`AnRbfT(7,{[1,1],[2,2]},t);`): elif nargs=1 and args[1]=ASrs then print(`ASrs(S,r,s): The set of permutations on the alphabet S such that pi[i+1]-pi[i]<>s for i=1..nops(S)-r, using dynamical programming`): elif nargs=1 and args[1]=ASrsv then print(`ASrsv(S,r,s,v): Inputs a set S of integers, pos. integers r and s and a vector of disctinct positive integers v or length r drawn from the members of S`): print(`outpts the set of permutations of S ending with the vecor v with the property that pi[i+r]-pi[r]<>s. Try:`): print(`ASrsv({1,2,3,4},2,2,[1,2]);`): elif nargs=1 and args[1]=B22rm then print(`B22rm(n): The image of Bnrs(n,2,2) under the Rintaro Matsuo transform. Try:`): print(`B22rm(7);`): elif nargs=1 and args[1]=B22rmBF then print(`B22rmBF(n): Same as B22rm(n) but by brute force, for checking only. UNDER CONSTRUCTION. DOES NOT YET WORK.Try:`): print(`evalb(B22rm(6)=B22rmBF(6));`): elif nargs=1 and args[1]=Bn then print(`Bn(n):The number of permutations of {1,..,n} using the closed-form expression of Ank(n,k), given by Bnk(n,k). Try:`): print(`[seq(Bn(i),i=1..10)];`): elif nargs=1 and args[1]=Bnk then print(`Bnk(n,k): the explicit expression of Ank(n,k), so the same thing, but computed via the formula. Needed for Bn(n)`): elif nargs=1 and args[1]=BnrsBF then print(`BnrsBF(n,r,s): The set of permutations of n such that |pi[i+r]-pi[i]|<>s for i=1..n-r, using brute force`): elif nargs=1 and args[1]=Coe then print(`Coe(f,x,p): The coefficient of the polynomial (in x[1],x[2],x[3],...) , f, of the monomial coming from the partition p`): elif nargs=1 and args[1]=Data22 then print(`Data22(x): The first 8 multi-variable generating functions of avoiding iB(i+2)`): elif nargs=1 and args[1]=Data22s then print(`Data22s(x): The TRIMED first 8 multi-variable generating functions of avoiding iB(i+2)`): elif nargs=1 and args[1]=DE1s then print(`DE1s(s,z,Dz,K): The differential operator in z, Dz, such that applying it to the generating function Sum(Vsr(1,s,i)*z^i,i=0..infinity) is a polynomial of low degree. `): print(`K is a guessing parameter.Try:`): print(`DE1s(2,z,Dz,50);`): elif nargs=1 and args[1]=f1n then print(`f1n(n,x): the weight-enumerator, accoring to the weight Weight(1^r)=x[r] of tilies of the form 1* of the board 1^n. Try:`): print(`f1n(10,x);`): elif nargs=1 and args[1]=f2n then print(`f2n(n,x): the weight-enumerator of tilings with tiles of the form (1X)*1 and 1, where weight(1)=x[1] and weight ((1X)^r1)=x[r+1] of the board 1^n. (X is blank) Try:`): print(`f2n(10,x);`): elif nargs=1 and args[1]=f2nR then print(`f2nR(n,x): The weight enumerator of subsets of {1,...,n-2} with weight(S):= x^|S|*2^NumberOf2-ComponentsOf[S}. Try:`): print(`f2nR(10,x);`): elif nargs=1 and args[1]=fsnR then print(`fsnR(s,n,x): The weight-enumerator of subsets of {1,...,n-2} with tiles {1},{1,1+s},{1,1+s,1+2*s},..., with weight(S):= x^|S|*2^NumberOfs-ComponentsOf[S}. Try:`): print(`fsnR(3,10,x);`): elif nargs=1 and args[1]=f2nr then print(`f2nr(n,r,x): the weight-enumerator of tilings with tiles of the form (1X)*1 and 1, where weight(1)=x[1] and weight( (1X)^r 1)=x[r+1] (X is blank) of the board 1^n(X1)^r. Note that f2nr(n,0,x) is f2n(n,x). Try:`): print(`f2nr(10,2,x);`): elif nargs=1 and args[1]=findrec then print(`findrec(f,DEGREE,ORDER,n,N): guesses a recurrence operator annihilating`): print(`the sequence f of degree DEGREE and order ORDER. The sequences starts at n=0.`): print(`For example, try: findrec([seq(i,i=0..10)],0,2,n,N);`): print(`For example, try: findrec([1,op(U22seqF(68)]),3,13,n,N);`): elif nargs=1 and args[1]=FindVi1 then print(`FindVi1(pi): Given a permutation pi finds the list of pairs [i,j] such that pi[i]=j and pi[i+1]=j+1. Try:`): print(`FindVi1([1,2,3,4]);`): elif nargs=1 and args[1]=FindVr then print(`FindVr(n,r): All the violations [i,j] in the members of AnrsRM(n,r,r) minus Anrs(n,1,1). Try:`): print(`FindVr(7,2);`): elif nargs=1 and args[1]=fsn then print(`fsn(s,n,x): the weight-enumerator of tilings with tiles of the form (1X^(s-1))*1 and 1, where weight(1)=x[1] and weight ((1X^(s-1))^r1)=x[s+1] of the board 1^n. (X is blank) Try:`): print(`fsn(2,10,x);`): elif nargs=1 and args[1]=IsBad then print(`IsBad(pi,S): does the permutation pi contain the schenario, for some i, pi[i+r]-pi[i]=s for some i?. and some [r,s] in S`): print(`Try IsBad([1,2,3,4,5],{[2,2]});`): elif nargs=1 and args[1]=IsBad1 then print(`IsBad1(pi,r,s): does the permutation pi contain the schenario, for some i, pi[i+r]-pi[i]=s for some i?.`): print(`Try IsBad1([1,2,3,4,5],2,2);`): elif nargs=1 and args[1]=IsRM2 then print(`IsRM2(pi): Is the permutation good in the sense of Rintaro Matsuo?. Try:`): print(`IsRM2([1,2,3,4]);`): elif nargs=1 and args[1]=Ker1s then print(`Ker1s(s,x,z,K): The rational function in x and z, let's call it F(x,z) such that Vrs(1,s,n) equals the coeff. of z^n in int(exp(-x)*F(x,z),x=0..infinity). `): print(`K is a guessing paramater. Try: `): print(`Ker1s(1,x,z,50);`): elif nargs=1 and args[1]=GO1s then print(`GO1s(s,n,N): Guesses an operator annihilationg Urs(1,s,n); Try:`): print(`GO1s(3,n,N);`): elif nargs=1 and args[1]=GG then print(`GG(n,S): given a positive integer n and set of patterns {[r,s]} outputs the set of permutations such that for each member of S pi[i+r]-pi[i] is never s`): elif nargs=1 and args[1]=matsuo then print(`matsuo(n,i,j,k,l,m): The quantity needed in Matsuo Rintaro's approach for computing U22(n), used in the more efficient procedure U22rm(n). Try:`): print(`matsuo(5,1,1,1,1,1);`): elif nargs=1 and args[1]=NuAnrs then print(`NuAnrs(n,r,s): The Number of permutations of n such that pi[i+r]-pi[i]<>s for i=1..n-r, using dynamical programming`): elif nargs=1 and args[1]=Ope21KK then print(`Ope21KK(n,N): the Kauers-Koutscan CONJETURED operator of order 8 and degree 7, annihilating {Vrs(1,2,n)}. It outputs the pair [ope,Initial conditions]. try:`): print(`Ope21KK(n,N);`): elif nargs=1 and args[1]=Ope21SBE then print(`Ope21SBE(n,N): the SBE rendition of the Kauers-Koutscan operator, Ope21KK, of order 10 and degree 3, for the (2,1)-case. It outputs the pair [ope,Initial conditions]. try:`): print(`Ope21SBE(n,N);`): elif nargs=1 and args[1]=Ope22KK then print(`Ope22KK(n,N): the CONJECTURED (but certain!) Kauers-Koutscan operator of order 8 and degree 11, for the sequence Urr(2,n) (what's called a_{22}(n) in the paper. It outputs the pair [ope,Initial conditions]. try:`): print(`ope:=Ope22KK(n,N): SeqFromRec(ope[1],n,N,ope[2],100);`): elif nargs=1 and args[1]=Ope22KKv then print(`Ope22KKv(n,N): the CONJECTURED (but certain!) Kauers-Koutscan operator of order 24 and degree 64, for the sequence Vrr(2,n), what's called b_{22}(n) in the paper. It outputs the pair [ope,Initial conditions]. try:`): print(`ope:=Ope22KKv(n,N):SeqFromRec(ope[1],n,N,ope[2],100);`): elif nargs=1 and args[1]=Ope22SBE then print(`Ope22SBE(n,N): the SBE rendition of the Kauers-Koutscan operator, of order 13 and degree 3, for the (2,2)-case. It outputs the pair [ope,Initial conditions]. try:`): print(`Ope22SBE(n,N);`): elif nargs=1 and args[1]=Paper1s then print(`Paper1s(s,n,K,M1,M2): An article about the sequence enumerating permutations of {1,..,n} such that |pi[i+1]-pi[i]|<>s. Try:`): print(`Paper1s(1,n,30,100,1000):`): elif nargs=1 and args[1]=REC1s then print(`REC1s(s,n,N,K): The recurrence operator in n, N, annihilating Vsr(1,s,n), followed by the initial conditions.`): print(`K is a guessing parameter.Try:`): print(` REC1s(2,n,N,50): `): elif nargs=1 and args[1]=RIN then print(`RIN(n,a,b): The number of permutations pi of {1,...,n} such that pi[i+1]-pi[i]<>1 except possibly when pi[i]=b and i=a using EXCLUSION EXCLUSION. Try:`): print(`RIN(10,5,5).`): elif nargs=1 and args[1]=RINv then print(`RINv(n,a,b): The number of permutations pi of {1,...,n} such that |pi[i+1]-pi[i]|<>1 except possibly when pi[i]=b and i=a using EXCLUSION EXCLUSION. Try:`): print(`RINv(10,5,5).`): elif nargs=1 and args[1]=RIN1 then print(`RIN1(n,a,b,k1,k2): Number of configurations of increasing components with a barrier between i=a and i=a+1 and whree you never have in any component, b followed immediately by b+1`): print(`and there are k1 components to the left of i=a and k2 components to the right of i=a. Needed to compute RIN(n,a,b) in the Inclusion-Exclusion argument. Try:`): print(`RIN1(10,3,3,2,2);`): elif nargs=1 and args[1]=RINv1 then print(`RINv1(n,a,b,k1,k2): Number of configurations of increasing components with a barrier between i=a and i=a+1 and whree you never have in any component, b next to b+1 (on either sides)`): print(`and there are k1 components to the left of i=a and k2 components to the right of i=a. Needed to compute RINv(n,a,b) in the Inclusion-Exclusion argument. Try:`): print(`RINv1(10,3,3,2,2);`): elif nargs=1 and args[1]=RM2 then print(`RM2(pi): The Rintaro Matsuo transform that takes members of Anrs(n,2,2) to The set of permutations avoiding [i,i+1] for i=1...n-1 except that [floor((n+1)/2),floor((n+1)/2)+1] is allowed and`): print(`it is ok to have a violation at the pair of locations [floor((n+1)/2),floor((n+1)/2)+1]. Try:`): print(`RM2([3,5,6,2,8,7,1,4]);`): elif nargs=1 and args[1]=RMg then print(`RMg(pi,g1,g2): The Generalized Rintaro Matsuo transform that takes members of Anrs(n,g1,g2) to The some set of permutations tbd`): print(`RMg(pi,2,2) is the same as RM2(pi). Try:`): print(`pi:=randperm(12); RMg(pi,3,3);`): elif nargs=1 and args[1]=RMperm then print(`RMperm(n,g): The Rintaro Matsuo's g-permutation of length n. Try:`): print(`RMperm(11,2);`): elif nargs=1 and args[1]=SeqFromRec then print(`SeqFromRec(ope,n,N,Ini,K): Given the first L`): print(`terms of the sequence Ini=[f(1), ..., f(L)]`): print(`satisfied by the recurrence ope(n,N)f(n)=0`): print(`extends it to the first K values. try`): print(`SeqFromRec(N-n-1,n,N,[1],10);`): elif nargs=1 and args[1]=SinsR then print(`SinsR(pi,R): Given a permutation pi and a set of restrictions R finds how many places have pi[i+r]-pi[i] in R for [r,s] in R. Try:`): print(`SinsR([1,2,3,4],{[1,1]});`): elif nargs=1 and args[1]=Til then print(`Til(R,T): Given a set of integers and a set of tiles, outputs all the tilings. Try:`): print(`Til({1,2,3,4,5,6},{{1,3},{1,3,5},{1}});`): elif nargs=1 and args[1]=U11 then print(`U11(n):The number of permutations of {1,..,n} such that pi[i+1]-pi[i]<>1 `): elif nargs=1 and args[1]=U11i then print(`U11i(n): Like U11(n) but using integration. Try: U11i(5); `): elif nargs=1 and args[1]=U12 then print(`U12(n):The number of permutations of {1,..,n} such that pi[i+1]-pi[i]<>2 `): elif nargs=1 and args[1]=U1sEN then print(`U1sEN(s,n): Same as Urs(1,s,n) but much faster, using Enrique Navarrete's formula`): print(`U1sEN(5,20);`): elif nargs=1 and args[1]=U1sF then print(`U1sF(s,n): Same as Urs(1,s,n) but much faster, using the recurrence. Try:`): print(`U1sF(5,20);`): elif nargs=1 and args[1]=U22 then print(`U22(n):The number of permutations of {1,..,n} such that pi[i+2]-pi[i]<>2. Try: `): print(`U22(10);`): elif nargs=1 and args[1]=U22rm then print(`U22rm(n): Same as U22(n) (q.v.), i.e.the number of permutations of {1,..,n} such that pi[i+2]-pi[i]<>2, BUT using Rintaro Matsuo's algorithm implemented in Maple. Try: `): print(`U22rm(10);`): elif nargs=1 and args[1]=U22i then print(`U22i(n): Same as U22(n), i.e. the number of permutations of {1,..,n} such that pi[i+2]-pi[i]<>2, BUT USING INTEGRATION OF KERNEL, hence slower, only of theoretical interest `): elif nargs=1 and args[1]=U22nr then print(`U22nr(n,r): Applying the umbra that gives U22(n) out of f2n(n,x) to the more general f2nr(n,r,x). Try:`): print(`U22nr(5,5);`): elif nargs=1 and args[1]=U22nri then print(`U22nri(n,r): Same as U22nr(n), BUT USING INTEGRATION `): print(`U22nri(5,5);`): elif nargs=1 and args[1]=V21seqF then print(`V21seqF(K): the first K terms of the sequence U22(n). USING THE CONJECTURE RECURRENCE.Try:`): print(`V21seqF(30);`): elif nargs=1 and args[1]=U22seq then print(`U22seq(K): the first K terms of the sequence U22(n). Try:`): print(`U22seq(30);`): elif nargs=1 and args[1]=U22seqRM then print(`U22seqRM(K): the first K terms of the sequence U22rm(n), i.e. same as U22seq(K), but using Rintaro Matsuo's algorithm. Try:`): print(`U22seqRM(30);`): elif nargs=1 and args[1]=U22seqRMie then print(`U22seqRMie(K): Same as U22seqRM(K), and hence the same as U22seq(K), but using the Inclusion-Exclusion approach, after applying Rintaro Matsuo's bijection. Try:`): print(`U22seqRMie(30);`): elif nargs=1 and args[1]=U22seqF then print(`U22seqF(K): the first K terms of the sequence U22(n). USING THE CONJECTURE RECURRENCE.Try:`): print(`U22seqF(30);`): elif nargs=1 and args[1]=Urs then print(`Urs(r,s,n): The number of permutations of {1,..,n} such that pi[i+r]-pi[i]<>s `): elif nargs=1 and args[1]=UrsT then print(`UrsT(r,s,n,t): The weight-enumerator of permutations of {1,..,n} according to the weight t^(#Number of occurences of pi[i+r]-pi[i]=s). Try: `): print(`UrsT(2.3,7,t);`): elif nargs=1 and args[1]=Urr then print(`Urr(r,n): The number of permutations of {1,..,n} such that pi[i+r]-pi[i]<>r. Same as Urs(r,r,n) but faster `): elif nargs=1 and args[1]=Ursi then print(`Ursi(r,s,n): Same as Urs(r,s,n), i.e., the number of permutations of {1,..,n} such that pi[i+r]-pi[i]<>s, BUT USING INTEGRATION, hence slower, only of theoretical interest. `): elif nargs=1 and args[1]=V11f then print(`V11f(n): Same output as Vrr(1,n), but using the formula in Roberto Tauraso's article INTEGERS 6 (2006) #A11. Try:`): print(`[seq(V11f(i),i=1..20)];`): elif nargs=1 and args[1]=V12R then print(`V12R(n): Same as Vrs(1,2,n) but using Dave Robbins' method. In other words the number of permutations of {1,...,n} such that |pi[i+1]-pi[i]|<>2. Try:`): print(`[seq(V12R(i),i=1..10)];`): elif nargs=1 and args[1]=V1sR then print(`V1sR(s,n): Same as Vrs(1,s,n) but using Dave Robbins' method. In other words the number of permutations of {1,...,n} such that |pi[i+1]-pi[i]|<>s. Try:`): print(`[seq(V1sR(3,i),i=1..10)];`): elif nargs=1 and args[1]=V21seq then print(`V21seq(K): the first K terms of the sequence Vrs(2,1,n).Try:`): print(`V21seq(20);`): elif nargs=1 and args[1]=V21seqF then print(`V21seqF(K): the first K terms of the sequence Vrs(2,1,n). USING THE CONJECTURE RECURRENCE.Try:`): print(`V21seqF(30);`): elif nargs=1 and args[1]=V22seqRMie then print(`V22seqRMie(K): Same as [seq(Vrr(2,i),i=1..N)],but using the Inclusion-Exclusion approach, after applying Rintaro Matsuo's bijection. Try:`): print(`V22seqRMie(30);`): elif nargs=1 and args[1]=Vrs then print(`Vrs(r,s,n): The number of permutations of {1,..,n} such that |pi[i+r]-pi[i]|<>s `): elif nargs=1 and args[1]=VrsT then print(`VrsT(r,s,n,t): The weight-enumerator of permutations of {1,..,n} according to the weight t^(#Number of occurences of |pi[i+r]-pi[i]|=s). Try: `): print(`VrsT(2.3,7,t);`): elif nargs=1 and args[1]=Vrr then print(`Vrr(r,n): The number of permutations of {1,..,n} such that |pi[i+r]-pi[i]|<>r. Same as Vrs(r,r,n) but faster `): elif nargs=1 and args[1]=WtTil then print(`WtTil(R,T,x): Given a set of integers and a set of tiles, outputs the weight-enumerator of all tilings. Try:`): print(`WtTil({1,2,3,4,5,6},{{1,3},{1,3,5},{1}},x);`): elif nargs=1 and args[1]=WtTilR then print(`WtTilR(R,T,x): Given a set of integers and a set of tiles T, outputs the setight-enumerator of all the set of all tilings with the weight`): print(`of a singleton being 1 and the weight of a larger tile t, being x^(nops(t)-1)*2. Needed for fsnR(s,n,x), and V1sR(s,n). Try:`): print(`WtTilR({1,2,3,4,5,6},{{1,3},{1,3,5},{1}},x);`): else print(`There is no such thing as`, args): fi: end: with(combinat): ################################Start From EKHAD.txt with(SolveTools): solve1:=Linear: pashetc:=proc(p,D) local gu,p1,i,mekh: p1:=normal(p): mekh:=denom(p1): p1:=numer(p1): p1:=expand(p1): gu:=0: for i from 0 to degree(p1,D) do gu:=gu+factor(coeff(p1,D,i))*D^i: od: RETURN(gu,mekh): end: goremc:=proc(D,ORDER) local i,gu: gu:=0: for i from 0 to ORDER do gu:=gu+(bpc[i])*D^i: od: RETURN(gu): end: gzor:=proc(f,x,n) local i,gu: gu:=f: for i from 1 to n do gu:=diff(gu,x): od: RETURN(gu): end: duisc:= proc(A,ORDER,y,x,D) local KAMA,gorem,conj, yakhas,lu1a,LU1,P,Q,R,S,j1,g,Q1,Q2,l1,l2,mumu, mekb1,fu,meka1,k1,gugu,eq,ia1,va1,va2,degg,i,shad: KAMA:=30: gorem:=goremc(D,ORDER): conj:=gorem: yakhas:=0: for i from 0 to degree(conj,D) do yakhas:=yakhas+normal(coeff(conj,D,i)*gzor(A,x,i)/A): yakhas:=normal(yakhas): od: lu1a:=yakhas: LU1:=numer(yakhas): yakhas:=1/denom(yakhas): yakhas:=normal(diff(yakhas,y)/yakhas+diff(A,y)/A): P:=LU1: Q:=numer(yakhas): R:=denom(yakhas): j1:=0: while j1 <= KAMA do g:=gcd(R,Q-j1*diff(R,y)): if g <> 1 then Q2:=(Q-j1*diff(R,y))/g: R:=normal(R/g): Q:=normal(Q2+j1*diff(R,y)): P:=P*g^j1: j1:=-1: fi: j1:=j1+1: od: P:=expand(P): R:=expand(R): Q:=expand(Q): Q1:=Q+diff(R,y): Q1:=expand(Q1): l1:=degree(R,y): l2:=degree(Q1,y): meka1:=coeff(R,y,l1): mekb1:=coeff(Q1,y,l2): if l1-1 <>l2 then k1:=degree(P,y)-max(l1-1,l2): else mumu:= -mekb1/meka1: if type(mumu,integer) and mumu > 0 then k1:=max(mumu, degree(P,y)-l2): else k1:=degree(P,y)-l2: fi: fi: fu:=0: if k1 < 0 then RETURN(0): fi: for ia1 from 0 to k1 do fu:=fu+apc[ia1]*y^ia1: od: gugu:=Q1*fu+R*diff(fu,y)-P: gugu:=expand(gugu): degg:=degree(gugu,y): for ia1 from 0 to degg do eq[ia1+1]:=coeff(gugu,y,ia1)=0: od: for ia1 from 0 to k1 do va1[ia1+1]:=apc[ia1]: od: for ia1 from 0 to ORDER do va2[ia1+1]:=bpc[ia1]: od: eq:=convert(eq,set): va1:=convert(va1,set): va2:=convert(va2,set): va1:=va1 union va2 : va1:=solve1(eq,va1): fu:=subs(va1,fu): gorem:=subs(va1,gorem): if fu=0 and gorem=0 then RETURN(0): fi: for ia1 from 0 to k1 do gorem:=subs(apc[ia1]=1,gorem): od: for ia1 from 0 to ORDER do gorem:=subs(bpc[ia1]=1,gorem): od: fu:=normal(fu): shad:=pashetc(gorem,D): S:=lu1a*R*fu*shad[2]/P: S:=subs(va1,S): S:=normal(S): S:=factor(S): for ia1 from 0 to k1 do S:=subs(apc[ia1]=1,S): od: for ia1 from 0 to ORDER do S:=subs(bpc[ia1]=1,S): od: shad[1],S: end: AZcI:=proc(A,y,x,D,a,b) local ORDER,gu,KAMA,ope,cert,yemin1: KAMA:=8: for ORDER from 1 to KAMA do gu:=duisc(A,ORDER,y,x,D): if gu<>0 and gu[1]<>0 then ope:=gu[1]: cert:=gu[2]: yemin1:=normal(subs(y=b, normal(cert*A))-subs(y=a, normal(cert*A))): RETURN([[ope,yemin1],cert]): fi: od: 0: end: ################################End EKHAD.txt #################Start from Cfinite.txt #GuessRec1(L,d): inputs a sequence L and tries to guess #a recurrence operator with constant cofficients of order d #satisfying it. It returns the initial d values and the operator # as a list. For example try: #GuessRec1([1,1,1,1,1,1],1); GuessRec1:=proc(L,d) local eq,var,a,i,n: if nops(L)<=2*d+2 then print(`The list must be of size >=`, 2*d+3 ): RETURN(FAIL): fi: var:={seq(a[i],i=1..d)}: eq:={seq(L[n]-add(a[i]*L[n-i],i=1..d),n=d+1..nops(L))}: var:=solve(eq,var): if var=NULL then RETURN(FAIL): else RETURN([[op(1..d,L)],[seq(subs(var,a[i]),i=1..d)]]): fi: end: #GuessRec(L): inputs a sequence L and tries to guess #a recurrence operator with constant cofficients #satisfying it. It returns the initial values and the operator # as a list. For example try: #GuessRec([1,1,1,1,1,1]); GuessRec:=proc(L) local gu,d: for d from 1 to trunc(nops(L)/2)-2 do gu:=GuessRec1(L,d): if gu<>FAIL then RETURN(gu): fi: od: FAIL: end: #SeqFromRecC(S,N): Inputs S=[INI,ope] #where INI is the list of initial conditions, a ope a list of #size L, say, and a recurrence operator ope, codes a list of #size L, finds the first N0 terms of the sequence satisfying #the recurrence f(n)=ope[1]*f(n-1)+...ope[L]*f(n-L). #For example, for the first 20 Fibonacci numbers, #try: SeqFromRecC([[1,1],[1,1]],20); SeqFromRecC:=proc(S,N) local gu,L,n,i,INI,ope: INI:=S[1]:ope:=S[2]: if not type(INI,list) or not type(ope,list) then print(`The first two arguments must be lists `): RETURN(FAIL): fi: L:=nops(INI): if nops(ope)<>L then print(`The first two arguments must be lists of the same size`): RETURN(FAIL): fi: if not type(N,integer) then print(`The third argument must be an integer`, L): RETURN(FAIL): fi: if Nexpand(SeqFromRecC(S,L+1)) then print([seq(coeff(f1,t,i),i=0..L)],SeqFromRecC(S,L+1)): RETURN(FAIL): else RETURN(f): fi: end: #################End from Cfinite.txt ###From FindRec.txt Yafe:=proc(ope,N) local i,ope1,coe1,L: if ope=0 then RETURN(1,0): fi: ope1:=expand(ope): L:=degree(ope1,N): coe1:=coeff(ope1,N,L): ope1:=normal(ope1/coe1): ope1:=normal(ope1): ope1:= convert( [seq(factor(coeff(ope1,N,i))*N^i,i=ldegree(ope1,N)..degree(ope1,N))],`+`): factor(coe1),ope1: end: #SeqFromRec(ope,n,N,Ini,K): Given the first L-1 #terms of the sequence Ini=[f(1), ..., f(L-1)] #satisfied by the recurrence ope(n,N)f(n)=0 #extends it to the first K values SeqFromRec:=proc(ope,n,N,Ini,K) local ope1,gu,L,n1,j1: ope1:=Yafe(ope,N)[2]: L:=degree(ope1,N): if nops(Ini)<>L then ERROR(`Ini should be of length`, L): fi: if K<=nops(Ini) then RETURN([op(1..K,Ini)]): fi: ope1:=expand(subs(n=n-L,ope1)/N^L): gu:=Ini: for n1 from nops(Ini)+1 to K do gu:=[op(gu), -add(gu[nops(gu)+1-j1]*subs(n=n1,coeff(ope1,N,-j1)), j1=1..L)]: od: gu: end: #findrec(f,DEGREE,ORDER,n,N): guesses a recurrence operator annihilating #the sequence f of degree DEGREE and order ORDER. The sequences starts at n=0. #For example, try: findrec([seq(i,i=0..10)],0,2,n,N); findrec:=proc(f,DEGREE,ORDER,n,N) local ope,var,eq,i,j,n0,kv,var1,eq1,mu,a: option remember: if (1+DEGREE)*(1+ORDER)+ORDER>nops(f) then ERROR(`Insufficient data for a recurrence of order`,ORDER, `degree`,DEGREE): fi: ope:=0: var:={}: for i from 0 to ORDER do for j from 0 to DEGREE do ope:=ope+a[i,j]*N^i*n^j: var:=var union {a[i,j]}: od: od: var:=var minus {a[ORDER,0]}: ope:=subs(a[ORDER,0]=1,ope): eq:={}: for n0 from 0 to nops(f)-ORDER-1 do eq1:=0: for i from 0 to ORDER do eq1:=eq1+subs(n=n0,coeff(ope,N,i))*op(n0+i+1,f): od: eq:= eq union {eq1}: od: var1:=solve(eq,var): if var1=NULL then RETURN(FAIL): fi: ope:=subs(var1,ope): add(factor(coeff(ope,N,i))*N^i,i=0..ORDER): end: #Findrec1(f,n,N,DEGREE): guesses a recurrence operator annihilating #the sequence f with as degree DEGREE. Try: #Findrec1([seq(i,i=0..10)],n,N,0); Findrec1:=proc(f,n,N,DEGREE) local gu,ORDER,ope: for ORDER from 0 while (1+DEGREE)*(1+ORDER)+ORDER<=nops(f)-3 do ope:=findrec(f,DEGREE,ORDER,n,N): if ope<>FAIL then RETURN(Yafe(ope,N)[2]): fi: od: FAIL: end: ###End From FindRec.txt #Ope21KK(n,N): the Kauers-Koutscan operator for the (2,1)-case. It outputs the pair [ope,Initial conditions] Ope21KK:=proc(n,N) local i: [ -(3+n)*(1+n)*(2*n^5+47*n^4+500*n^3+2927*n^2+9224*n+11992)+(-2*n^7-59*n^6-806*n^ 5-6440*n^4-31672*n^3-94045*n^2-154400*n-107168)*N+(8*n^7+222*n^6+2803*n^5+20415 *n^4+90129*n^3+231363*n^2+302284*n+140888)*N^2+(8*n^6+224*n^5+2670*n^4+17018*n^ 3+63518*n^2+134550*n+121692)*N^3+(-6*n^7-209*n^6-3164*n^5-27636*n^4-149934*n^3-\ 499311*n^2-927548*n-727632)*N^4+(2*n^7+71*n^6+1092*n^5+9694*n^4+54270*n^3+ 190635*n^2+377204*n+310424)*N^5+(3+n)*(10*n^5+237*n^4+2480*n^3+14549*n^2+46400* n+67088)*N^6+(-4*n^6-118*n^5-1484*n^4-10760*n^3-46772*n^2-113610*n-114900)*N^7+ (4*n^5+74*n^4+664*n^3+3378*n^2+9384*n+10480)*N^8, [1, 2, 2, 8, 28, 152, 952, 7208]]: end: #Ope21SBE(n,N): the Shalosh operator for the (2,1)-case. It outputs the pair [ope,Initial conditions] Ope21SBE:=proc(n,N) local i: [ (3+n)*(980453*n+3990186)*(1+n)+(1522666*n^3+16519562*n^2+57292242*n+62716844)*N+(-4089147*n^3-38912424*n^2-123637984*n-129513299)*N^2 +(-2878400*n^3-46484816*n^2-232547312*n-377085582)*N^3+(5779551*n^3+99906908*n^2+497904198*n+743421267)*N^4+(646186*n^3+18951230*n^2+ 154547618*n+377930244)*N^5+(-2670857*n^3-73541872*n^2-590915284*n-1425749491)*N^6+(709548*n^3+18941288*n^2+165440372*n+443190910)*N^7 +(4632166*n^2+72558479*n+270723095)*N^8+(-1419096*n^2-29114552*n-128057688)*N^9+(1419096*n+9581878)*N^10, [1, 2, 2, 8, 28, 152, 952, 7208, 62296, 605864]] end: #Ope22SBE(n,N): the Shalosh operator for the (2,2)-case. It outputs the pair [ope,Initial conditions] Ope22SBE:=proc(n,N) local i: [1/4*n*(n-1)^2+1/16*n*(11*n^2+38*n-16)*N+(139/16*n+29/8*n^2+3/16*n^3+3/2)*N^2+(-51/4*n-79/16*n^2-5/8*n^3-21/4)*N^3+(-1/8*n+5/4*n^2+1/8*n^3-15/2)*N^4+(307/4*n+49/4*n^2+11/16*n^3+603/4)*N^5+(-533/16*n-49/8*n^2-5/16*n^3-41)*N^6 +(-303/16*n^2-3/4*n^3-161*n-911/2)*N^7+(-1/4*n^3-417/4*n-37/4*n^2-363)*N^8+(-53*n-11/4*n^2-993/4)*N^9+(-93/4*n-n^2-130)*N^10+(-2*n-71/4)*N^11+(-n-10)*N^12+N^13, [1, 2, 5, 18, 75, 410, 2729, 20906, 181499, 1763490, 18943701, 222822578, 2847624899]]: end: #Ope22KK(n,N): the Kauers-Koutscan operator for the (2,2)-case. It outputs the pair [ope,Initial conditions] Ope22KK:=proc(n,N) local i: [ (262711*n + 1387742*n^2 - 824875*n^3 - 1855253*n^4 - 111530*n^5 + 680983*n^6 + 364242*n^7 + 84992*n^8 + 10332*n^9 + 640*n^10 + 16*n^11) + (-1050844*n - 9705192*n^2 - 7414683*n^3 + 3536494*n^4 + 6459004*n^5 + 3326393*n^6 + 903534*n^7 + 144684*n^8 + 13756*n^9 + 720*n^10 + 16*n^11)*N + (3492344 - 2212342*n - 8507169*n^2 - 11544227*n^3 - 12034116*n^4 - 8216995*n^5 - 3442049*n^6 - 890050*n^7 - 142300*n^8 - 13660*n^9 - 720*n^10 - 16*n^11)*N^2 + (19817984 + 45323852*n + 825228*n^2 - 57004661*n^3 - 57059306*n^4 - 28077270*n^5 - 8398637*n^6 - 1631510*n^7 - 207980*n^8 - 16828*n^9 - 784*n^10 - 16*n^11)*N^3 + (9586160 + 6680237*n - 13772613*n^2 - 27689586*n^3 - 22162455*n^4 - 9855085*n^5 - 2629562*n^6 - 427656*n^7 - 41332*n^8 - 2176*n^9 - 48*n^10)*N^4 + (22192864 + 44710768*n - 2924668*n^2 - 52385912*n^3 - 45161616*n^4 - 18784740*n^5 - 4549208*n^6 - 674256*n^7 - 60400*n^8 - 3008*n^9 - 64*n^10)*N^5 + (557152 - 2032472*n - 2937392*n^2 - 1594200*n^3 - 517688*n^4 - 122032*n^5 - 19856*n^6 - 1792*n^7 - 64*n^8)*N^6 + (3786960 + 7105324*n - 1191064*n^2 - 8059160*n^3 - 5938996*n^4 - 2073752*n^5 - 402736*n^6 - 44528*n^7 - 2624*n^8 - 64*n^9)*N^7 + (-598208 - 943004*n + 414196*n^2 + 1213772*n^3 + 728648*n^4 + 203584*n^5 + 29616*n^6 + 2176*n^7 + 64*n^8)*N^8, [1, 2, 5, 18, 75, 410, 2729, 20906]]: end: Data22s:=proc(x) [-1/(-1+x[1]), -1/(-1+x[1]+x[2]), -1/(-x[1]*x[2]*x[3]+x[1]+x[2]+x[3]-1), -1/(2*x[1]*x[2]*x[3]*x[4]-x[1]*x[2]*x[3]-x[1]*x[2]*x[4]-x[1]*x[3]*x[4]-x[2]*x[3]*x[4]+x[1]+x[2]+x[3]+x[4]-1), (-x[1]*x[3]*x[5]+1)/(-x[1]*x[2]*x[3]*x[4]*x[5]-2*x[1]*x[2]*x[3]*x[4]+x[1]*x[2]*x[3]*x[5]+x [1]*x[3]*x[4]*x[5]-2*x[2]*x[3]*x[4]*x[5]+x[1]*x[2]*x[3]+x[1]*x[2]*x[4]+x[1]*x[3]*x[4]+x[1]*x[3]*x[5]+x[2]*x[3]*x[4]+x[2]*x[3]*x[5]+x[2]*x[4]*x[5]+x[3]*x[4]*x[5]-x[1]-x[2]-x[3]-x[4]-x[5]+1), (-x[1]*x[3]*x[5]-x[2]*x[4]*x[6]+1)/(2*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]-x[1]*x[2]*x[3]* x[4]*x[5]-x[1]*x[2]*x[3]*x[4]*x[6]-2*x[1]*x[2]*x[3]*x[5]*x[6]-2*x[1]*x[2]*x[4]*x[5]*x[6]-x[1]*x[3]*x[4]*x[5]*x[6]-x[2]*x[3]*x[4]*x[5]*x[6]-2*x[1]*x[2]*x[3]*x[4]+x[1]*x[2]*x[3]*x[5]+x[1]*x[2]*x[4]*x[6]+x[1]*x[3]*x[4]*x[5]-2*x[1]*x[3]*x[4]*x[6]+x[1]*x[3]*x[5]*x[6]-2*x[2]*x[3 ]*x[4]*x[5]+x[2]*x[3]*x[4]*x[6]+x[2]*x[4]*x[5]*x[6]-2*x[3]*x[4]*x[5]*x[6]+x[1]*x[2]*x[3]+x[1]*x[2]*x[4]+x[1]*x[3]*x[4]+x[1]*x[3]*x[5]+x[1]*x[3]*x[6]+x[1]*x[4]*x[6]+x[2]*x[3]*x[4]+x[2]*x[3]*x[5]+x[2]*x[4]*x[5]+x[2]*x[4]*x[6]+x[3]*x[4]*x[5]+x[3]*x[4]*x[6]+x[3]*x[5]*x[6]+x[4] *x[5]*x[6]-x[1]-x[2]-x[3]-x[4]-x[5]-x[6]+1), (x[1]*x[2]*x[3]*x[5]*x[7]+x[1]*x[3]*x[4]*x[5]*x[7]+x[1]*x[3]*x[5]*x[6]*x[7]-x[1]*x[3]*x[5]-x[2]*x[4]*x[6]-x[3]*x[5]*x[7]+1)/(3*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]*x[7]+2*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]+x[1]*x[2]*x[3]*x[4]*x[6]*x[7]-2*x[ 1]*x[2]*x[3]*x[5]*x[6]*x[7]+x[1]*x[2]*x[4]*x[5]*x[6]*x[7]+2*x[2]*x[3]*x[4]*x[5]*x[6]*x[7]-x[1]*x[2]*x[3]*x[4]*x[5]-x[1]*x[2]*x[3]*x[4]*x[6]-2*x[1]*x[2]*x[3]*x[5]*x[6]-3*x[1]*x[2]*x[3]*x[5]*x[7]-2*x[1]*x[2]*x[4]*x[5]*x[6]-2*x[1]*x[2]*x[4]*x[6]*x[7]-x[1]*x[3]*x[4]*x[5]*x[6]-\ 3*x[1]*x[3]*x[4]*x[5]*x[7]-3*x[1]*x[3]*x[5]*x[6]*x[7]-x[2]*x[3]*x[4]*x[5]*x[6]-x[2]*x[3]*x[4]*x[5]*x[7]-2*x[2]*x[3]*x[4]*x[6]*x[7]-2*x[2]*x[3]*x[5]*x[6]*x[7]-x[2]*x[4]*x[5]*x[6]*x[7]-x[3]*x[4]*x[5]*x[6]*x[7]-2*x[1]*x[2]*x[3]*x[4]+x[1]*x[2]*x[3]*x[5]+x[1]*x[2]*x[4]*x[6]+x[1 ]*x[3]*x[4]*x[5]-2*x[1]*x[3]*x[4]*x[6]+x[1]*x[3]*x[5]*x[6]-2*x[2]*x[3]*x[4]*x[5]+x[2]*x[3]*x[4]*x[6]+x[2]*x[3]*x[5]*x[7]+x[2]*x[4]*x[5]*x[6]-2*x[2]*x[4]*x[5]*x[7]+x[2]*x[4]*x[6]*x[7]-2*x[3]*x[4]*x[5]*x[6]+x[3]*x[4]*x[5]*x[7]+x[3]*x[5]*x[6]*x[7]-2*x[4]*x[5]*x[6]*x[7]+x[1]*x [2]*x[3]+x[1]*x[2]*x[4]+x[1]*x[3]*x[4]+x[1]*x[3]*x[5]+x[1]*x[3]*x[6]+x[1]*x[3]*x[7]+x[1]*x[4]*x[6]+x[1]*x[5]*x[7]+x[2]*x[3]*x[4]+x[2]*x[3]*x[5]+x[2]*x[4]*x[5]+x[2]*x[4]*x[6]+x[2]*x[4]*x[7]+x[2]*x[5]*x[7]+x[3]*x[4]*x[5]+x[3]*x[4]*x[6]+x[3]*x[5]*x[6]+x[3]*x[5]*x[7]+x[4]*x[5] *x[6]+x[4]*x[5]*x[7]+x[4]*x[6]*x[7]+x[5]*x[6]*x[7]-x[1]-x[2]-x[3]-x[4]-x[5]-x[6]-x[7]+1), (-2*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]*x[7]*x[8]+x[1]*x[2]*x[3]*x[5]*x[7]+x[1]*x[2]*x[4]*x[6]*x[8]+x[1]*x[3]*x[4]*x[5]*x[7]+x[1]*x[3]*x[5]*x[6]*x[7]+x[1]*x[3]*x[5]*x[7]*x[8]+x[2]*x[3]*x[4] *x[6]*x[8]+x[2]*x[4]*x[5]*x[6]*x[8]+x[2]*x[4]*x[6]*x[7]*x[8]-x[1]*x[3]*x[5]-x[2]*x[4]*x[6]-x[3]*x[5]*x[7]-x[4]*x[6]*x[8]+1)/(3*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]*x[7]+3*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]*x[8]+4*x[1]*x[2]*x[3]*x[4]*x[5]*x[7]*x[8]+4*x[1]*x[2]*x[3]*x[4]*x[6]*x[7]*x[8]+ 4*x[1]*x[2]*x[3]*x[5]*x[6]*x[7]*x[8]+4*x[1]*x[2]*x[4]*x[5]*x[6]*x[7]*x[8]+3*x[1]*x[3]*x[4]*x[5]*x[6]*x[7]*x[8]+3*x[2]*x[3]*x[4]*x[5]*x[6]*x[7]*x[8]+2*x[1]*x[2]*x[3]*x[4]*x[5]*x[6]+x[1]*x[2]*x[3]*x[4]*x[5]*x[8]+x[1]*x[2]*x[3]*x[4]*x[6]*x[7]-2*x[1]*x[2]*x[3]*x[5]*x[6]*x[7]+x [1]*x[2]*x[3]*x[5]*x[6]*x[8]-2*x[1]*x[2]*x[3]*x[5]*x[7]*x[8]+x[1]*x[2]*x[4]*x[5]*x[6]*x[7]-2*x[1]*x[2]*x[4]*x[5]*x[6]*x[8]-2*x[1]*x[2]*x[4]*x[6]*x[7]*x[8]+2*x[1]*x[3]*x[4]*x[5]*x[6]*x[8]-2*x[1]*x[3]*x[4]*x[5]*x[7]*x[8]+x[1]*x[3]*x[4]*x[6]*x[7]*x[8]+x[1]*x[4]*x[5]*x[6]*x[7] *x[8]+2*x[2]*x[3]*x[4]*x[5]*x[6]*x[7]+x[2]*x[3]*x[4]*x[5]*x[7]*x[8]-2*x[2]*x[3]*x[4]*x[6]*x[7]*x[8]+x[2]*x[3]*x[5]*x[6]*x[7]*x[8]+2*x[3]*x[4]*x[5]*x[6]*x[7]*x[8]-x[1]*x[2]*x[3]*x[4]*x[5]-x[1]*x[2]*x[3]*x[4]*x[6]-2*x[1]*x[2]*x[3]*x[5]*x[6]-3*x[1]*x[2]*x[3]*x[5]*x[7]-2*x[1]* x[2]*x[3]*x[5]*x[8]-2*x[1]*x[2]*x[4]*x[5]*x[6]-2*x[1]*x[2]*x[4]*x[6]*x[7]-3*x[1]*x[2]*x[4]*x[6]*x[8]-x[1]*x[3]*x[4]*x[5]*x[6]-3*x[1]*x[3]*x[4]*x[5]*x[7]-2*x[1]*x[3]*x[4]*x[5]*x[8]-x[1]*x[3]*x[4]*x[6]*x[8]-3*x[1]*x[3]*x[5]*x[6]*x[7]-x[1]*x[3]*x[5]*x[6]*x[8]-3*x[1]*x[3]*x[5] *x[7]*x[8]-2*x[1]*x[4]*x[5]*x[6]*x[8]-2*x[1]*x[4]*x[6]*x[7]*x[8]-x[2]*x[3]*x[4]*x[5]*x[6]-x[2]*x[3]*x[4]*x[5]*x[7]-2*x[2]*x[3]*x[4]*x[6]*x[7]-3*x[2]*x[3]*x[4]*x[6]*x[8]-2*x[2]*x[3]*x[5]*x[6]*x[7]-2*x[2]*x[3]*x[5]*x[7]*x[8]-x[2]*x[4]*x[5]*x[6]*x[7]-3*x[2]*x[4]*x[5]*x[6]*x[8 ]-3*x[2]*x[4]*x[6]*x[7]*x[8]-x[3]*x[4]*x[5]*x[6]*x[7]-x[3]*x[4]*x[5]*x[6]*x[8]-2*x[3]*x[4]*x[5]*x[7]*x[8]-2*x[3]*x[4]*x[6]*x[7]*x[8]-x[3]*x[5]*x[6]*x[7]*x[8]-x[4]*x[5]*x[6]*x[7]*x[8]-2*x[1]*x[2]*x[3]*x[4]+x[1]*x[2]*x[3]*x[5]+x[1]*x[2]*x[4]*x[6]+x[1]*x[3]*x[4]*x[5]-2*x[1]*x [3]*x[4]*x[6]+x[1]*x[3]*x[5]*x[6]+x[1]*x[3]*x[5]*x[8]-2*x[1]*x[3]*x[6]*x[8]+x[1]*x[4]*x[6]*x[8]-2*x[2]*x[3]*x[4]*x[5]+x[2]*x[3]*x[4]*x[6]+x[2]*x[3]*x[5]*x[7]+x[2]*x[4]*x[5]*x[6]-2*x[2]*x[4]*x[5]*x[7]+x[2]*x[4]*x[6]*x[7]-2*x[3]*x[4]*x[5]*x[6]+x[3]*x[4]*x[5]*x[7]+x[3]*x[4]*x [6]*x[8]+x[3]*x[5]*x[6]*x[7]-2*x[3]*x[5]*x[6]*x[8]+x[3]*x[5]*x[7]*x[8]-2*x[4]*x[5]*x[6]*x[7]+x[4]*x[5]*x[6]*x[8]+x[4]*x[6]*x[7]*x[8]-2*x[5]*x[6]*x[7]*x[8]+x[1]*x[2]*x[3]+x[1]*x[2]*x[4]+x[1]*x[3]*x[4]+x[1]*x[3]*x[5]+x[1]*x[3]*x[6]+x[1]*x[3]*x[7]+x[1]*x[3]*x[8]+x[1]*x[4]*x[6 ]+x[1]*x[5]*x[7]+x[1]*x[6]*x[8]+x[2]*x[3]*x[4]+x[2]*x[3]*x[5]+x[2]*x[4]*x[5]+x[2]*x[4]*x[6]+x[2]*x[4]*x[7]+x[2]*x[4]*x[8]+x[2]*x[5]*x[7]+x[2]*x[6]*x[8]+x[3]*x[4]*x[5]+x[3]*x[4]*x[6]+x[3]*x[5]*x[6]+x[3]*x[5]*x[7]+x[3]*x[5]*x[8]+x[3]*x[6]*x[8]+x[4]*x[5]*x[6]+x[4]*x[5]*x[7]+x [4]*x[6]*x[7]+x[4]*x[6]*x[8]+x[5]*x[6]*x[7]+x[5]*x[6]*x[8]+x[5]*x[7]*x[8]+x[6]*x[7]*x[8]-x[1]-x[2]-x[3]-x[4]-x[5]-x[6]-x[7]-x[8]+1)]: end: #IsBad1(pi,r,s): does the permutation pi contain the schenario, for some i, pi[i+r]-pi[i]=s for some i?. #Try IsBad1([1,2,3,4,5],2,2); IsBad1:=proc(pi,r,s) local n,i: n:=nops(pi): if r=0 then RETURN(FAIL): fi: if r>0 then for i from 1 to n-r do if pi[i+r]-pi[i]=s then RETURN(true): fi: od: RETURN(false): fi: if r<0 then for i from -r+1 to n do if pi[i+r]-pi[i]=s then RETURN(true): fi: od: RETURN(false): fi: end: #IsBad(pi,S): does the permutation pi contain the schenario, for some i, pi[i+r]-pi[i]=s for some i?. and some [r,s] in S #Try IsBad([1,2,3,4,5],{[2,2]}); IsBad:=proc(pi,S) local t: for t in S do if IsBad1(pi,t[1],t[2]) then RETURN(true): fi: od: false: end: #GG(n,S): given a positive integer n and set of patterns {[r,s]} outputs the set of permutations such that for each member of S pi[i+r]-pi[i] is never s GG:=proc(n,S) local gu,mu,pi: mu:=permute(n): gu:={}: for pi in mu do if not IsBad(pi,S) then gu:=gu union {pi}: fi: od: gu: end: #GGg(T,S): given a set T and set of patterns {[r,s]} outputs the set of permutations of T such that for each member of S pi[i+r]-pi[i] is never s GGg:=proc(T,S) local gu,mu,pi: mu:=convert(permute(T),set): gu:={}: for pi in mu do if not IsBad(pi,S) then gu:=gu union {pi}: fi: od: gu: end: #Anrs(n,r,s): The set of permutations of n such that pi[i+r]-pi[i]<>s for i=1..n-r, using dynamical programming Anrs:=proc(n,r,s) local i1 : ASrs({seq(i1,i1=1..n)},r,s): end: #ASrs(S,r,s): The set of permutations on the alphabet S such that pi[i+r]-pi[i]<>s for i=1..nops(S)-r, using dynamical programming ASrs:=proc(S,r,s) local i1,Vecs,v: if nops(S)s. Try: #ASrsv({1,2,3,4},2,2,[1,2]); ASrsv:=proc(S,r,s,v) local S1,Y1,u,i1,Vecs,G,lu,ru,ru1,pi: option remember: if S={op(v)} then RETURN({v}): fi: if nops(S)<2*r then ru:=GGg(S,{[r,s]}): ru1:={}: for pi in ru do if [op(nops(pi)-r+1..nops(pi),pi)]=v then ru1:=ru1 union {pi}: fi: od: RETURN(ru1): fi: S1:=S minus {op(v)}: Vecs:=choose(convert(S1,list),r): Vecs:={seq(op(permute(Vecs[i1])),i1=1..nops(Vecs))}: Y1:={}: for u in Vecs do if not member(s,{seq(v[i1]-u[i1],i1=1..r)}) then Y1:=Y1 union {u}: fi: od: G:={}: for u in Y1 do lu:=ASrsv(S1,r,s,u): G:=G union {seq([op(lu[i1]),op(v)],i1=1..nops(lu))}: od: G: end: #NuASrsv(S,r,s,v): Inputs a set S of integers, pos. integers r and s and a vector of disctinct positive integers v or length r drawn from the members of S #outpts the NUMBER of permutations of S ending with the vecor v with the property that pi[i+r]-pi[r]<>s. Try: #NuASrsv({1,2,3,4},2,2,[1,2]); NuASrsv:=proc(S,r,s,v) local S1,Y1,u,i1,Vecs,G,ru,ru1,pi: option remember: if S={op(v)} then RETURN(1): fi: if nops(S)<2*r then ru:=GGg(S,{[r,s]}): ru1:={}: for pi in ru do if [op(nops(pi)-r+1..nops(pi),pi)]=v then ru1:=ru1 union {pi}: fi: od: RETURN(nops(ru1)): fi: S1:=S minus {op(v)}: Vecs:=choose(convert(S1,list),r): Vecs:={seq(op(permute(Vecs[i1])),i1=1..nops(Vecs))}: Y1:={}: for u in Vecs do if not member(s,{seq(v[i1]-u[i1],i1=1..r)}) then Y1:=Y1 union {u}: fi: od: G:=0: for u in Y1 do G:=G+NuASrsv(S1,r,s,u): od: G: end: #NuAnrs(n,r,s): The Number of permutations of n such that pi[i+r]-pi[i]<>s for i=1..n-r, using dynamical programming NuAnrs:=proc(n,r,s) local i1 : NuASrs({seq(i1,i1=1..n)},r,s): end: #NuASrs(S,r,s): The Number of permutations on the alphabet S such that pi[i+r]-pi[i]<>s for i=1..nops(S)-r, using dynamical programming NuASrs:=proc(S,r,s) local i1,Vecs,v: if nops(S)2 #f2n(n,x): the weight-enumerator of tilings of the form (1X)* f2n:=proc(n,x) option remember: f2nr(n,0,x): end: #f2nr(n,r,x): The weight enumerator of tilings of the board 1^n(X1)^r f2nr:=proc(n,r,x) local i: option remember: #if n<0 or r<0 then # RETURN(0): #fi: if n=0 and r=0 then RETURN(1): fi: expand(add(x[i]*f2nr(n,r-i,x),i=1..r)+ add(x[r+i]*f2nr(n+1-2*i,i-1,x),i=1..trunc((n+1)/2))): end: #U22oldold(n): U22oldold:=proc(n) local x,f,i,P,c,j,s,p1: f:=f2n(n,x): P:=partition(n): s:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c:=f: for j from 1 to nops(p1) do c:=coeff(c,x[p1[j][1]],p1[j][2]): od: s:=s+c^2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: s: end: #U12(n): U12:=proc(n) local x,f1,f2,i,P,c1,c2,j,s,p1: f1:=f1n(n,x): f2:=f2n(n,x): P:=partition(n): s:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: s:=s+c1*c2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: s: end: #Coe(f,x,p): The coefficient in f of the partition coming from p1 Coe:=proc(f,x,p) local c,j,p1: p1:=ParToF(p): c:=f: for j from 1 to nops(p1) do c:=coeff(c,x[p1[j][1]],p1[j][2]): od: c: end: #Ank(n,k): the weight-enumeartor of objects with clusters, using the recurrence Ank:=proc(n,k) option remember: if k=1 then RETURN((-1)^(n+1)): elif k>n or k<0 then RETURN(0): else RETURN(k*Ank(n-1,k-1)-Ank(n-1,k)): fi: end: #An(n):The number of permutations of {1,..,n} such that pi[i+1]-pi[i]<>1 for i=2..n, using Inclusion-Exclusion An:=proc(n) local k: add(Ank(n,k),k=1..n): end: Bnk:=proc(n,k) : (-1)^(n+k)*k*(n-1)!/(n-k)!:end: Bn:=proc(n) option remember: if n=1 then 1: else ((n+1)*(n-1)*Bn(n-1)+(-1)^(n-1))/n: fi: end: #start tiling #Til(R,T): Given a set of integers and a set of tiles, outputs all the tilings. Try: #Til({1,2,3,4,5,6},{{1,3},{1,3,5},{1})}; Til:=proc(R,T) local i,S,t,R1,s,t1,S1: option remember: if R={} then RETURN({{}}): fi: i:=min(op(R)): S:={}: for t in T do if {seq(i-1+t1,t1 in t)} minus R={} then R1:=R minus {seq(i-1+t1,t1 in t)}: S1:=Til(R1,T): S:=S union {seq(s union {[i,t]}, s in S1)}: fi: od: S: end: #Wt(t,x): The weight of the tiling t Wt:=proc(t,x) local t1: mul(x[nops(t1[2])],t1 in t): end: #WtTil(R,T,X): Given a set of integers and a set of tiles, outputs the setight-enumerator of all the set of all tilings. Try: #WtTil({1,2,3,4,5,6},{{1,3},{1,3,5},{1}},X); WtTil:=proc(R,T,X) local i,S,t,R1,s,t1: option remember: if R={} then RETURN(1): fi: i:=min(op(R)): S:=0: for t in T do if {seq(i-1+t1,t1 in t)} minus R={} then R1:=R minus {seq(i-1+t1,t1 in t)}: S:=expand(S+X[nops(t)]*WtTil(R1,T,X)): fi: od: S: end: #fsn(s,n,x): fsn(n,x) done directly fsn:=proc(s,n,x) local R,T,i,gu,r: if s=1 then RETURN(f1n(n,x)): elif s=2 then RETURN(f2n(n,x)): else R:={seq(i,i=1..n)}: T:={}: for r from 0 to trunc((n-1)/s) do T:=T union {{seq(1+s*i,i=0..r)}}: od: WtTil(R,T,x): fi: end: #end tiling #Urs(r,s,n): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s Urs:=proc(r,s,n) local x,f1,f2,i,P,j,su,p1,c1,c2: f1:=fsn(r,n,x): f2:=fsn(s,n,x): P:=partition(n): su:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: su:=su+c1*c2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: su: end: a11Old:=proc(n) local x,t,gu,i: if n>15 then RETURN(FAIL): fi: gu:=Data11(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: a22Old:=proc(n) local x,t,gu,i: if n>8 then RETURN(FAIL): fi: gu:=Data22(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: #e11n(n,x): e11n:=proc(n,x) local i,j: option remember: if n=1 then RETURN(1-x[1]): fi: e11n(n-1,x)+add((-1)^(i-1)*mul(x[n-j],j=0..i),i=0..n-1): end: #a11(n): Like U11(n) but via the multi-generating function (much slower) a11:=proc(n) local x,t,gu,i: gu:=1/e11n(n,x): gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: a11R:=proc(n,R) local x,t,gu,i: gu:=1/e11n(n,x): gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,R*n+1): gu:=expand(coeff(gu,t,R*n)): for i from 1 to n do gu:=coeff(gu,x[i],R): od: gu: end: a22R:=proc(n,R) local x,t,gu,i: if n>8 then RETURN(FAIL): fi: gu:=Data22(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,R*n+1): gu:=expand(coeff(gu,t,R*n)): for i from 1 to n do gu:=coeff(gu,x[i],R): od: gu: end: #Trim1(P,x,n): Trims the polynomial in x[1],...,x[n] by removing larger powerss than 1 Trim1:=proc(P,x,n) local P1,i: P1:=P: for i from 1 to n do P1:=expand(coeff(P1,x[i],0)+coeff(P1,x[i],1)*x[i]): od: P1: end: Trim:=proc(P,x,n) Trim1(numer(P),x,n)/Trim1(denom(P),x,n): end: a22:=proc(n) local x,t,gu,i: if n>8 then RETURN(FAIL): fi: gu:=Data22s(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: #U22i(n): Like U22(n) but using integration U22i:=proc(n) local x,gu,z,f,i: f:=f2n(n,x): gu:=(-1)^n*subs({seq(x[i]=-x[i],i=1..n)},f)*subs({seq(x[i]=z[i]/x[i],i=1..n)},f): for i from 1 to n do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n)): for i from 1 to n do gu:=int(gu,z[i]=0..infinity): od: gu: end: #Ursi(r,s,n): Like Urs(r,s,n) but using integration Ursi:=proc(r,s,n) local x,gu,z,f1,f2,i: f1:=fsn(r,n,x): f2:=fsn(s,n,x): gu:=(-1)^n*subs({seq(x[i]=-x[i],i=1..n)},f1)*subs({seq(x[i]=z[i]/x[i],i=1..n)},f2): for i from 1 to n do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n)): for i from 1 to n do gu:=int(gu,z[i]=0..infinity): od: gu: end: #U22nr(n,r): Applying the umbra that gives U22(n) out of f2n(n,x) to the more general f2nr(n,r,x). Try: #U22nr(5,5); U22nr:=proc(n,r) local x,f,i,P,c,j,s,p1: f:=f2nr(n,r,x): P:=partition(n+r): s:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c:=f: for j from 1 to nops(p1) do c:=coeff(c,x[p1[j][1]],p1[j][2]): od: s:=s+c^2*(-1)^(n+r)*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: s: end: #U22nri(n,r): U22nri:=proc(n,r) local x,gu,z,f,i: f:=f2nr(n,r,x): gu:=(-1)^(n+r)*subs({seq(x[i]=-x[i],i=1..n+r)},f)*subs({seq(x[i]=z[i]/x[i],i=1..n+r)},f): for i from 1 to n+r do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n+r)): for i from 1 to n+r do gu:=int(gu,z[i]=0..infinity): od: gu: end: #U22old(n): U22old:=proc(n) local x,f,i,i1,f1,j,s,c,p1: if n=1 then RETURN(1): elif n=2 then RETURN(2): else f:=f2n(n,x): s:=0: for i from 1 to nops(f) do f1:=op(i,f): p1:=[seq(degree(f1,x[i1]),i1=1..n)]: c:=normal(f1/mul(x[i1]^p1[i1],i1=1..n)): s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*s): fi: end: #U22(n): U22:=proc(n) local x,f,i,i1,f1,j,s,c,p1: if n=1 then RETURN(1): elif n=2 then RETURN(2): else f:=f2n(n,x): s:=0: for i from 1 to nops(f) do f1:=op(i,f): c:=op(1,f1): if not type(c,integer) then c:=1: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*s): fi: end: #U22seq(N): the first K terms of the sequence U22(n). Try: #U22seq(30); U22seq:=proc(K) local i: [seq(U22(i),i=1..K)]: end: #U22seqRM(N): the first K terms of the sequence U22rm(n) (i.e. using Rintaro Matsuo's algorithm). Try: #U22seqRM(30); U22seqRM:=proc(K) local i: [seq(U22rm(i),i=1..K)]: end: #U22seqF(K): the first K terms of the sequence U22(n), by using the CONJECTURED recurrence. Try: #U22seqF(200); U22seqF:=proc(K) local n,N,ope: ope:=Ope22SBE(n,N): SeqFromRec(ope[1],n,N,ope[2],K): end: #V21seq(N): the first K terms of the sequence V21(n). Try: #V21seq(30); V21seq:=proc(K) local i: [seq(Vrs(2,1,i),i=1..K)]: end: #V21seqF(K): the first K terms of the sequence Vrs(2,1,n), by using the CONJECTURED recurrence. Try: #V21seqF(200); V21seqF:=proc(K) local n,N,ope: ope:=Ope21SBE(n,N): SeqFromRec(ope[1],n,N,ope[2],K): end: IsGoodrs:=proc(pi,r,s) local i: for i from 1 to nops(pi)-r do if pi[i+r]-pi[i]=s then RETURN(false): fi: od: true: end: #AnrsBF(n,r,s): Like Anrs(n,r,s) but by brute-force. try: #AnrsBF(6,1,1); AnrsBF:=proc(n,r,s) local mu,mu1,gu: mu:=permute(n): gu:={}: for mu1 in mu do if IsGoodrs(mu1,r,s) then gu:=gu union {mu1}: fi: od: gu: end: IsGoodrsA:=proc(pi,r,s) local i: for i from 1 to nops(pi)-r do if abs(pi[i+r]-pi[i])=s then RETURN(false): fi: od: true: end: #BnrsBF(n,r,s): #BnrsBF(6,1,1); BnrsBF:=proc(n,r,s) local mu,mu1,gu: mu:=permute(n): gu:={}: for mu1 in mu do if IsGoodrsA(mu1,r,s) then gu:=gu union {mu1}: fi: od: gu: end: #IsGoodR(pi,R): is the permutation pi such that for each pair [r,s] in R and s in S pi[i+r]-pi[i]<>s. Try: #IsGoodR([1,2,3,4],{[2,2],[1,1]}); IsGoodR:=proc(pi,R) local pa: for pa in R do if not IsGoodrs(pi,pa[1],pa[2]) then RETURN(false): fi: od: true: end: #GO1s(s,n,N): Guesses an operator annihilationg Urs(1,s,n); Try: #GO1s(3,n,N); GO1s:=proc(s, n, N) local n1: subs(n = n - s, findrec([seq(Urs(1, s, n1), n1 = s .. s+10)], 1, 2, n, N)); end: #U1sF(s,n): Same as Urs(1,s,n) but much faster, using the recurrence. Try: #U1sF(5,20); U1sF:=proc(s,n) option remember: if n<0 then RETURN(0): fi: if (n<=s and n>=0) then RETURN(n!): fi: (n-1)*U1sF(s,n-1)+(n-1-s)*U1sF(s,n-2): end: #U1sEN(s,n): Same as Urs(1,s,n) but much faster, using Enrique Navarrete's formula #U1sEN(5,20); U1sEN:=proc(s,n) local j: add((-1)^j*binomial(n-s,j)*(n-j)!,j=0..n-s): end: #U11i(n): Like U11(n) but using integration. Try: U11(7); U11i:=proc(n) local x,gu,z,f,i: f:=f1n(n,x): gu:=(-1)^(n)*subs({seq(x[i]=-x[i],i=1..n)},f)*subs({seq(x[i]=z[i]/x[i],i=1..n)},f): for i from 1 to n do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n)): for i from 1 to n do gu:=int(gu,z[i]=0..infinity): od: gu: end: #AnRbf(n,R): Given a set of pairs and a positive integer n, outputs the set of permutations such that for all the pairs [r,s] in R #all s in R pi[i+r]-pi[i]<>s. Try: #AnRbf(7,{[1,1],[2,2]}); AnRbf:=proc(n,R) local mu,mu1,gu,pair: mu:=permute(n): gu:={}: for mu1 in mu do if IsGoodR(mu1,R) then gu:=gu union {mu1}: fi: od: gu: end: ###New stuff #Vrs(r,s,n): The number of permutations of {1,...,n} such that |pi[i+r]-pi[r]|<>s Vrs:=proc(r,s,n) local x,f1,f2,i,P,j,su,p1,c1,c2,i1: f1:=fsn(r,n,x): f1:=subs({seq(x[i1]=2*x[i1],i1=2..n)},f1): f2:=fsn(s,n,x): P:=partition(n): su:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: su:=su+c1*c2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: su: end: #Vrr(r,n): Same as Vrs(r,s,n) with s=r, but hopefully faster/ Vrr:=proc(r,n) local x,f,i,i1,f1,j,s,c,p1: f:=fsn(r,n,x): if not type(f,`+`) then RETURN(Vrs(r,r,n)): fi: s:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c:=1: else c:=op(1,f1): if not type(c,integer) then c:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1))*2^add(p1[i1],i1=2..n): od: RETURN((-1)^n*s): end: #V11f(n): Same output as Vrr(1,n), but using the formula in Roberto Tauraso's article INTEGERS 6 (2006) #A11. Try: #[seq(V11f(i),i=1..20)]; V11f:=proc(n) local r,c: add((-1)^r*(n-r)!*add(2^c*binomial(r-1,c-1)*binomial(n-r,c),c=0..r),r=0..n-1):end: ##Start Dynamical programming for set of restictions #ASrsv(S,r,s,v): Inputs a set S of integers, pos. integers r and s and a vector of disctinct positive integers v or length r drawn from the members of S #outpts the set of permutations of S ending with the vecor v with the property that pi[i+r]-pi[r]<>s. Try: #ASrsv({1,2,3,4},2,2,[1,2]); ASrsv:=proc(S,r,s,v) local S1,Y1,u,i1,Vecs,G,lu,ru,ru1,pi: option remember: if S={op(v)} then RETURN({v}): fi: if nops(S)<2*r then ru:=GGg(S,{[r,s]}): ru1:={}: for pi in ru do if [op(nops(pi)-r+1..nops(pi),pi)]=v then ru1:=ru1 union {pi}: fi: od: RETURN(ru1): fi: S1:=S minus {op(v)}: Vecs:=choose(convert(S1,list),r): Vecs:={seq(op(permute(Vecs[i1])),i1=1..nops(Vecs))}: Y1:={}: for u in Vecs do if not member(s,{seq(v[i1]-u[i1],i1=1..r)}) then Y1:=Y1 union {u}: fi: od: G:={}: for u in Y1 do lu:=ASrsv(S1,r,s,u): G:=G union {seq([op(lu[i1]),op(v)],i1=1..nops(lu))}: od: G: end: #NuASrsv(S,r,s,v): Inputs a set S of integers, pos. integers r and s and a vector of disctinct positive integers v or length r drawn from the members of S #outpts the NUMBER of permutations of S ending with the vecor v with the property that pi[i+r]-pi[r]<>s. Try: #NuASrsv({1,2,3,4},2,2,[1,2]); NuASrsv:=proc(S,r,s,v) local S1,Y1,u,i1,Vecs,G,ru,ru1,pi: option remember: if S={op(v)} then RETURN(1): fi: if nops(S)<2*r then ru:=GGg(S,{[r,s]}): ru1:={}: for pi in ru do if [op(nops(pi)-r+1..nops(pi),pi)]=v then ru1:=ru1 union {pi}: fi: od: RETURN(nops(ru1)): fi: S1:=S minus {op(v)}: Vecs:=choose(convert(S1,list),r): Vecs:={seq(op(permute(Vecs[i1])),i1=1..nops(Vecs))}: Y1:={}: for u in Vecs do if not member(s,{seq(v[i1]-u[i1],i1=1..r)}) then Y1:=Y1 union {u}: fi: od: G:=0: for u in Y1 do G:=G+NuASrsv(S1,r,s,u): od: G: end: #NuAnrs(n,r,s): The Number of permutations of n such that pi[i+r]-pi[i]<>s for i=1..n-r, using dynamical programming NuAnrs:=proc(n,r,s) local i1 : NuASrs({seq(i1,i1=1..n)},r,s): end: #NuASrs(S,r,s): The Number of permutations on the alphabet S such that pi[i+r]-pi[i]<>s for i=1..nops(S)-r, using dynamical programming NuASrs:=proc(S,r,s) local i1,Vecs,v: if nops(S)2 #f2n(n,x): the weight-enumerator of tilings of the form (1X)* f2n:=proc(n,x) option remember: f2nr(n,0,x): end: #f2nr(n,r,x): The weight enumerator of tilings of the board 1^n(X1)^r f2nr:=proc(n,r,x) local i: option remember: if n=0 and r=0 then RETURN(1): fi: expand(add(x[i]*f2nr(n,r-i,x),i=1..r)+ add(x[r+i]*f2nr(n+1-2*i,i-1,x),i=1..trunc((n+1)/2))): end: #U22oldold(n): U22oldold:=proc(n) local x,f,i,P,c,j,s,p1: f:=f2n(n,x): P:=partition(n): s:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c:=f: for j from 1 to nops(p1) do c:=coeff(c,x[p1[j][1]],p1[j][2]): od: s:=s+c^2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: s: end: #U12(n): U12:=proc(n) local x,f1,f2,i,P,c1,c2,j,s,p1: f1:=f1n(n,x): f2:=f2n(n,x): P:=partition(n): s:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: s:=s+c1*c2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: s: end: #Coe(f,x,p): The coefficient in f of the partition coming from p1 Coe:=proc(f,x,p) local c,j,p1: p1:=ParToF(p): c:=f: for j from 1 to nops(p1) do c:=coeff(c,x[p1[j][1]],p1[j][2]): od: c: end: #Ank(n,k): the weight-enumeartor of objects with clusters, using the recurrence Ank:=proc(n,k) option remember: if k=1 then RETURN((-1)^(n+1)): elif k>n or k<0 then RETURN(0): else RETURN(k*Ank(n-1,k-1)-Ank(n-1,k)): fi: end: #An(n):The number of permutations of {1,..,n} such that pi[i+1]-pi[i]<>1 for i=2..n, using Inclusion-Exclusion An:=proc(n) local k: add(Ank(n,k),k=1..n): end: Bnk:=proc(n,k) : (-1)^(n+k)*k*(n-1)!/(n-k)!:end: Bn:=proc(n) option remember: if n=1 then 1: else ((n+1)*(n-1)*Bn(n-1)+(-1)^(n-1))/n: fi: end: #start tiling #Til(R,T): Given a set of integers and a set of tiles, outputs all the tilings. Try: #Til({1,2,3,4,5,6},{{1,3},{1,3,5},{1})}; Til:=proc(R,T) local i,S,t,R1,s,t1,S1: option remember: if R={} then RETURN({{}}): fi: i:=min(op(R)): S:={}: for t in T do if {seq(i-1+t1,t1 in t)} minus R={} then R1:=R minus {seq(i-1+t1,t1 in t)}: S1:=Til(R1,T): S:=S union {seq(s union {[i,t]}, s in S1)}: fi: od: S: end: #Wt(t,x): The weight of the tiling t Wt:=proc(t,x) local t1: mul(x[nops(t1[2])],t1 in t): end: #WtTil(R,T,X): Given a set of integers and a set of tiles, outputs the setight-enumerator of all the set of all tilings. Try: #WtTil({1,2,3,4,5,6},{{1,3},{1,3,5},{1}},X); WtTil:=proc(R,T,X) local i,S,t,R1,s,t1: option remember: if R={} then RETURN(1): fi: i:=min(op(R)): S:=0: for t in T do if {seq(i-1+t1,t1 in t)} minus R={} then R1:=R minus {seq(i-1+t1,t1 in t)}: S:=expand(S+X[nops(t)]*WtTil(R1,T,X)): fi: od: S: end: #fsn(s,n,x): fsn(n,x) done directly fsn:=proc(s,n,x) local R,T,i,gu,r: if s=1 then RETURN(f1n(n,x)): elif s=2 then RETURN(f2n(n,x)): else R:={seq(i,i=1..n)}: T:={}: for r from 0 to trunc((n-1)/s) do T:=T union {{seq(1+s*i,i=0..r)}}: od: WtTil(R,T,x): fi: end: #end tiling #UrsSlow(r,s,n): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s UrsSlow:=proc(r,s,n) local x,f1,f2,i,P,j,su,p1,c1,c2: f1:=fsn(r,n,x): f2:=fsn(s,n,x): P:=partition(n): su:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: su:=su+c1*c2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: su: end: a11Old:=proc(n) local x,t,gu,i: if n>15 then RETURN(FAIL): fi: gu:=Data11(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: a22Old:=proc(n) local x,t,gu,i: if n>8 then RETURN(FAIL): fi: gu:=Data22(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: #e11n(n,x): e11n:=proc(n,x) local i,j: option remember: if n=1 then RETURN(1-x[1]): fi: e11n(n-1,x)+add((-1)^(i-1)*mul(x[n-j],j=0..i),i=0..n-1): end: #a11(n): Like U11(n) but via the multi-generating function (much slower) a11:=proc(n) local x,t,gu,i: gu:=1/e11n(n,x): gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: a11R:=proc(n,R) local x,t,gu,i: gu:=1/e11n(n,x): gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,R*n+1): gu:=expand(coeff(gu,t,R*n)): for i from 1 to n do gu:=coeff(gu,x[i],R): od: gu: end: a22R:=proc(n,R) local x,t,gu,i: if n>8 then RETURN(FAIL): fi: gu:=Data22(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,R*n+1): gu:=expand(coeff(gu,t,R*n)): for i from 1 to n do gu:=coeff(gu,x[i],R): od: gu: end: #Trim1(P,x,n): Trims the polynomial in x[1],...,x[n] by removing larger powerss than 1 Trim1:=proc(P,x,n) local P1,i: P1:=P: for i from 1 to n do P1:=expand(coeff(P1,x[i],0)+coeff(P1,x[i],1)*x[i]): od: P1: end: Trim:=proc(P,x,n) Trim1(numer(P),x,n)/Trim1(denom(P),x,n): end: a22:=proc(n) local x,t,gu,i: if n>8 then RETURN(FAIL): fi: gu:=Data22s(x)[n]: gu:=subs({seq(x[i]=x[i]*t,i=1..n)},gu): gu:=taylor(gu,t=0,n+1): gu:=expand(coeff(gu,t,n)): for i from 1 to n do gu:=coeff(gu,x[i],1): od: gu: end: #U22i(n): Like U22(n) but using integration U22i:=proc(n) local x,gu,z,f,i: f:=f2n(n,x): gu:=(-1)^n*subs({seq(x[i]=-x[i],i=1..n)},f)*subs({seq(x[i]=z[i]/x[i],i=1..n)},f): for i from 1 to n do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n)): for i from 1 to n do gu:=int(gu,z[i]=0..infinity): od: gu: end: #Ursi(r,s,n): Like Urs(r,s,n) but using integration Ursi:=proc(r,s,n) local x,gu,z,f1,f2,i: f1:=fsn(r,n,x): f2:=fsn(s,n,x): gu:=(-1)^n*subs({seq(x[i]=-x[i],i=1..n)},f1)*subs({seq(x[i]=z[i]/x[i],i=1..n)},f2): for i from 1 to n do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n)): for i from 1 to n do gu:=int(gu,z[i]=0..infinity): od: gu: end: #U22nr(n,r): Applying the umbra that gives U22(n) out of f2n(n,x) to the more general f2nr(n,r,x). Try: #U22nr(5,5); U22nr:=proc(n,r) local x,f,i,P,c,j,s,p1: f:=f2nr(n,r,x): P:=partition(n+r): s:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c:=f: for j from 1 to nops(p1) do c:=coeff(c,x[p1[j][1]],p1[j][2]): od: s:=s+c^2*(-1)^(n+r)*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: s: end: #U22nri(n,r): U22nri:=proc(n,r) local x,gu,z,f,i: f:=f2nr(n,r,x): gu:=(-1)^(n+r)*subs({seq(x[i]=-x[i],i=1..n+r)},f)*subs({seq(x[i]=z[i]/x[i],i=1..n+r)},f): for i from 1 to n+r do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n+r)): for i from 1 to n+r do gu:=int(gu,z[i]=0..infinity): od: gu: end: #U22old(n): U22old:=proc(n) local x,f,i,i1,f1,j,s,c,p1: if n=1 then RETURN(1): elif n=2 then RETURN(2): else f:=f2n(n,x): s:=0: for i from 1 to nops(f) do f1:=op(i,f): p1:=[seq(degree(f1,x[i1]),i1=1..n)]: c:=normal(f1/mul(x[i1]^p1[i1],i1=1..n)): s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*s): fi: end: #U22(n): U22:=proc(n) local x,f,i,i1,f1,j,s,c,p1: if n=1 then RETURN(1): elif n=2 then RETURN(2): else f:=f2n(n,x): s:=0: for i from 1 to nops(f) do f1:=op(i,f): c:=op(1,f1): if not type(c,integer) then c:=1: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*s): fi: end: #Urr(r,n): Same as Urs(r,s,n) with s=r, but hopefully faster/ Urr:=proc(r,n) local x,f,i,i1,f1,j,s,c,p1: f:=fsn(r,n,x): if not type(f,`+`) then RETURN(Urs(r,r,n)): fi: s:=0: for i from 1 to nops(f) do f1:=op(i,f): c:=op(1,f1): if not type(c,integer) then c:=1: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*s): end: #U22seq(N): the first K terms of the sequence U22(n). Try: #U22seq(30); U22seq:=proc(K) local i: [seq(U22(i),i=1..K)]: end: #U22seqF(K): the first K terms of the sequence U22(n), by using the CONJECTURED recurrence. Try: #U22seqF(200); U22seqF:=proc(K) local n,N,ope: ope:=Ope22SBE(n,N): SeqFromRec(ope[1],n,N,ope[2],K): end: IsGoodrs:=proc(pi,r,s) local i: for i from 1 to nops(pi)-r do if pi[i+r]-pi[i]=s then RETURN(false): fi: od: true: end: #AnrsBF(n,r,s): Like Anrs(n,r,s) but by brute-force. try: #AnrsBF(6,1,1); AnrsBF:=proc(n,r,s) local mu,mu1,gu: mu:=permute(n): gu:={}: for mu1 in mu do if IsGoodrs(mu1,r,s) then gu:=gu union {mu1}: fi: od: gu: end: #IsGoodR(pi,R): is the permutation pi such that for each pair [r,s] in R and s in S pi[i+r]-pi[i]<>s. Try: #IsGoodR([1,2,3,4],{[2,2],[1,1]}); IsGoodR:=proc(pi,R) local pa: for pa in R do if not IsGoodrs(pi,pa[1],pa[2]) then RETURN(false): fi: od: true: end: #GO1s(s,n,N): Guesses an operator annihilationg Urs(1,s,n); Try: #GO1s(3,n,N); GO1s:=proc(s, n, N) local n1: subs(n = n - s, findrec([seq(Urs(1, s, n1), n1 = s .. s+10)], 1, 2, n, N)); end: #U1sF(s,n): Same as Urs(1,s,n) but much faster, using the recurrence. Try: #U1sF(5,20); U1sF:=proc(s,n) option remember: if n<0 then RETURN(0): fi: if (n<=s and n>=0) then RETURN(n!): fi: (n-1)*U1sF(s,n-1)+(n-1-s)*U1sF(s,n-2): end: #U1sEN(s,n): Same as Urs(1,s,n) but much faster, using Enrique Navarrete's formula #U1sEN(5,20); U1sEN:=proc(s,n) local j: add((-1)^j*binomial(n-s,j)*(n-j)!,j=0..n-s): end: #U11i(n): Like U11(n) but using integration. Try: U11(7); U11i:=proc(n) local x,gu,z,f,i: f:=f1n(n,x): gu:=(-1)^(n)*subs({seq(x[i]=-x[i],i=1..n)},f)*subs({seq(x[i]=z[i]/x[i],i=1..n)},f): for i from 1 to n do gu:=coeff(gu,x[i],0): od: gu:=gu*exp(-add(z[i],i=1..n)): for i from 1 to n do gu:=int(gu,z[i]=0..infinity): od: gu: end: #AnRbf(n,R): Given a set of pairs and a positive integer n, outputs the set of permutations such that for all the pairs [r,s] in R #all s in R pi[i+r]-pi[i]<>s. Try: #AnRbf(7,{[1,1],[2,2]}); AnRbf:=proc(n,R) local mu,mu1,gu,pair: mu:=permute(n): gu:={}: for mu1 in mu do if IsGoodR(mu1,R) then gu:=gu union {mu1}: fi: od: gu: end: ###New stuff #VrsSlow(r,s,n): The number of permutations of {1,...,n} such that |pi[i+r]-pi[r]|<>s VrsSlow:=proc(r,s,n) local x,f1,f2,i,P,j,su,p1,c1,c2,i1: f1:=fsn(r,n,x): f1:=subs({seq(x[i1]=2*x[i1],i1=2..n)},f1): f2:=fsn(s,n,x): P:=partition(n): su:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: su:=su+c1*c2*(-1)^n*(-1)^nops(P[i])*mul(p1[j][2]!,j=1..nops(p1)): od: su: end: #Vrr(r,n): Same as Vrs(r,s,n) with s=r, but hopefully faster/ Vrr:=proc(r,n) local x,f,i,i1,f1,j,s,c,p1: f:=fsn(r,n,x): if not type(f,`+`) then RETURN(Vrs(r,r,n)): fi: s:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c:=1: else c:=op(1,f1): if not type(c,integer) then c:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1))*2^add(p1[i1],i1=2..n): od: RETURN((-1)^n*s): end: #V11f(n): Same output as Vrr(1,n), but using the formula in Roberto Tauraso's article INTEGERS 6 (2006) #A11. Try: #[seq(V11f(i),i=1..20)]; V11f:=proc(n) local r,c: add((-1)^r*(n-r)!*add(2^c*binomial(r-1,c-1)*binomial(n-r,c),c=0..r),r=0..n-1):end: ###### #Urr(r,n): Same as Vrs(r,s,n) with s=r, but hopefully faster/ Urr:=proc(r,n) local x,f,i,i1,f1,j,s,c,p1: f:=fsn(r,n,x): if not type(f,`+`) then RETURN(Urs(r,r,n)): fi: s:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c:=1: else c:=op(1,f1): if not type(c,integer) then c:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: s:=s+c^2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*s): end: #Urs(r,s,n): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s Urs:=proc(r,s,n) local x,f1,i,j,su,p1,c1,c2,f,g,i1: g:=fsn(r,n,x): f:=fsn(s,n,x): if not type(f,`+`) then RETURN(UrsSlow(r,s,n)): fi: su:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c1:=1: else c1:=op(1,f1): if not type(c1,integer) then c1:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: c2:=g: for i1 from 1 to n do c2:=coeff(c2,x[i1],p1[i1]): od: su:=su+c1*c2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1)): od: RETURN((-1)^n*su): end: #Vrs(r,s,n): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s Vrs:=proc(r,s,n) local x,f1,i,j,su,p1,c1,c2,f,g,i1: g:=fsn(r,n,x): f:=fsn(s,n,x): if not type(f,`+`) then RETURN(VrsSlow(r,s,n)): fi: su:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c1:=1: else c1:=op(1,f1): if not type(c1,integer) then c1:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: c2:=g: for i1 from 1 to n do c2:=coeff(c2,x[i1],p1[i1]): od: su:=su+c1*c2*(-1)^(convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1))*2^add(p1[i1],i1=2..n): od: RETURN((-1)^n*su): end: #####START Generalized inclusion-exclusion with t #SinsR(pi,R): Given a permutation pi and a set of restrictions R finds how many places have pi[i+r]-pi[i] in R for [r,s] in R SinsR:=proc(pi,R) local n,i,co,r: n:=nops(pi): co:=0: for i from 1 to n do for r in R do if i+r[1]<=n and pi[i+r[1]]-pi[i]=r[2] then co:=co+1: fi: od: od: co: end: #AnRbfT(n,R,t): Given a set of pairs and a positive integer n, outputs the weight-enumerator according to the number of violations #Try: #AnRbfT(7,{[1,1],[2,2]},t); AnRbfT:=proc(n,R,t) local mu,mu1: mu:=permute(n): add(t^SinsR(mu1,R),mu1 in mu): end: #UrsSlowT(r,s,n,t): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s UrsSlowT:=proc(r,s,n,t) local x,f1,f2,i,P,j,su,p1,c1,c2: f1:=fsn(r,n,x): f2:=fsn(s,n,x): P:=partition(n): su:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: su:=expand(su+c1*c2*(t-1)^(n-nops(P[i]))*mul(p1[j][2]!,j=1..nops(p1))): od: su: end: #UrsT(r,s,n,t): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s UrsT:=proc(r,s,n,t) local x,f1,i,j,su,p1,c1,c2,f,g,i1: g:=fsn(r,n,x): f:=fsn(s,n,x): if not type(f,`+`) then RETURN(UrsSlowT(r,s,n,t)): fi: su:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c1:=1: else c1:=op(1,f1): if not type(c1,integer) then c1:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: c2:=g: for i1 from 1 to n do c2:=coeff(c2,x[i1],p1[i1]): od: su:=expand(su+c1*c2*(t-1)^(n-convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1))): od: su: end: #VrsSlowT(r,s,n,t): The number of permutations of {1,...,n} such that |pi[i+r]-pi[r]|<>s VrsSlowT:=proc(r,s,n,t) local x,f1,f2,i,P,j,su,p1,c1,c2,i1: f1:=fsn(r,n,x): f1:=subs({seq(x[i1]=2*x[i1],i1=2..n)},f1): f2:=fsn(s,n,x): P:=partition(n): su:=0: for i from 1 to nops(P) do p1:=ParToF(P[i]): c1:=f1: for j from 1 to nops(p1) do c1:=coeff(c1,x[p1[j][1]],p1[j][2]): od: c2:=f2: for j from 1 to nops(p1) do c2:=coeff(c2,x[p1[j][1]],p1[j][2]): od: su:=expand(su+c1*c2*(t-1)^(n-nops(P[i]))*mul(p1[j][2]!,j=1..nops(p1))): od: su: end: #VrsT(r,s,n,t): The number of permutations of {1,...,n} such that pi[i+r]-pi[r]<>s VrsT:=proc(r,s,n,t) local x,f1,i,j,su,p1,c1,c2,f,g,i1: g:=fsn(r,n,x): f:=fsn(s,n,x): if not type(f,`+`) then RETURN(VrsSlow(r,s,n)): fi: su:=0: for i from 1 to nops(f) do f1:=op(i,f): if not type(f1,`*`) then c1:=1: else c1:=op(1,f1): if not type(c1,integer) then c1:=1: fi: fi: p1:=[seq(degree(f1,x[i1]),i1=1..n)]: c2:=g: for i1 from 1 to n do c2:=coeff(c2,x[i1],p1[i1]): od: su:=expand(su+c1*c2*(t-1)^(n-convert(p1,`+`))*mul(p1[j]!,j=1..nops(p1))*2^add(p1[i1],i1=2..n)): od: RETURN(su): end: #a11G(n): Same as V11f(n) or Vrr(1,n) but using George Spahn's scheme. Try: #[seq(a11G(n),n=1..30)]; a11G:=proc(n) option remember: if n=0 then 1 elif n=1 then 1: else (n-2)*a11G(n-1)+b11G(n-1)-c11G(n-1): fi: end: b11G:=proc(n) option remember: if n=0 then 0 elif n=1 then 0: else 2*(n-1)*a11G(n-1)+2*b11G(n-1)+b11G(n-2): fi: end: c11G:=proc(n) option remember: if n=0 then 0 elif n=1 then 0: else 2*a11G(n-1)+c11G(n-1): fi: end: ###start Robbins approach #f2nR(n,x): the weight-enumerator of tilings of the form (1X)* f2nR:=proc(n,x) option remember: f2nrR(n,0,x): end: #f2nrR(n,r,x): The weight enumerator of tilings of the board 1^n(X1)^r f2nrR:=proc(n,r,x) local i: option remember: if n=0 and r=0 then RETURN(1): fi: if r>0 then RETURN(expand(f2nrR(n,r-1,x)+add(2*x^(i-1)*f2nrR(n,r-i,x),i=2..r)+ add(2*x^(r+i-1)*f2nrR(n+1-2*i,i-1,x),i=1..trunc((n+1)/2)))): elif r=0 then RETURN(expand(f2nrR(n-1,0,x)+ add(2*x^(i-1)*f2nrR(n+1-2*i,i-1,x),i=2..trunc((n+1)/2)))): else RETURN(FAIL): fi: end: #f2nRslow(n,x): The weight enumerator of subsets of {1,...,n-2} with weight(S):= x^|S|*2^NumberOf2-ComponentsOf[S}. Try: #f2nRslow(10,x); f2nRslow:=proc(n,x) local gu,i: gu:=f2n(n,x): expand(subs({x[1]=1,seq(x[i]=x^(i-1)*2,i=2..n)},gu)): end: #V12Rslow(n): Same as Vrs(1,2,n) but using Dave Robbins' method. In other words the number of permutations of {1,...,n} such that |pi[i+1]-pi[i]|<>2. Try: #[seq(V12Rslow(i),i=1..10)]; V12Rslow:=proc(n) local gu,x,i: gu:=f2nRslow(n,x): gu:=subs(x=-x,gu): add(coeff(gu,x,i)*(n-i)!,i=0..degree(gu,x)): end: #V12R(n): Same as Vrs(1,2,n) but using Dave Robbins' method. In other words the number of permutations of {1,...,n} such that |pi[i+1]-pi[i]|<>2. Try: #[seq(V12R(i),i=1..10)]; V12R:=proc(n) local gu,x,i: gu:=f2nR(n,x): gu:=subs(x=-x,gu): add(coeff(gu,x,i)*(n-i)!,i=0..degree(gu,x)): end: #fsnRslow(s,n,x): The weight enumerator of subsets of {1,...,n-2} with weight(S):= x^|S|*2^NumberOfComponentsOf[S}. Try: #fsnRslow(3,10,x); fsnRslow:=proc(s,n,x) local gu,i: gu:=fsn(s,n,x): expand(subs({x[1]=1,seq(x[i]=x^(i-1)*2,i=2..n)},gu)): end: #V1sRslow(s,n): Same as Vrs(1,s,n) but using Dave Robbins' method. In other words the number of permutations of {1,...,n} such that |pi[i+1]-pi[i]|<>s. Try: #[seq(V1sRslow(i),i=1..10)]; V1sRslow:=proc(s,n) local gu,x,i: gu:=fsnRslow(s,n,x): gu:=subs(x=-x,gu): add(coeff(gu,x,i)*(n-i)!,i=0..degree(gu,x)): end: #WtTilR(R,T,x): Given a set of integers and a set of tiles, outputs the setight-enumerator of all the set of all tilings with the weight #of a singleton being 1 and the weight of a larger tile being x^(nops(t)-1)*y. Needed for V1sR(s,n). Try: #WtTilR({1,2,3,4,5,6},{{1,3},{1,3,5},{1}},x,y); WtTilR:=proc(R,T,x) local i,S,t,R1,t1: option remember: if R={} then RETURN(1): fi: i:=min(op(R)): S:=0: for t in T do if {seq(i-1+t1,t1 in t)} minus R={} then R1:=R minus {seq(i-1+t1,t1 in t)}: if nops(t)=1 then S:=expand(S+WtTilR(R1,T,x)): else S:=expand(S+x^(nops(t)-1)*2*WtTilR(R1,T,x)): fi: fi: od: S: end: #fsnR(s,n,x): fsn(n,x) for V1sR(s,n): fsnR:=proc(s,n,x) local R,T,i,r: R:={seq(i,i=1..n)}: T:={}: for r from 0 to trunc((n-1)/s) do T:=T union {{seq(1+s*i,i=0..r)}}: od: WtTilR(R,T,x): end: #V1sR(s,n): Same as Vrs(1,s,n) but using Dave Robbins' method. In other words the number of permutations of {1,...,n} such that |pi[i+1]-pi[i]|<>s. Try: #[seq(V1sR(3,i),i=1..10)]; V1sR:=proc(s,n) local gu,x,i: gu:=fsnR(s,n,x): gu:=subs(x=-x,gu): add(coeff(gu,x,i)*(n-i)!,i=0..degree(gu,x)): end: #Ker1s(s,x,z,K): The rational function in x and z, let's call it F(x,z) such that Vrs(1,s,n) equals the coeff. of z^n in int(exp(-x)*F(x,z),x=0..infinity). #K is a guessing paramater/ #Try: #Ker1s(1,x,z,50); Ker1s:=proc(s,x,z,K) local K1,gu,i,lu,mu: gu:=[seq(fsnR(s,i,x),i=0..20)]: for K1 from 30 by 10 to K do lu:=GuessRec(gu): if lu<>FAIL then mu:=CtoR(lu,z): mu:=subs(x=-x,mu): mu:=factor(subs({x=1/x,z=z*x},mu)): RETURN(mu): else gu:=[seq(fsnR(s,i,x),i=0..K1)]: fi: od: FAIL: end: #CheckKer1s(s,K,N1): Checks Ker1s(s,x,z,K). Try: #CheckKer1s(3,50,20); CheckKer1s:=proc(s,K,N1) local z,x,F,F1,mu,i,gu: F:=Ker1s(s,x,z,K): if F=FAIL then print(`Make `, K, `bigger `): RETURN(FAIL): fi: F1:=taylor(F,z=0,N1+2): mu:=expand([seq(coeff(F1,z,i),i=1..N1)]): mu:=[seq(int(exp(-x)*mu[i],x=0..infinity),i=1..N1)]: gu:=[seq(Vrs(1,s,i),i=1..N1)]: evalb(mu=gu): end: #DE1s(s,z,Dz,K): The differential operator in z, Dz, such that applying it to the generating function Sum(Vsr(1,s,i)*z^i,i=0..infinity) is a polynomial of low degree. #K is a guessing parameter.Try: #DE1s(2,z,Dz,50): DE1s:=proc(s,z,Dz,K) local F,x,A,ope,yemin,lu: F:=Ker1s(s,x,z,K): if F=FAIL then print(`Make `, K, ` bigger `): RETURN(FAIL): fi: lu:=AZcI(exp(-x)*F,x,z,Dz,0,A)[1]: ope:=lu[1]: yemin:=limit(lu[2],A=infinity): [ope,yemin]: end: #REC1s(s,n,N,K): The recurrence operator in n, N, annihilating Vsr(1,s,n) #K is a guessing parameter.Try: #REC1s(2,n,N,50): REC1s:=proc(s,n,N,K) local z,Dz,gu,gu1,gu11,i,j,ope,d,j1,INI: gu:=DE1s(s,z,Dz,K): if degree(gu[2],z)=FAIL then print(`Not yet implemented `): RETURN(FAIL): fi: gu:=gu[1]: ope:=0: for i from 0 to degree(gu,z) do gu1:=coeff(gu,z,i): for j from 0 to degree(gu1,Dz) do gu11:=coeff(gu1,Dz,j): ope:=ope+ gu11*mul(n-i+j1,j1=1..j)*N^(j-i): od: od: ope:=expand(ope): d:=ldegree(ope,N): ope:=expand(subs(n=n-d,ope))/N^d: ope:=add(factor(coeff(ope,N,i))*N^i,i=0..degree(ope,N)): ope:=Yafe(ope,N)[2]: ope:=add(factor(coeff(ope,N,i))*N^i,i=0..degree(ope,N)): ope:=collect(ope,N): INI:=[seq(V1sR(s,i),i=1..degree(ope,N))]: if SeqFromRec(ope,n,N,INI,nops(INI)+10)<>[seq(V1sR(s,i),i=1..nops(INI)+10)] then print(`Something went wrong `): RETURN(FAIL): fi: [ope,INI]: end: #Paper1s(s,n,K,M1,M2): An article about the sequence enumerating permutations of {1,..,n} such that |pi[i+1]-pi[i]|<>s. Try: #Paper1s(1,n,20,100,1000): Paper1s:=proc(s,n,K,M1,M2) local F,x,z,F1,lu,G,ku,N,ope,INI,a,Dz,i,t0: t0:=time(): F:=Ker1s(s,x,z,K): if F=FAIL then RETURN(FAIL): fi: print(``): print(`-----------------------------------------------------`): print(``): print(`Efficient computation of the sequence enumerating permutations of {1,...,n} such that |pi[i+1]-pi[i]| is never`, s): print(``): print(`By Shalosh B. Ekhad `): print(``): print(`Let a(n) be the number of permutations pi of {1,...,n} such that |pi[i+1]-pi[i]| is never`, s): print(``): print(`Then `): print(``): print(Sum(a(n)*z^n,n=0..infinity)=Int(exp(-x)*F,x=0..infinity)): print(``): print(`and in Maple notation `): print(``): lprint(Sum(a(n)*z^n,n=0..infinity)=Int(exp(-x)*F,x=0..infinity)): print(``): F1:=taylor(F,z=0,M1+1): print(``): print(`This can be used for fast computation of this sequence. Let's do it from n=1 to n=`, M1): print(``): lprint([seq(int(expand(coeff(F1,z,i))*exp(-x),x=0..infinity),i=1..M1)]): print(``): lu:=DE1s(s,z,Dz,K): print(``): print(`Using the Continuous Almkvist-Zeilberger algorithm it follows that the generating function,let's call it G(z), namely`): print(``): print(G(z)=Sum(a(n)*z^n,n=0..infinity)): print(``): print(`satisdies the following differential equation `): print(``): print(coeff(lu[1],Dz,0)*G(z)+add(factor(coeff(lu[1],Dz,i))*diff(G(z),z$i),i=1..degree(lu[1],Dz))=lu[2]): print(``): print(`and in Maple notation `): print(``): lprint(coeff(lu[1],Dz,0)*G(z)+add(factor(coeff(lu[1],Dz,i))*diff(G(z),z$i),i=1..degree(lu[1],Dz))=lu[2]): print(``): print(`This translates to a linear recurrence equation with polynomial coefficients for the sequence`): print(``): if s>2 then print(`Since the right side is not a polynomial, we decided not to compute the resulting (inhomog.) recurrence, that of course implies a higher-order homog. recurrence`): elif s=1 or s=2 then ku:=REC1s(s,n,N,K): print(``): ope:=ku[1]: INI:=ku[2]: print(add(coeff(ope,N,i)*a(n+i),i=0..degree(ope,N))=0 ): print(``): print(`and in Maple notation `): print(``): lprint(add(coeff(ope,N,i)*a(n+i),i=0..degree(ope,N))=0 ): print(``): print(`Subject to the initial conditions `): print(``): lprint(seq(a(i)=INI[i],i=1..nops(INI))): print(``): print(`Using this recurrence we can compute many terms. Just for fun we have `): print(``): print(a(M2)=SeqFromRec(ope,n,N,INI,M2)[M2]): print(``): fi: print(``): print(`-----------------------------------------------`): print(``): print(`This ends this paper that took`, time()-t0, `seconds to generate. `): print(``): end: #Start Maple implementation of Rintaro Matsuo's algoritm #U22rm(n): Same as U22(n) (q.v.) but using Rintaro Matsuo's brilliant O(n^4) algorithm. Try: #U22rm(10); U22rm:=proc(n) matsuo(n,n,floor((n+1)/2),0,0,1) + matsuo(n,n,floor((n+1)/2),0,0,2) + matsuo(n,n,floor((n+1)/2),0,0,3) end: #matsuo(n,i,j,k,l,m): The quantity needed in Matsuo Rintaro's approach for computing U22(n), used in the more efficient procedure U22rm(n). Try: #matsuo(5,1,1,1,1,1); matsuo:=proc(n, i, j, k, l, m) local half, iiminus1ok, su: option remember: if i = 0 or j < 0 or j > i or k < 0 or (k >= j and k > 0) or l < 0 or (l >= (i - j) and l > 0) then #print(n,i,j,k,l,m, "suzero"): return(0): fi: if i = 1 and j=0 and k=0 and l=0 and m=3 then #print(n,i,j,k,l,m, "su", 1): return(1): fi: if i = 1 and j=1 and k=0 and l=0 and m=2 then #print(n,i,j,k,l,m, "su", 1): return(1): fi: half := floor((n+1)/2): iiminus1ok := 0: if i-1 = half then iiminus1ok := 1: fi: su := 0: if m = 1 then su := su + matsuo(n, i-1, j-1, k, l, 1)* (j-1-k -1+iiminus1ok) : su := su + matsuo(n, i-1, j-1, k, l, 2)* (j-1-k): su := su + matsuo(n, i-1, j-1, k, l, 3)* (j-1-k): su := su + matsuo(n, i-1, j-1, k-1, l, 1)* (1-iiminus1ok) : su := su + matsuo(n, i-1, j-1, k+1, l, 1)* (k+1) : su := su + matsuo(n, i-1, j-1, k+1, l, 2)* (k+1): su := su + matsuo(n, i-1, j-1, k+1, l, 3)* (k+1): fi: if m = 2 then su := su + matsuo(n, i-1, j-1, k, l, 1)* 1 : su := su + matsuo(n, i-1, j-1, k, l, 2)* (iiminus1ok): su := su + matsuo(n, i-1, j-1, k, l, 3)* 1 : su := su + matsuo(n, i-1, j-1, k-1, l, 2)* (1-iiminus1ok) : fi: if m = 3 then su := su + matsuo(n, i-1, j, k, l, 1)* (i-j-l) : su := su + matsuo(n, i-1, j, k, l, 2)* (i-j-l): su := su + matsuo(n, i-1, j, k, l, 3)* (i-j-l -1+iiminus1ok): su := su + matsuo(n, i-1, j, k, l-1, 3)* (1-iiminus1ok) : su := su + matsuo(n, i-1, j, k, l+1, 1)* (l+1) : su := su + matsuo(n, i-1, j, k, l+1, 2)* (l+1): su := su + matsuo(n, i-1, j, k, l+1, 3)* (l+1): fi: #print(n,i,j,k,l,m, "su", su): su: end: #End Maple implementation of Rintaro Matsuo's algoritm #RM2(pi): The Rintaro Matsuo transform that takes members of Anrs(n,2,2) to The set of permutations avoiding [i,i+1] for i=1...n-1 except that [floor((n+1)/2),floor((n+1)/2)+1] is allowed and #it is ok to have a violation at the pair of locations [floor((n+1)/2),floor((n+1)/2)+1]. Try: #RM2([1,4,3,5,2]); RM2:=proc(pi) local n,i,R,m,Rinv: n:=nops(pi): if n mod 2=0 then m:=n/2: R:=[seq( op([i,i+m]),i=1..n/2)]: else m:=(n+1)/2: R:=[seq( op([i,i+m]),i=1..m-1),m]: fi: for i from 1 to n do Rinv[R[i]]:=i: od: [seq(R[pi[Rinv[i]]],i=1..n)]: end: #A22rm(n): The image of Anrs(n,2,2) under the Rintaro Matsuo transform. Try: #A22rm(7); A22rm:=proc(n) local gu,gu1: gu:=Anrs(n,2,2): {seq(RM2(gu1),gu1 in gu)}: end: #IsRM2(pi): Is the permutation good in the sense of Rintaro Matsuo? IsRM2:=proc(pi) local n,i,m: n:=nops(pi): m:=trunc((n+1)/2): for i from 1 to n-1 do if i<>m then if (pi[i]<>m and pi[i+1]-pi[i]=1) then RETURN(false): fi: fi: od: true: end: #A22rmBF(n): Same as A22rm(n) but by brute force, for checking only. Try: #evalb(A22rm(6)=A22rmBF(6)); A22rmBF:=proc(n) local mu,mu1,gu: mu:=permute(n): gu:={}: for mu1 in mu do if IsRM2(mu1) then gu:=gu union {mu1}: fi: od: gu: end: #B22rm(n): The image of Bnrs(n,2,2) under the Rintaro Matsuo transform. Try: #B22rm(7); B22rm:=proc(n) local gu,gu1: gu:=BnrsBF(n,2,2): {seq(RM2(gu1),gu1 in gu)}: end: #IsRM2a(pi): Is the permutation good in the sense of Rintaro Matsuo, and with absolute value? IsRM2a:=proc(pi) local n,i,m: n:=nops(pi): m:=trunc((n+1)/2): for i from 1 to n-1 do if i<>m then if (pi[i]<>m and abs(pi[i+1]-pi[i])=1) then RETURN(false): fi: fi: od: true: end: #B22rmBF(n): Same as B22rm(n) but by brute force, for checking only. Try: #evalb(B22rm(6)=B22rmBF(6)); B22rmBF:=proc(n) local mu,mu1,gu: print(`B22rmBF(n): UNDER CONSTRUCTION. DOES NOT YET WORK`): mu:=permute(n): gu:={}: for mu1 in mu do if IsRM2a(mu1) then gu:=gu union {mu1}: fi: od: gu: end: #RMperm(n,g): The Rintaro Matsuo's g-permutation. Try: #RMperm(11,2); RMperm:=proc(n,g) local i,cu,SE,R: R[1]:=1: SE:={seq(i,i=2..n)}: cu:=1: for i from 2 while SE<>{} do cu:=cu+g: if not member(cu,SE) then cu:=min(op(SE)): fi: R[cu]:=i: SE:=SE minus {cu}: od: [seq(R[i],i=1..n)]: end: #inv1(pi): The inverse of the permutation pi inv1:=proc(pi) local T,i,n: n:=nops(pi): for i from 1 to n do T[pi[i]]:=i: od: [seq(T[i],i=1..n)]: end: #RMg(pi,g1,g2): The Generalized Rintaro Matsuo transform that takes members of Anrs(n,g1,g2) to The some set of permutations tbd #RMg(pi,2,2) is the same as RM2(pi). Try: #pi:=randperm(12); RMg(pi,3,3); RMg:=proc(pi,g1,g2) local n,R,Rinv,i: n:=nops(pi): R:=RMperm(n,g1): Rinv:=inv1(RMperm(n,g2)): [seq(R[pi[Rinv[i]]],i=1..n)]: end: #AnrsRM(n,r,s): The image of Anrs(n,r,s) under the Rintaro Matsuo transform pi->RMg(pi,r,s). Try: #AnrsRM(7,2,2); It should be the same as A22rm(7); AnrsRM:=proc(n,r,s) local gu,gu1: gu:=Anrs(n,r,s): {seq(RMg(gu1,r,s),gu1 in gu)}: end: #FindVi1(pi): Given a permutation pi finds the list of pairs [i,j] such that pi[i]=j and pi[i+1]=j+1. Try: #FindVi1([1,2,3,4]); FindVi1:=proc(pi) local gu,i: gu:={}: for i from 1 to nops(pi)-1 do if pi[i+1]=pi[i]+1 then gu:=gu union {[i,pi[i]]}: fi: od: gu: end: #FindVr(n,r): All the violations [i,j] in the members of AnrsRM(n,r,r) minus Anrs(n,1,1). Try: #FindV(7,2); FindVr:=proc(n,r) local gu,gu1: gu:=AnrsRM(n,r,r) minus Anrs(n,1,1): {seq(op(FindVi1(gu1)), gu1 in gu)}: end: #RIN(n,a,b): The number of permutations pi of {1,...,n} such that pi[i+1]-pi[i]<>1 except possibly when pi[i]=b and i=a using EXCLUSION EXCLUSION. Try: #RIN(10,5,5). RIN:=proc(n,a,b) local k1,k2: option remember: if a>0 and a=0 and a<=n and b>=1 and b<=n then RETURN(FAIL): fi: if a>0 and a0 then RETURN(0): fi: if b0 then RETURN(0): fi: if b1 except possibly when pi[i]=b and i=a using EXCLUSION EXCLUSION. Try: #RINv(10,5,5). RINv:=proc(n,a,b) local k1,k2: option remember: if a>0 and a=0 and a<=n and b>=1 and b<=n then RETURN(FAIL): fi: if a>0 and a0 then RETURN(0): fi: if b0 then RETURN(0): fi: if b i or k < 0 or (k >= j and k > 0) or l < 0 or (l >= (i - j) and l > 0) then #print(n,i,j,k,l,m, "sumzero"): return(0): fi: if i = 1 and j=0 and k=0 and l=0 and m <= 2 then #print(n,i,j,k,l,m, "sum", 1): return(1): fi: if i = 1 and j=1 and k=0 and l=0 and m >2 then #print(n,i,j,k,l,m, "sum", 1): return(1): fi: if i = 2 and j = 0 then if k <> 0 then return(0): fi: if m <> 3 then return(0): fi: if l=1 then return(2): fi: if l=0 then return(0): fi: fi: if i = 2 and j = 1 then if k <> 0 or l <> 0 then return(0): fi: if m=1 or m = 3 then return(0): fi: return(1): fi: if i = 2 and j = 2 then if l <> 0 then return(0): fi: if m <> 1 then return(0): fi: if k=1 then return(2): fi: if k=0 then return(0): fi: fi: half := floor((n+1)/2): iim1ok := false: if i-1 = half then iim1ok := true: fi: im1im2ok := false: if i-2 = half then im1im2ok := true: fi: sum := 0: if iim1ok then #cannot introduce any new mistakes if m=1 then sum := sum + v22rm1 (n, i-1, j-1, k, l, 1)* (1) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 2)* (2) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 1)* (1) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 2)* (0) : elif m=2 then sum := sum + v22rm1 (n, i-1, j-1, k, l, 1)* ( max(j-1-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 2)* ( max(j-2-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 3)* ( max(j-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 4)* ( max(j-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 1)* ( k ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 2)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 3)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 4)* ( k+1 ) : elif m=3 then sum := sum + v22rm1 (n, i-1, j, k, l, 3)* (1) : sum := sum + v22rm1 (n, i-1, j, k, l, 4)* (2) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 3)* (1) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 4)* (0) : else sum := sum + v22rm1 (n, i-1, j, k, l, 1)* ( max(i-j-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 2)* ( max(i-j-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 3)* ( max(i-j-1-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 4)* ( max(i-j-2-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 1)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 2)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 3)* ( l ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 4)* ( l+1 ) : fi: elif im1im2ok then # can introduce new mistake, cannot remove one specific mistake if m=1 then sum := sum + v22rm1 (n, i-1, j-1, k-1, l, 1)* (2) : sum := sum + v22rm1 (n, i-1, j-1, k-1, l, 2)* (2) : elif m=2 then sum := sum + v22rm1 (n, i-1, j-1, k, l, 1)* ( max(j-2-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 2)* ( max(j-2-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 3)* ( max(j-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 4)* ( max(j-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 1)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 2)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 3)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 4)* ( k+1 ) : elif m=3 then sum := sum + v22rm1 (n, i-1, j, k, l-1, 3)* (2) : \ sum := sum + v22rm1 (n, i-1, j, k, l-1, 4)* (2) : else sum := sum + v22rm1 (n, i-1, j, k, l, 1)* ( max(i-j-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 2)* ( max(i-j-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 3)* ( max(i-j-2-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 4)* ( max(i-j-2-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 1)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 2)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 3)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 4)* ( l+1 ) : fi: else if m=1 then sum := sum + v22rm1 (n, i-1, j-1, k, l, 1)* (1) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 2)* (0) : sum := sum + v22rm1 (n, i-1, j-1, k-1, l, 1)* (1) : sum := sum + v22rm1 (n, i-1, j-1, k-1, l, 2)* (2) : elif m=2 then sum := sum + v22rm1 (n, i-1, j-1, k, l, 1)* ( max(j-1-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 2)* ( max(j-2-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 3)* ( max(j-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k, l, 4)* ( max(j-k,0) ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 1)* ( k ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 2)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 3)* ( k+1 ) : sum := sum + v22rm1 (n, i-1, j-1, k+1, l, 4)* ( k+1 ) : elif m=3 then sum := sum + v22rm1 (n, i-1, j, k, l, 3)* (1) : sum := sum + v22rm1 (n, i-1, j, k, l, 4)* (0) : sum := sum + v22rm1 (n, i-1, j, k, l-1, 3)* (1) : sum := sum + v22rm1 (n, i-1, j, k, l-1, 4)* (2) : else sum := sum + v22rm1 (n, i-1, j, k, l, 1)* ( max(i-j-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 2)* ( max(i-j-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 3)* ( max(i-j-1-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l, 4)* ( max(i-j-2-l,0) ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 1)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 2)* ( l+1 ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 3)* ( l ) : sum := sum + v22rm1 (n, i-1, j, k, l+1, 4)* ( l+1 ) : fi: fi: sum: end: #V22seqF(N): Same as [seq(Vrr(2,n),n=1..N)], but using the CONJECTURED recurrence, hence this is "semi-rigorous". Try: #V22seqF(20); V22seqF:=proc(N) local ope,n,Sn: ope:=Ope22KKv(n,Sn): SeqFromRec(ope[1],n,Sn,ope[2],N): end: