###################################################################### ## SMCramsey: Save this file as SMCramsey. To use it, stay in the # ## same directory, get into Maple (by typing: maple ) # ## and then type: read SMCramsey : # ## Then follow the instructions given there # ## # ## Written by Doron Zeilberger, Rutgers University , # # zeilberg at math dot rutgers dot edu # ###################################################################### #Created: Jan. 2004 #This version: Nov. 2005 #SMCramsey: A Maple package to implement Doron Zeilberger's #Symbolic Moment Calculus applied to Graph Colorings #Please report bugs to zeilberg at math dot rutgers dot edu print(`Created: Jan. 2004`): print(`This version: Nov. 2005`): print(`This is SMramsey, a Maple package to implement`): print(`Doron Zeilberger's Symbolic Moment Calculus`): print(`applied to Graph Colorings `): print(`as it is described in his article`): print(` Symbolic Moment Calculus: II. Why is Ramsey Theory Soooo`): print(` Eeeeenormoulsy Hard?`): print(``): print(`Written by Doron Zeilberger, zeilberg at math dot rutgers dot edu`): lprint(``): print(`Please report bugs to zeilberg at math dot rutgersd dot edu`): lprint(``): print(`The most current version of this package and paper`): print(` are available from`): print(`http://www.math.rutgers.edu/~zeilberg/`): print(`For a list of the MAIN procedures type ezra();, for help with`): print(`a specific procedure, type ezra(procedure_name);`): print(``): with(combinat): ezra1:=proc(): if args=NULL then print(`All the procedures are:`): print(` CheckkthMoment, CheckkthMomentS,FkL, GrL `): print(`kthMoment, kthMomentE, kthMomentS , Kurtosis`): print(` Mishkal , ThirdMoment, ThirdMomentT, Variance`): else ezra(args): fi: end: ezra:=proc() if args=NULL then print(`kthMoment, kthMomentE, kthMomentS , Kurtosis`): fi: if nargs=1 and args[1]=Bonf1 then print(`Bonf1(k,r,n): The first Bonferonni Expression for the number`): print(`of monochromatic K_k in the r-coloring of the edges of K_n`): fi: if nargs=1 and args[1]=Bonf3 then print(`Bonf3(k,r,n): The third Bonferonni Expression for the number`): print(`of monochromatic K_k in the r-coloring of the edges of K_n`): print(`k is a numeric integer, but r and n are symbolic (or numeric)`): fi: if nargs=1 and args[1]=CheckkthMoment then print(`CheckkthMoment(a,k,n0): Check the procedure kthMoment by`): print(`for the special case of 2 colors, and n=n0, for the k-th`): print(`moment of the r.v. number of monochromatic K_a `): print(`For example, type: CheckkthMoment(3,2,4)`): fi: if nargs=1 and args[1]=CheckkthMomentS then print(`CheckkthMomentS(a0,k,r,N,T): checks procedure kthMomentS`): print(`by comparing it with numerical instances (a=a0) of`): print(` kthMoment(a,k,r,N). `): print(`For example, try: CheckkthMomentS(4,2,r,N,3);`): fi: if nargs=1 and args[1]=FkL then print(`FkL(k,L) : The performance of Bonf1 vs. Bonf3 for #colors 2 `): print(`through L with monochromatic K_k (k numeric)`): fi: if nargs=1 and args[1]=GrL then print(`GrL(r,L) : The performance of Bonf1 vs. Bonf3 for r colors `): print(`with monochromatic K_k (k 3 through L`): fi: if nargs=1 and args[1]=kthMoment then print(`kthMoment(a,k,r,N): The k^th moment (about the mean) of `): print(`the r.v. #monochromatic K_{a} in edge-r-colorings of K_N`): print(`For example, try kthMoment(3,3,r,N);`): fi: if nargs=1 and args[1]=kthMomentE then print(`kthMomentE(a,k,n0) : If a and k , and n0, are positive intgers,`): print(`computes empirically the kth moment, about the mean)`): print(` of the r.v. "number of mono-chromatic `): print(`K_a over S_{n_0} in two-colorings of all the labelled graphs on n0`): print(`vertices"; For example, try: kthMomentE(3,2,4);`): fi: if nargs=1 and args[1]=kthMomentS then print(`kthMomentS(a,k,r,N,T) : The first T+1 terms in the asymp. expression`): print(`for the k-th moment (about the mean) of the r.v. "Number of mono-`): print(`chromatic K_a in r-edge-colorings of the K_N`): print(`divided by ((N^a/a!/r^(binomial(a,2)-1))^k `): print(` a, r, N are symbolic, k and T are numeric integers. `): print(`For example, try: kthMomentS(a,2,r,N,3);`): fi: if nargs=1 and args[1]=Kurtosis then print(`Kurtosis(k,r,N,T): The first T terms in the asymptotic expression`): print(`in N (with k fixed but symbolic)`): print(`for the Kurtosis of the r.v. #of mono-chrmomatic K_k in`): print(`r-colorings of the edges of K_N (r also symbolic)`): fi: if nargs=1 and args[1]=Mishkal then print(`Mishkal(L,r,n): Given a list of integers, L`): print(`of length k, say, finds`): print(`E[X_L[1]*...*X_L[k]]`): print(`where X_i is the number of monochromatic r-edge-colorings K_i`): print(`For example, try Mishkal([3,3,3],r,n);`): fi: if nargs=1 and args[1]=ThirdMoment then print(`ThirdMoment(k,r,n):The third moment (about the mean)`): print(`of of the r.v. #number of mono-chromatic K_k `): print(`for r-colorings of the edges of K_n (k numeric, r,n symbolic)`): print(`Try ThirdMoment(3,r,n);`): fi: if nargs=1 and args[1]=ThirdMomentT then print(`ThirdMomentT(a,r,n,T) The first T leading terms (in n)`): print(`of ThirdMoment(a,r,n) `): print(` r, and n symbolic, a numeric) For example, `): print(` try ThirdMomentT(3,r,n,1);`): fi: if nargs=1 and args[1]=Variance then print(`Variance(k,r,n):The second moment (about the mean)`): print(`of of the r.v. #number of mono-chromatic K_k `): print(`for r-colorings of the edges of K_n (k numeric, r,n symbolic)`): print(`Try Variance(3,r,n);`): fi: end: Tupe:=proc(Resh,S) local gu,i,GU,lu: GU:={}: gu:=Tup(Resh,S): for i from 1 to nops(gu) do lu:=gu[i]: lu:=Yofi(lu): GU:=GU union {lu}: od: GU: end: Yofi:=proc(Resh) local i: [seq(sort(convert(Resh[i],list)),i=1..nops(Resh))]: end: #Tup(Resh,S): Given a list of integers Resh, let's call #[k_1, ..., k_r], and a set S #finds all tuples [S_1, ..., S_r], #where S_1, ..., S_r are subsets of S and their union #equals S and nops(S_i)=k_i Tup:=proc(Resh,S) local i,gu,r,mu,lu,j,Resh1,gu1,k: option remember: r:=nops(Resh): if r=1 then if nops(S)=Resh[1] then RETURN({[S]}): else RETURN({}) fi: fi: Resh1:=[op(1..r-1,Resh)]: lu:=Subsets1(S,Resh[r]): gu:={}: for i from 1 to nops(lu) do mu:=powerset(lu[i]): for j from 1 to nops(mu) do gu1:= Tup(Resh1,S minus mu[j]): gu:=gu union {seq([op(gu1[k]),lu[i]],k=1..nops(gu1))}: od: od: gu: end: #Subsets1(S,k): all the k-subsets of S Subsets1:=proc(S,k) local i,gu1,gu2,S1,N: option remember: if k=0 then RETURN({{}}): fi: if nops(S)=0 or nops(S)1 GraphE:=proc(ListS) local i1,i2,G: G:={}: for i1 from 1 to nops(ListS) do for i2 from i1+1 to nops(ListS) do if nops(ListS[i1] intersect ListS[i2])>1 then G:=G union {{i1,i2}}: fi: od: od: G: end: #Neighs1(G,v): all the neighbors of the vertex v in G (including v) #where V is the set of edges Neighs1:=proc(G,V,v) local gu,i1,V1: gu:={v}: V1:= V minus {v}: for i1 from 1 to nops(V1) do if member({V1[i1],v},G) then gu:=gu union {V1[i1]}: fi: od: gu: end: #Neighs(G,V,SetV): the collective neigbors of the vertices in SetV Neighs:=proc(G,V,SetV) local i1: {seq(op(Neighs1(G,V,SetV[i1])),i1=1..nops(SetV))}: end: #ConComp1(G,V,v): The connected component of vertex v ConComp1:=proc(G,V,v) local gu,gu1: gu:={v}: gu1:=Neighs(G,V,gu): while gu1<>gu do gu:=gu1: gu1:=Neighs(G,V,gu1): od: gu1: end: #ConComps(G,V): the connected components of the graph G,V ConComps:=proc(G,V) local i1: {seq(ConComp1(G,V,V[i1]),i1=1..nops(V))}: end: #Wt(ListS,r): The weight of the list of sets ListS, in edge-r-coloring Wt:=proc(ListS,r) local G,i1: G:=GraphE(ListS): r^(nops(ConComps(G,{seq(i1,i1=1..nops(ListS))}))-nops(Edges(ListS))): end: #Mishkal(L,r,n): Given a list of integers, L #of length k, say, finds #E[X_L[1]*...*X_L[k]] #where X_i is the number of monochromatic r-edge-colorings K_i #For example, try Mishkal([3,3,3],r,n); MishkalSlow:=proc(mu,r,n) local gu,i,t,j,mish,mish1,katan,gadol: option remember: katan:=max(op(mu)): gadol:=convert(mu,`+`): mish:=0: for t from katan to gadol do gu:=Tup(mu,{seq(i,i=1..t)}): mish1:=0: for j from 1 to nops(gu) do mish1:=mish1+Wt(gu[j],r): od: mish:=expand(mish+expand(binomial(n,t))*mish1): od: factor(expand(mish)): end: #MishkalT0(Resh1): A quick way to compute #the (coeff. of the) leading contribution to Mishkal(Resh1,n) #i.e. Mishkal(Resh1,n,sum(Resh1[i],i=1..nops(Resh1))) #For example, try MishkalT0([2,2,2]); MishkalT0:=proc(Resh1): MTC(op(Resh1)):end: MTC:=proc() local gu,i: gu:=[args]: for i from 1 to nops(gu) do if gu[i]<0 then ERROR(`All arguments must be non-negative `): fi: od: convert(gu,`+`)!/ convert([seq(gu[i]!,i=1..nops(gu))],`*`): end: #MishkalT1(Resh1): A quick way to compute #the (coeff. of the) second-to-leading contribution to Mishkal(Resh1,n) #i.e. Mishkal(Resh1,n,sum(Resh1[i],i=1..nops(Resh1))-1) #For example, try MishkalT0([2,2,2]); MishkalT1:=proc(Resh1) local i,j,gu,gu1: gu:=0: for i from 1 to nops(Resh1) do for j from i+1 to nops(Resh1) do gu1:=[op(1..i-1,Resh1),Resh1[i]-1,op(i+1..j-1,Resh1), Resh1[j]-1, op(j+1..nops(Resh1),Resh1),1]: gu:=gu+MTC(op(gu1)): od: od: gu: end: #Mishkal(L,r,n): Given a list of integers, L #of length k, say, finds #E[X_L[1]*...*X_L[k]] #where X_i is the number of monochromatic r-edge-colorings K_i #For example, try Mishkal([3,3,3],r,n); Mishkal1:=proc(mu,r,n) local gu,i,t,j,mish,mish1,katan,gadol: option remember: katan:=max(op(mu)): gadol:=convert(mu,`+`): mish1:= r^(nops(mu)-convert([seq(binomial(mu[i],2),i=1..nops(mu))],`+`)): mish:=mish1*MTC(op(mu))*binomial(n,gadol): for t from katan to gadol-1 do gu:=Tup(mu,{seq(i,i=1..t)}): mish1:=0: for j from 1 to nops(gu) do mish1:=mish1+Wt(gu[j],r): od: mish:=expand(mish+expand(binomial(n,t))*mish1): od: factor(expand(mish)): end: #Mishkal(L,r,n): Given a list of integers, L #of length k, say, finds #E[X_L[1]*...*X_L[k]] #where X_i is the number of monochromatic r-edge-colorings K_i #For example, try Mishkal([3,3,3],r,n); Mishkal:=proc(mu,r,n) local gu,i,t,j,mish,mish1,katan,gadol: option remember: katan:=max(op(mu)): gadol:=convert(mu,`+`): mish1:= r^(nops(mu)-convert([seq(binomial(mu[i],2),i=1..nops(mu))],`+`)): mish:=mish1*MTC(op(mu))*binomial(n,gadol): mish1:= r^(nops(mu)-convert([seq(binomial(mu[i],2),i=1..nops(mu))],`+`)): mish:=mish+mish1*MishkalT1(mu)*binomial(n,gadol-1): for t from katan to gadol-2 do gu:=Tup(mu,{seq(i,i=1..t)}): mish1:=0: for j from 1 to nops(gu) do mish1:=mish1+Wt(gu[j],r): od: mish:=expand(mish+expand(binomial(n,t))*mish1): od: factor(expand(mish)): end: #kthMomentV(a,k,r,n): The k^th moment of the r.v. product of #the r.v. #monochromatic K_{a} in edge-r-colorings of K_n #For example, try kthMomentV(3,3,r,n); kthMomentV:=proc(a,k,r,n) if k=0 then 1: else Mishkal([a$k],r,n): fi: end: #kthMoment(a,k,r,n): The k^th moment about the mean #of the r.v. product of #the r.v.s #monochromatic K_a in edge-r-colorings of K_n #For example, try kthMoment(3,3,r,n); kthMoment:=proc(a,k,r,n) local m,k1: m:=expand(binomial(n,a))/r^(binomial(a,2)-1): factor(expand(convert([seq((-1)^(k-k1)*binomial(k,k1)* m^(k-k1)*kthMomentV(a,k1,r,n),k1=0..k)],`+`))): end: ####Emprical Part########### #LabGraphs(n): all the labelled graphs on {1, .., n}. For example, #try: LabGraphs(4); LabGraphs:=proc(n) local i,j,gu1,gu,lu: option remember: if n=0 or n=1 then RETURN({{}}): fi: gu1:=LabGraphs(n-1): lu:={seq({i,n},i=1..n-1)}: lu:=powerset(lu): gu:={}: for i from 1 to nops(gu1) do for j from 1 to nops(lu) do gu:=gu union {gu1[i] union lu[j]}: od: od: gu: end: #IsClique(G,S): Given a graph G , and a subset of the #set of vertices, S, decides whether S is a clique IsClique:=proc(G,S) local i,j,gu: gu:={seq(seq({S[i],S[j]},j=i+1..nops(S)),i=1..nops(S))}: if gu minus G={} then true:else false: fi: end: #IsAClique(G,S): Given a graph G and a subset of the #set of vertices, S, decides whether S is an anti-clique IsAClique:=proc(G,S) local i,j,gu: gu:={seq(seq({S[i],S[j]},j=i+1..nops(S)),i=1..nops(S))}: if gu intersect G={} then true:else false: fi: end: #NumCliques(G,k,n): Given a graph G on {1, ... , n}, outputs #the number of cliques of size k NumCliques:=proc(G,k,n) local lu,i,co: lu:=choose({seq(i,i=1..n)},k): co:=0: for i from 1 to nops(lu) do if IsClique(G,lu[i]) then co:=co+1: fi: od: co: end: #NumCliques(G,k,n): Given a graph G on {1, ... , n}, outputs #the number of cliques of size k NumCliques:=proc(G,k,n) local lu,i,co: lu:=choose({seq(i,i=1..n)},k): co:=0: for i from 1 to nops(lu) do if IsClique(G,lu[i]) then co:=co+1: fi: od: co: end: #NumACliques(G,k,n): Given a graph G on {1, ... , n}, outputs #the number of anti-cliques of size k NumACliques:=proc(G,k,n) local lu,i,co: lu:=choose({seq(i,i=1..n)},k): co:=0: for i from 1 to nops(lu) do if IsAClique(G,lu[i]) then co:=co+1: fi: od: co: end: #kthMomentVE(a,k,n0) : If a and k , and n0, are positive intgers, #computes empirically the kth moment of the r.v. `number of mono-chromatic' #K_a over S_{n_0} in two-colorings of all the labelled graphs on n0 #vertices; For example, try: kthMomentVE(3,2,4); kthMomentVE:=proc(a,k,n0) local gu,i,co,G: gu:=LabGraphs(n0): co:=0: for i from 1 to nops(gu) do G:=gu[i]: co:=co+(NumACliques(G,a,n0)+NumCliques(G,a,n0))^k: od: co/nops(gu): end: #kthMomentE(a,k,n0) : If a and k , and n0, are positive intgers, #computes empirically the kth moment, (about the mean) # of the r.v. "number of mono-chromatic #K_a over S_{n_0} in two-colorings of all the labelled graphs on n0 #vertices" For example, try: kthMomentE(3,2,4); kthMomentE:=proc(a,k,n0) local gu,i,co,G,m: gu:=LabGraphs(n0): m:=binomial(n0,a)/2^(binomial(a,2)-1): co:=0: for i from 1 to nops(gu) do G:=gu[i]: co:=co+(NumACliques(G,a,n0)+NumCliques(G,a,n0)-m)^k: od: co/nops(gu): end: ####End Emprical Part FirstMomentV:=proc(k,r,n): binomial(n,k)/r^(binomial(k,2)-1): end: #SecondMomentV(k,r,n):The second moment of #number of mono-chromatic K_k #for r-colorings of the edges of K_n (k numeric, r,n symbolic) #Try SecondMomomentV(3,r,n); SecondMomentV:=proc(k,r,n) local gu,k1: gu:=0: for k1 from 0 to 1 do gu:=gu+expand(MTC(k1,k-k1,k-k1)/r^(2*binomial(k,2)-binomial(k1,2)-2)* binomial(n,2*k-k1)): od: for k1 from 2 to k do gu:=gu+expand(MTC(k1,k-k1,k-k1)/r^(2*binomial(k,2)-binomial(k1,2)-1)* binomial(n,2*k-k1)): od: factor(gu): end: #Variance(k,r,n):The second moment (about the mean) #of of the r.v. #number of mono-chromatic K_k #for r-colorings of the edges of K_n (k numeric, r,n symbolic) #Try Variance(3,r,n); Variance:=proc(k,r,n): factor(SecondMomentV(k,r,n)-(expand(binomial(n,k))/r^(binomial(k,2)-1))^2): end: #Skd(k,d): The set of vectors of integers [a12,a13,a23] that #are feasible, see the paper Skd:=proc(k,d) local a12,a13,a23,gu: gu:={}: for a12 from max(0,d-k) to d do for a13 from max(0,d-k) to d do for a23 from max(0,d-k) to d do if a12+a13<=d and a12+a23<=d and a13+a23<=d and a12+a13+a23-d>=0 then gu:=gu union {[a12,a13,a23]}: fi: od: od: od: gu: end: #WtVec3(vec): the weight attached to a feasible vector vec=[a12,a13,a23] #for the K_k problem with the union of S1, S2, S3 having 3k-d elements #For example, try WtVec3([0,1,0],3,1,r); WtVec3:=proc(vec,k,d,r) local NuE,NuC,lu,a12,a13,a23,a123,i: a12:=vec[1]:a13:=vec[2]:a23:=vec[3]: a123:=a12+a13+a23-d: NuE:=3*binomial(k,2)-binomial(a12,2)-binomial(a13,2)-binomial(a23,2) +binomial(a123,2): if a12<=1 and a13<=1 and a23<=1 then NuC:=3: elif (a12>=2 and a13<=1 and a23<=1) or (a12<=1 and a13>=2 and a23<=1) or (a12<=1 and a13<=1 and a23>=2) then NuC:=2: else NuC:=1: fi: NuE,NuC: lu:=[k+a12-d,k+a13-d,k+a23-d,d-a12-a23,d-a12-a13,d-a13-a23,a123]: if {seq(evalb(lu[i]>=0),i=1..nops(lu))}<>{true} then ERROR(`Bad input`): fi: MTC(op(lu))*r^(NuC-NuE): end: #Aadr the sum of the weights correspoding to all config. #[S1,S2,S3] with nops(S1 union S2 union S3)=3*a-d Aadr:=proc(a,d,r) local gu,i: gu:=Skd(a,d): convert([seq(WtVec3(gu[i],a,d,r),i=1..nops(gu))],`+`):end: #ThirdMomentV(a,r,n):The (straight) third moment of the random variable #Number of monochromatic K_a in r-edge colorings of K_n #(a an integer, r, and n symbolic) For example, try ThirdMomentV(3,r,n); ThirdMomentV:=proc(a,r,n) local d: factor(convert([seq(expand(binomial(n,3*a-d))*Aadr(a,d,r),d=0..2*a)],`+`)): end: #ThirdMoment(a,r,n):The third moment about the mean of the random variable #Number of monochromatic K_a in r-edge colorings of K_n #(a an integer, r, and n symbolic) For example, try ThirdMoment(3,r,n); ThirdMoment:=proc(a,r,n) local m: m:=expand(binomial(n,a))/r^(binomial(a,2)-1): factor(expand(ThirdMomentV(a,r,n)-3*m*SecondMomentV(a,r,n)+2*m^3)): end: #Bonf1(k,r,n): The first Bonferonni Expression for the number #of monochromatic K_k in the r-coloring of the edges of K_n Bonf1:=proc(k,r,n): expand(1-FirstMomentV(k,r,n)): end: #Bonf2(k,r,n): The second Bonferonni Expression for the number #of monochromatic K_k in the r-coloring of the edges of K_n Bonf2:=proc(k,r,n) local P,X: P:=expand(1-X+binomial(X,2)): expand(coeff(P,X,0)+coeff(P,X,1)*FirstMomentV(k,r,n) +coeff(P,X,2)*SecondMomentV(k,r,n)): end: #Bonf3(k,r,n): The third Bonferonni Expression for the number #of monochromatic K_k in the r-coloring of the edges of K_n Bonf3:=proc(k,r,n) local P,X,gu: P:=expand(1-X+binomial(X,2)-binomial(X,3)): gu:= coeff(P,X,0)+coeff(P,X,1)*FirstMomentV(k,r,n) +coeff(P,X,2)*SecondMomentV(k,r,n)+coeff(P,X,3)*ThirdMomentV(k,r,n): gu:=expand(expand(gu)): end: #LastPos1(P1,n): Given an expression P1, in n, finds the last n0 #such that P1(n0)>0 LastPos1Old:=proc(P1,n) local n0: for n0 from 1 do if subs(n=n0,P1)<=0 then RETURN(n0-1): fi: od: end: #LastPos(P,n,r,L): Given an expression P and n and r #and a positive integer L, returns the sequence a[i] (i=L1..L2) #where n=a[i] is the last integer for which P(i,n) is positive LastPos:=proc(P,n,r,L) local r0 : [seq(LastPos1(subs(r=r0,P),n),r0=1..L)]: end: #LastPos1(P,n): Given an expression P, in n, finds the last n0 #such that P(n0)>0 LastPos1:=proc(P,n) local n0,P1,L1,L2,L3,L4: Digits:=50: P1:=evalf(P): for n0 from 0 while subs(n=n0*1000,P1)>0 do od: L4:=n0-1: for n0 from 0 while subs(n=L4*1000+n0*100,P1)>0 do od: L3:=n0-1: for n0 from 0 while subs(n=L4*1000+L3*100+10*n0,P1)>0 do od: L2:=n0-1: for n0 from 0 while subs(n=L4*1000+L3*100+10*L2+n0,P1)>0 do od: L1:=n0-1: L4*1000+L3*100+L2*10+L1: end: #FkL(k,L) : The performance of Bonf1 vs. Bonf3 for #colors 2 through L #with monochromatic K_k (k numeric) FkL:=proc(k,L) local n,r0,r,lu1,lu3: lu1:=Bonf1(k,r,n): lu3:=Bonf3(k,r,n): seq([LastPos1(subs(r=r0,lu1), n),LastPos1(subs(r=r0,lu3), n)],r0=2..L): end: #GrL(r,L) : The performance of Bonf1 vs. Bonf3 for r colors #with monochromatic K_k (k 3 through L GrL:=proc(r,L) local n,k0: seq([LastPos1(Bonf1(k0,r,n), n),LastPos1(Bonf3(k0,r,n), n)],k0=3..L): end: #ThirdMomentTV(a,r,n,T):The first T leading terms (in n) #of ThirdMomentV(a,r,n): #(r, and n are symbolic, a numeric, T is numeric) #For example, try ThirdMomentTV(3,r,n,3); ThirdMomentTV:=proc(a,r,n,T) local d,deg,gu,i: gu:=convert([seq(expand(binomial(n,3*a-d))*Aadr(a,d,r),d=0..T)],`+`): deg:=degree(gu,n): gu:=convert([seq(coeff(gu,n,deg-i)*n^(deg-i),i=0..T-1)],`+`): end: #ThirdMomentT(a,r,n,T) The first T leading terms (in n) #of ThirdMomentT(a,r,n) #( r, and n symbolic, a numeric) For example, try ThirdMomentT(3,r,n,1); ThirdMomentT:=proc(a,r,n,T) local m,i,deg,gu: m:=expand(expand(binomial(n,a))/r^(binomial(a,2)-1)): gu:=expand(expand(ThirdMomentTV(a,r,n,T+5)-3*m*SecondMomentV(a,r,n)+2*m^3)): deg:=degree(gu,n): convert([seq(coeff(gu,n,deg-i)*n^(deg-i),i=0..T-1)],`+`): end: #erk(r,k), if r is a positive integer and k is a symbolic integer #computes the polynomial in k: e_r(1,..., k), where #e_r is the r^th elementary symmetric function in k variables erk:=proc(r,k) local gu,K: if r=0 then RETURN(1): fi: gu:=k*subs(k=k-1,erk(r-1,k)): gu:=sum(gu,k=r-1..K): expand(subs(K=k,gu)): end: #binomialT1(n,k,T): The first T leading terms of binomial(n,k) #(k symbolic!). For example, try binomialT(n,k,1); (it should be n^k/k!) binomialT1:=proc(n,k,T) local r: convert([seq(n^(k-r)*(-1)^r*expand(subs(k=k-1,erk(r,k)))/k!,r=0..T-1)],`+`): end: #binomialT(n,k,T): The first T leading terms of binomial(n,k) #(k symbolic!). For example, try binomialT(n,k,1); (it should be n^k/k!) binomialT:=proc(n,k,T) local k1:subs(k1=k,binomialT1(n,k1,T)):end: #ThirdMomentTVS(a,r,n,T):The first T leading terms (in n) #of ThirdMomentV(a,r,n): (as an expression in a, r, and n) #(a, r, and n are symbolic, T is numeric) #For example, try ThirdMomentTV(3,r,n,3); ThirdMomentTVS:=proc(a,r,n,T) local d,deg,gu,i: gu:=convert([seq(expand(binomialT(n,3*a-d,T))*Aadr(a,d,r),d=0..T)],`+`): deg:=degree(gu,n): gu:=convert([seq(coeff(gu,n,deg-i)*n^(deg-i),i=0..T-1)],`+`): end: #SdS(d): The set of vectors of integers [a12,a13,a23] that #are feasible for symbolic k), see the paper SdS:=proc(d) local a12,a13,a23,gu: gu:={}: for a12 from 0 to d do for a13 from 0 to d do for a23 from 0 to d do if a12+a13<=d and a12+a23<=d and a13+a23<=d and a12+a13+a23-d>=0 then gu:=gu union {[a12,a13,a23]}: fi: od: od: od: gu: end: #AadrS the sum of the weights correspoding to all config. #[S1,S2,S3] with nops(S1 union S2 union S3)=3*a-d (symbolic a) AadrS:=proc(a,d,r) local gu,i: gu:=SdS(d): convert([seq(WtVec3S(gu[i],a,d,r),i=1..nops(gu))],`+`):end: #WtVec3S(vec,k,d,r): the weight attached to a feasible vector vec=[a12,a13,a23] #for the K_k problem with the union of S1, S2, S3 having 3k-d elements #symbolic k #For example, try WtVec3([0,1,0],3,1,r); WtVec3S:=proc(vec,k,d,r) local NuE,NuC,lu,a12,a13,a23,a123: a12:=vec[1]:a13:=vec[2]:a23:=vec[3]: a123:=a12+a13+a23-d: NuE:=expand(3*binomial(k,2)-binomial(a12,2)-binomial(a13,2)-binomial(a23,2) +binomial(a123,2)): if a12<=1 and a13<=1 and a23<=1 then NuC:=3: elif (a12>=2 and a13<=1 and a23<=1) or (a12<=1 and a13>=2 and a23<=1) or (a12<=1 and a13<=1 and a23>=2) then NuC:=2: else NuC:=1: fi: NuE,NuC: lu:=[k+a12-d,k+a13-d,k+a23-d,d-a12-a23,d-a12-a13,d-a13-a23,a123]: MTCS(op(lu))*r^(NuC-NuE): end: MTCS:=proc() local gu,i: gu:=[args]: for i from 1 to nops(gu) do if type(gu[i],integer) and gu[i]<0 then RETURN(0): fi: od: convert(gu,`+`)!/ convert([seq(gu[i]!,i=1..nops(gu))],`*`): end: #ThirdMomentTVS(a,r,n,T):The first T leading terms (in n) #of ThirdMomentV(a,r,n): (as an expression in a, r, and n) #(a, r, and n are symbolic, T is numeric) #For example, try ThirdMomentTV(3,r,n,3); ThirdMomentTVS:=proc(a,r,n,T) local d,gu,i: print(`No Good`): gu:=convert([seq(expand(binomialT(n,3*a-d,T+3))*AadrS(a,d,r),d=0..T)],`+`): gu:=expand(simplify(a!^3*gu/n^(3*a))): gu:=convert([seq(coeff(gu,n,-i)*n^(-i),i=0..T-1)],`+`): gu:=simplify(expand(gu)): n^(3*a)*gu/a!^3: end: #ThirdMomentTS(a,r,n,T) The first T leading terms (in n) #of ThirdMomentT(a,r,n) #( a, r, and n symbolic) For example, try ThirdMomentTS(k,r,n,1); ThirdMomentTS:=proc(a,r,n,T) local m,i,deg,gu: m:=expand(expand(binomialT(n,a,T+3))/r^(binomial(a,2)-1))/n^a: m:=expand(m): #HERE gu:=expand(expand(ThirdMomentTV(a,r,n,T+5)-3*m*SecondMomentV(a,r,n)+2*m^3)): deg:=degree(gu,n): convert([seq(coeff(gu,n,deg-i)*n^(deg-i),i=0..T-1)],`+`): end: #SecondMomentTV(k,r,n,T):The T leading tens #of second moment of #number of mono-chromatic K_k #for r-colorings of the edges of K_n (k numeric, r,n symbolic) #Try SecondMomomentTV(3,r,n); SecondMomentTV:=proc(k,r,n) local gu,k1: gu:=0: for k1 from 0 to 1 do gu:=gu+expand(MTC(k1,k-k1,k-k1)/r^(2*binomial(k,2)-binomial(k1,2)-2)* binomial(n,2*k-k1)): od: for k1 from 2 to k do gu:=gu+expand(MTC(k1,k-k1,k-k1)/r^(2*binomial(k,2)-binomial(k1,2)-1)* binomial(n,2*k-k1)): od: factor(gu): end: NorThirdMomentV:=proc(a,r,N) local gu,lu : gu:=ThirdMomentV(a,r,N): lu:=N^(a*3)/a!^3/r^(3*binomial(a,2)-3): expand(gu/lu): end: NorThirdMoment:=proc(a,r,N) local gu,lu : gu:=ThirdMoment(a,r,N): lu:=N^(a*3)/a!^3/r^(3*binomial(a,2)-3): expand(gu/lu): end: NorVariance:=proc(k,r,n): expand(Variance(k,r,n)/(N^k/k!/r^(binomial(k,2)-1))^2): end: #CheckkthMoment(a,k,n0): #Checks the procedure kthMoment by doing it empirically #for 2-colors, K_n0, kth-monent. number of mon-chromatic a CheckkthMoment:=proc(a,k,n0): evalb(kthMomentE(a,k,n0)= subs({r=2,n=n0},kthMoment(a,k,r,n))): end: ####The Symbolic Part (in k)#### ez:=proc():print(`WComp(L,s), Fnks(n,k,s), Merge1(S1,S2), Fns(n,s) `): print(`AlphVec(Vec1,n), WtVec1(Vec1,n,k)`): print(`WtVec(Vec1,n,k,s,N,T), kthMomentSV(a,k,r,N,T), ConVec1(Vec1,n)`): print(`EdgeDef(Vec1,n), CCIG(Vec1,n) , kthMomentS(a,k,r,N,T)`): print(`Kurtosis(a,r,N,T)`): end: with(combinat): #WComp(L,s): given a list L, of length n, say, #and an integer s , finds the #set of vectors of non-negative integers [a1, ..., an] #such that a[1]*L[1]+ ... + a[n]*L[n]=s WComp:=proc(L,s) local gu,an,n,L1,s1,gu1,i: n:=nops(L): if s<0 then RETURN({}): fi: if n=0 then if s=0 then RETURN({[]}): else RETURN({}): fi: fi: gu:={}: for an from 0 to trunc(s/L[n]) do L1:=[op(1..n-1,L)]: s1:=s-an*L[n]: gu1:=WComp(L1,s1): gu:=gu union {seq([op(gu1[i]),an],i=1..nops(gu1))}: od: gu: end: #Fnks(n,k,s): all the functions f: choose(n,k)->N that sum up to s Fnks:=proc(n,k,s) local gu,L1,lu,i,j: gu:=choose(n,k): L1:=[1$nops(gu)]: lu:=WComp(L1,s): {seq({seq([gu[j],lu[i][j]],j=1..nops(gu))},i=1..nops(lu))}: end: #Merge1(S1,S2): Given two sets S1 and S2, outputs the #set {S1[i] union S2[j]} Merge1:=proc(S1,S2) local i,j,gu: gu:={}: for i from 1 to nops(S1) do for j from 1 to nops(S2) do gu:=gu union {S1[i] union S2[j]}: od: od: gu: end: #Fns(n,s) #All eligible funcions a[S] (|S|>=2), S subsets of {1,..n} that add up to #s Fns:=proc(n,s) local mu,L,k,gu,mu1,gu1,i: L:=[seq(k-1,k=2..n)]: mu:=WComp(L,s): gu:={}: for i from 1 to nops(mu) do mu1:=mu[i]: gu1:=Fnks(n,2,mu1[1]): for k from 3 to n do gu1:=Merge1(gu1,Fnks(n,k,mu1[k-1])): od: gu:=gu union gu1: od: gu: end: #AlphVec(Vec1,n): Given a function f on integer-valued functions #on powerset({1, ..., n}) minus {{},seq({i},i=1..n)} #returns the vector a[i]=sum(f(a[i union S]),S) where #S ranges over all non-empty subsets of {1, ...,i-1,i+1, ..n}. AlphVec:=proc(Vec1,n) local T,i,gu,i1,lu,su,ku: if nops(Vec1)<>2^n-n-1 then ERROR(`The size of the first argument should be`, 2^n-n-1 ): fi: for i from 1 to nops(Vec1) do T[Vec1[i][1]]:=Vec1[i][2]: od: gu:=[]: for i from 1 to n do lu:=powerset({seq(i1,i1=1..i-1),seq(i1,i1=i+1..n)}) minus {{}}: su:=0: for i1 from 1 to nops(lu) do ku:=lu[i1] union {i}: ku:=sort(convert(ku,list)): su:=su+T[ku]: od: gu:=[op(gu),su]: od: gu: end: #WtVec1(Vec1,n,k): The weight of a function Vec1 in one of the Fns(n,s) #(k is the symbol denoting the K_k in the monochromaic K_k) #For example, WtVec1(Fns(3,2)[1],3,k); WtVec1:=proc(Vec1,n,k) local lu,i: lu:=AlphVec(Vec1,n): convert([seq(expand(simplify(k!/(k-lu[i])!)),i=1..nops(lu))],`*`)* convert([seq(1/(Vec1[i][2])!,i=1..nops(Vec1))],`*`): end: #erk(r,k), if r is a positive integer and k is a symbolic integer #computes the polynomial in k: e_r(1,..., k), where #e_r is the r^th elementary symmetric function in k variables erk:=proc(r,k) local gu,K: if r=0 then RETURN(1): fi: gu:=k*subs(k=k-1,erk(r-1,k)): gu:=sum(gu,k=r-1..K): expand(subs(K=k,gu)): end: #FF1(n,k,T): The first T+1 leading terms of n(n-1) ... (n-k+1)/n^k #(k symbolic!). For example, try FF1(n,k,0); (it should be 1) FF1:=proc(n,k,T) local r: convert([seq(n^(-r)*(-1)^r*expand(subs(k=k-1,erk(r,k))),r=0..T)],`+`): end: #FF(n,k,T): The first T leading terms of n(n-1) ... (n-k+1)/n^k #(k symbolic!). For example, try FF(n,k,1); (it should be 1) FF:=proc(n,k,T) local k1:subs(k1=k,FF1(n,k1,T)):end: #WtVec(Vec1,n,k,s,N,T): The first T terms of the vector Vec1 #that came from Fns(n,s) WtVec:=proc(Vec1,n,k,s,N,T) local gu: gu:=WtVec1(Vec1,n,k): gu:=expand((gu/N^s)*FF(N,n*k-s,T-s)): end: #kthMomentSV(a,k,r,N,T): The explicit expression for the #first T+1 leading terms for the (straight) n-th moment #of the r.v. #monochmromatic K_a in edge-coloring of #K_N with r colorings as a (symbolic) expression in #N and a. N , a and r are symbols, k and T are (numeric) integers #For example, try kthMomentSV(a,2,r,N,2); kthMomentSV:=proc(a,k,r,N,T) local gu,su,i,s,ku: if k=0 then RETURN(1): fi: if k=1 then RETURN(FF(N,a,T)): fi: su:=0: for s from 0 to T do gu:=Fns(k,s): for i from 1 to nops(gu) do ku:=WtVec(gu[i],k,a,s,N,T)*WtVec2(gu[i],k,r): su:=expand(su+ku): od: od: su: end: #ConVec1(Vec1,n): Given a function from sets to non-negative integers #(a typical member of Fns(n,s) for some s) outpus #the function g(T)=sum(f(S)) over all S supersets of T ConVec1:=proc(Vec1,n) local T,gu,i,GU,i1,TU,su,mu,tu: gu:={}: for i from 1 to nops(Vec1) do gu:=gu union {Vec1[i][1]}: T[Vec1[i][1]]:=Vec1[i][2]: od: GU:={}: for i from 1 to nops(gu) do mu:=gu[i]: tu:={seq(i1,i1=1..n)} minus convert(mu,set): TU:=powerset(tu): su:=0: for i1 from 1 to nops(TU) do su:=su+T[sort(convert(TU[i1] union convert(mu,set),list))]: od: GU:=GU union {[mu,su]}: od: GU: end: #The number of edges missing from n binomial(k,2) #in a config. [S1, ... , Sn] from the profile Vec1 #For example EdgeDef(Fns(3,3)[1],3); EdgeDef:=proc(Vec1,n) local gu,i,su: gu:=ConVec1(Vec1,n): su:=0: for i from 1 to nops(gu) do su:=su+(-1)^nops(gu[i][1])*binomial(gu[i][2],2): od: su: end: #CCIG(Vec1,n): Number of connected componets in the intersection #graph ####Begin onnected components of a graph #Neighs(G,V,SetV): the collective neigbors of the vertices in SetV Neighs:=proc(G,V,SetV) local i1: {seq(op(Neighs1(G,V,SetV[i1])),i1=1..nops(SetV))}: end: #ConComp1(G,V,v): The connected component of vertex v ConComp1:=proc(G,V,v) local gu,gu1: gu:={v}: gu1:=Neighs(G,V,gu): while gu1<>gu do gu:=gu1: gu1:=Neighs(G,V,gu1): od: gu1: end: #Neighs1(G,v): all the neighbors of the vertex v in G (including v) #where V is the set of edges Neighs1:=proc(G,V,v) local gu,i1,V1: gu:={v}: V1:= V minus {v}: for i1 from 1 to nops(V1) do if member({V1[i1],v},G) then gu:=gu union {V1[i1]}: fi: od: gu: end: #NuCC(G,V): the connected components of the graph G,V NuCC:=proc(G,V) local i1: nops({seq(ConComp1(G,V,V[i1]),i1=1..nops(V))}): end: ####Ends Connected components of a graph CCIG:=proc(Vec1,n) local i,G,V,gu: V:={seq(i,i=1..n)}: gu:=ConVec1(Vec1,n): G:={}: for i from 1 to nops(gu) do if nops(gu[i][1])=2 and gu[i][2]>=2 then G:=G union {convert(gu[i][1],set)}: fi: od: NuCC(G,V): end: WtVec2:=proc(Vec1,n,r): r^(EdgeDef(Vec1,n)+CCIG(Vec1,n)-n): end: kthMomentS:=proc(a,k,r,N,T) local k1,gu,T1,m,i: T1:=T+k: m:=kthMomentSV(a,1,r,N,T1): gu:=0: for k1 from 0 to k do gu:=gu+binomial(k,k1)*kthMomentSV(a,k1,r,N,T1)*(-m)^(k-k1): gu:=expand(gu): gu:=convert([seq(coeff(gu,N,-i)/N^i,i=0..T1)],`+`): od: gu:=expand(gu): sort( convert([seq(factor(coeff(gu,N,i))*N^i,i=ldegree(gu,N)..degree(gu,N))],`+`) ): end: #Kurtosis(k,r,N,T): The first T terms in the asymptotic expression #in N (with k fixed but symbolic) #for the Kurtosis of the r.v. #of mono-chrmomatic K_k in #r-colorings of the edges of K_N (r also symbolic) Kurtosis:=proc(a,r,N,T) local gu4,gu2,x,gu,i: gu2:=kthMomentS(a,2,r,N,T+2): gu4:=kthMomentS(a,4,r,N,T+2): gu:=normal(gu4/gu2^2): gu:=normal(subs(N=1/x,gu)): gu:=taylor(gu,x=0,T+2): convert([seq(factor(coeff(gu,x,i))/N^i,i=0..T)],`+`): end: #kthMomentBV(n,p,k): the straight k-th moment of tbe Binomial distribution #Binomial(n,p) kthMomentBV:=proc(n,p,k) local gu,i,t: gu:=(p*t+1-p)^n: for i from 1 to k do gu:=t*diff(gu,t): od: factor(subs(t=1,gu)): end: kthMomentB:=proc(n,p,k) local k1,gu,m: m:=kthMomentBV(n,p,1): gu:=0: for k1 from 0 to k do gu:=gu+binomial(k,k1)*kthMomentBV(n,p,k1)*(-m)^(k-k1): gu:=expand(gu): od: gu: end: #CheckkthMomentS(a0,k,r,N,T): checks procedure kthMomentS #by comparing it with numerical instances (a=a0) of kthMoment(a,k,r,N): #For example, try: CheckkthMomentS(4,2,r,N,3); CheckkthMomentS:=proc(a0,k,r,N,T) local lu,mu,a,i,mu1: lu:=subs(a=a0,kthMomentS(a,k,r,N,T)): mu:=kthMoment(a0,k,r,N): mu:=mu/(N^a0/a0!/r^(binomial(a0,2)-1))^k: mu1:=convert([seq(coeff(mu,N,i)*N^i,i=degree(mu,N)-T+1..degree(mu,N))],`+`): lu,mu1,expand(lu-mu1): end: ####End of the Symbolic Part (in k)