###################################################################### ##HANS: Save this file as HANS # ## To use it, stay in the # ##same directory, get into Maple (by typing: maple ) # ##and then type: read REK # ##Then follow the instructions given there # # # ##Written by Andrew V. Sills, Georgia Southern University, and # ##Doron Zeilberger, Rutgers University , # #Please report bugs to: # #zeilberg at math dot rutgers dot edu # ###################################################################### #Created: Aug-Oct., 2011 #First posted, Sept. 29, 2011 # Version of October 19, 2011. with(numtheory): print(`Created: Aug.-Oct., 2011`): print(` This is HANS `): print(`It accompanies the article `): print(`"Rademacher's Infinite Partial Fraction Conjecture is (almost certainly) False"`): print(`by Andrew V. Sills and Doron Zeilberger`): print(`and also available from Sills' and Zeilberger's websites.`): print(``): print(`Please report bugs to zeilberg at math dot rutgers dot edu`): print(``): print(`The most current version of this package and paper`): print(` are available from`): print(`http://www.math.rutgers.edu/~zeilberg/ .`): print(`For a list of the procedures type ezra();, for help with`): print(`a specific procedure, type ezra(procedure_name); .`): print(``): with(combinat): ezra1:=proc() if args=NULL then print(` The supporting procedures are: Ana1, AsyB`): print(` C01d, C1hklN, ChkFormulaList, Chkgf, L32, pmngf, RecEqn`): print(` SeriesExpansionCoeffs, shk, whk, Zinn `): else ezra(args): fi: end: ezra:=proc() if args=NULL then print(`The main procedures are: Ana, ChkFormula, ChklInfinity, ChklN, `): print(` ChklNdirect, ChklNfloat,ChklSeq, ChklSeqF, CloseEncounters `): print(`ClosestN, CNi, CNid, `): print(` C01, C01f, C01Seq, C01Seqf, C01Table, C01TableF, C01TableMinus `): print(` HansTopDownAutoPaper, PaV, PNz, SipurAsy01, `): elif nops([args])=1 and op(1,[args])=Ana then print(`Ana(L,n): inputs a sequence of real numbers L, and outputs`): print(`conjenctured asymptotics in terms of a list of length`): print(`P, where P is the "period" and the ith-term of the`): print(`output list gives an estimate for the asymptotic behavior of`): print(`L[P*n+i], for i from 1 to P-1, and the P-th term an`): print(`for L[P*n]. `): print(`It also returns the list [max,min] where`): print(`max is the integer mod P where the local maxima occur`): print(`and min is the integer mod P where the local minima occur`): print(`For example, try:`): print(`read HANSC01: Ana(evalf(C01r[1]),n);`): elif nops([args])=1 and op(1,[args])=Ana1 then print(`Ana1(L): analyzes the pattern of differences between peaks and valleys`): print(`of the sequence L, for example, once you have the data file HANS800, and uploaded it by typing`): print(`read HANS800: `): print(` try `): print(`Ana1(HANS800[1]);`): elif nops([args])=1 and op(1,[args])=AsyB then print(`AsyB:=proc(L): Given a list of positive numbers, tries`): print(`to find an asympototic expression for them as C*mu^n*n^theta`): print(`it returns the approx. conjenctured asymptotics `): print(`in the form C*mu^n*n^theta. `): print(`AsyB([seq(5*2^i*i^2,i=1..10)],n);`): elif nops([args])=1 and op(1,[args])=C01Table then print(`C01Table(K,L): A table of C_{01l}(N) for l from 1 to L`): print(`and N from 1 to K, it is given as a list of lists`): print(`such that T[l][N] is C[0,1,l](N). For example, try:`): print(`C01Table(50,10);`): elif nops([args])=1 and op(1,[args])=C01TableMinus then print(`C01TableMinus(K,L): A table of C_{01l}(N) for l from 0 to -L`): print(`and N from 1 to K, it is given as a list of lists`): print(`such that T[l+1][N] is C[0,1,-l](N). For example, try:`): print(`C01TableMinus(50,10);`): elif nops([args])=1 and op(1,[args])=C01TableF then print(`C01TableF(K,L,Dig): A table, in floating-point`): print(`with Dig digits`): print(`of C_{01l}(N) for l from 1 to L`): print(`and N from 1 to K, it is given as a list of lists`): print(`such that T[l][N] is C[0,1,l](N). For example, try:`): print(`C01TableF(50,10,40);`): elif nops([args])=1 and op(1,[args])=ChklN then print(`ChklN(h,k,l,N): the coeff. of 1/(x-exp(2*Pi*h/k))^l in `): print(`1/(1-x)/(1-x^2)/.../(1-x^N). computed by the recurrence`): print(`For example, try:`): print(`ChklN(0,1,1,5);`): elif nops([args])=1 and op(1,[args])=ChkFormula then print(`ChkFormula(h,k,r,n,R): returns a formula, in symbolic n for`): print(`the Rademacher coefficient C[h,k,n-R](k*n+r) ` ): print(`For example, try:`): print(`ChkFormula(1,2,1,n,2)`): elif nops([args])=1 and op(1,[args])=ChkFormulaList then print(`ChkFormulaList(h,k,r,n,R): returns a list of formulas for`): print(`the Rademacher coefficient C[h,k,n-j](k*n+r) for j from 0 to R-1`): print(`For example try:`): print(`ChkFormulaList(1,2,1,n,3)`): elif nops([args])=1 and op(1,[args])=Chkgf then print(`Chkgf(N,x,h,k): returns the Taylor series, up to`); print(` O( (x-exp(2*Pi*I*h/k)^N ) of the function`): print(` x^(floor(N/k))/( (1-x)*(1-x^2)* . . . * (1-x^N) )`): print(` the first floor(N/k) coefficients of which are `): print(` Rademacher's coefficients C[h,k,n-j](k*n+r), j=0..n-1 `): print(` For example, try:`): print(`Chkgf(6,x,1,4)`): elif nops([args])=1 and op(1,[args])=ChklInfinity then print(`ChklInfinity(h,k,l): Rademacher's definition of C[h,k,l](infinity)`): print(`Eq. (130.6) in p. 302`): print(`of his book "Topics in Analytic Number Theory"`): print(`(Springer, 1973, Grundlehren series Band 169)`): print(`For example, try:`): print(`ChklInfinity(0,1,1);`): elif nops([args])=1 and op(1,[args])=ChklNfloat then print(`ChklNfloat(h,k,l,N,Dig): the coeff. of 1/(x-exp(2*Pi*h/k))^l in `): print(`1/(1-x)/(1-x^2)/.../(1-x^N). computed by the recurrence`): print(`but using floating point with Dig digits.`): print(`For example, try:`): print(`ChklNfloat(0,1,1,5,30);`): elif nops([args])=1 and op(1,[args])=C1hklN then print(`C1hklN(h,k,l,N): the coeff. of 1/(x-exp(2*Pi*h/k))^(trunc(N/k)-l) in`): print(`1/(1-x)/(1-x^2)/.../(1-x^N). computed by the recurrence`): print(`For example, try:`): print(`C1hklN(0,1,1,5);`): elif nops([args])=1 and op(1,[args])=ChklNdirect then print(`ChklNdirect(h,k,l,N): the coeff. of 1/(x-exp(2*Pi*h/k))^l in`): print(`1/(1-x)/(1-x^2)/.../(1-x^N). computed directly`): print(`For example, try:`): print(`ChklNdirect(0,1,1,5);`): elif nops([args])=1 and op(1,[args])=ChklSeq then print(`ChklSeq(h,k,l,N0):The first N0 terms of the Rademacher sequence`): print(`C_{hkl}(0,1,1,60). Try: `): print(`ChklSeq(0,1,1,100);`): elif nops([args])=1 and op(1,[args])=ChklSeqF then print(`ChklSeqF(h,k,l,N0,Dig):The first N0 terms of the Rademacher sequence`): print(`C_{hkl}(N). In floating-point with Dig Digits. Try: `): print(`ChklSeqF(0,1,1,100,20);`): elif nops([args])=1 and op(1,[args])=CloseEncounters then print(`CloseEncounters(K1,L1,N0): The list of the closest`): print(`approach to Rademacher's C_{hkl}(infinity)`): print(`together with the deviations for all 0<=h<=k<=K1`): print(` 1<=l<=L1 for N<=N0. For example try:` ): print(`CloseEncounters(3,4,50);`): elif nops([args])=1 and op(1,[args])=ClosestN then print(`ClosestN(h,k,l,K): The integer N between 1 and K that comes`): print(`closest to ChklInfinity(h,k,l), together with the difference`): print(`For example, try:`): print(`ClosestN(0,1,1,100);`): elif nops([args])=1 and op(1,[args])=CNi then print(`CNi(N,i): the coeff. of z^i in z^(N-1)/((z+1)^2-1))/((z+1)^3-1)...`): print(`by the recurrence. For example, try:`): print(`CNi(5,4);`): elif nops([args])=1 and op(1,[args])=CNid then print(`CNid(N,i): the coeff. of z^i in z^(N-1)/((z+1)^2-1))/((z+1)^3-1)...`): print(`computed directly. For example, try:`): print(`CNid(5,4);`): elif nops([args])=1 and op(1,[args])=C01 then print(`C01(j,N): The Rademacher coeff. C_{01j}(N)`): print(`computed by the efficient recurrence. For example, try:`): print(`C01(1,8);`): elif nops([args])=1 and op(1,[args])=C01d then print(`C01d(j,N): The Rademacher coeff. C_{01j}(N)`): print(`computed directly. For example, try:`): print(`C01(1,8);`): elif nops([args])=1 and op(1,[args])=C01f then print(`C01f(j,N,Dig): The Rademacher coeff. C_{01j}(N) in FLOATING point`): print(`using Dig Digits`): print(`could give wrong (and often nonsense) answers if Dig is not`): print(`large enough. Try:`): print(`C01f(1,8,20);`): elif nops([args])=1 and op(1,[args])=C01Seq then print(`C01Seq(j,N0):The first N0 terms of the Rademacher sequence`): print(`C_{01j}(N). Try: `): print(`C01Seq(1,100);`): elif nops([args])=1 and op(1,[args])=C01Seqf then print(`C01Seqf(j,N0):The first N0 terms of the Rademacher sequence`): print(`C_{01j}(N), in floating-point. Try: `): print(`C01Seqf(1,100);`): elif nops([args])=1 and op(1,[args])=HansTopDownAutoPaper then print(`HansTopDownAutoPaper(h,k,r,n,C,R): `): print(`automatically generates a paper which states and proves`): print(`formulas for Rademacher coefficients C[h,k,n-j](k*n+r) `): print(`for j from 0 to R-1`): print(`For example, try:`): print(`HansTopDownAutoPaper(3,4,0,n,C,3): `): elif nops([args])=1 and op(1,[args])=L32 then print(`L32(y1): L_{3/2} in Rademacher's formula on p.302`): print(`of his book "Topics in Analytic Number Theory"`): print(`(Springer, 1973, Grundlehren series Band 169)`): print(`For example, try`): print(`L32(-4);`): elif nops([args])=1 and op(1,[args])=PaV then print(`PaV(L) the locations and values of the local max and min of `): print(`a sequence. For example, try:`): print(`PaV([3,2,1,5,7,9,5,4]);`): elif nops([args])=1 and op(1,[args])=pmngf then print(`pmngf(m,x): the generating function for the number of partitions`): print(`of n into at most m parts, factored into cyclotomic polynomials`): print(`For example, try`): print(`pmngf(4,x)`): elif nops([args])=1 and op(1,[args])=PNz then print(`PNz(N,z): The rational generating function`): print(` Sum(CNi(N,i)*z^i,i=0..infinity)`): print(`For example, try:`): print(`PNz(2,z);`): elif nops([args])=1 and op(1,[args])=QuasiZinn then print(`QuasiZinn(L,p): given a list of numbers L, and a pos. integer p`): print(`tries to guess, for each n mod p, the appx. asymptotics`): print(`in the form C*mu^n*n^theta, the output is a list of triples`): print(`[C,theta,mu] for each residue class`): print(`For example, try:`): print(`read HANS500: QuasiZinn(lu[1],32);`): elif nops([args])=1 and op(1,[args])=RecEqn then print(`ReqEqn(h,k,r,n,j,C): returns the recurrence equation that results`): print(`from comparing the coefficients of`, (x-exp(2*Pi*I*h/k))^j ): print(`on either side of the generating function recurrence equation`): print(`(1-x^(k*n+r))*(1-x^(k*n+r-1))*...*(1-x^(k*n+r-k+1))* Chkgf(k*n+r,x,h,k) = (x-exp(2*Pi*I*h/k)) * Chkgf(k*n+r-k,x,h,k)`): elif nops([args])=1 and op(1,[args])=SeriesExpansionCoeffs then print(`SeriesExpansionCoeffs(h,k,r,n,s): `): print(`the coefficient of (x- exp(2*Pi*I*h/k))^s in the Taylor expansion`): print(`of (1-x^(k*n+r))*(1-x^(k*n+r-1))*...* (1-x^(k*n+r-k+1))`): print(`about x = exp(2*Pi*I*h/k) `): print(`For example, try`): print(`SeriesExpansionCoeffs(1,2,0,n,3)`): elif nops([args])=1 and op(1,[args])=SipurAsy01 then print(`SipurAsy01(n,N): The guessed asymptotics of C01l(32n+r) for l from 1 to 15 and r from 0 to 31.`): print(`Type:`): print(`SipurAsy01(n,N):`): elif nops([args])=1 and op(1,[args])=shk then print(`shk(h,k): The Dedekind sum Eq. (5.2.6) in George Andrews's`): print(`The Theory of Partitions" (Chap. 5.2, p. 72)`): print(`For example, try:`): print(`shk(11,5);`): elif nops([args])=1 and op(1,[args])=whk then print(`whk(h,k): The 24-th root of unity of Rademacher`): elif nops([args])=1 and op(1,[args])=Zinn then print(`Zinn(resh): Zinn-Justin's method to estimate`): print(`the mu, and theta such that`): print(`resh[i] is appx. Const*mu^i*i^theta. `): print(`For example, try:`): print(`Zinn([seq(5*i*2^i,i=1..30)]);`): else print(`There is no ezra for`,args): fi: end: ez:=proc(): print(`C01(j,N), C01d(j,N), CNid(N,i), CNi(N,i) `): end: #C01d(j,N): The Rademacher coeff. C_{01j}(N) #computed directly C01d:=proc(j,N) local gu,y,i: gu:=normal((y-1)^N/mul(1-y^i,i=1..N)): gu:=subs(y=1+y,gu): gu:=taylor(gu,y=0,N-j+1): coeff(gu,y,N-j): end: #CNid(N,i): the coeff. of z^i in 1/*(1+(1+z))*(1+(1+z)+(1+z+z^2)*...)) #computed directly CNid:=proc(N,j) local gu,y,i: gu:=normal((y-1)^N/mul(1-y^i,i=1..N)): gu:=subs(y=1+y,gu): gu:=taylor(gu,y=0,j+1): coeff(gu,y,j): end: #CNi(N,i): the coeff. of z^i in z^(N-1)/((z+1)^2-1))/((z+1)^3-1)... #by the recurrence. For example, try: #CNi(5,4); CNi:=proc(N,i) local r: option remember: if i<0 then RETURN(0): fi: if N<1 then RETURN(FAIL): fi: if N=1 then if i=0 then RETURN(-1): else RETURN(0): fi: fi: (-CNi(N-1,i)-add(binomial(N,r+1)*CNi(N,i-r),r=1..N-1))/N: end: #C01(j,N): The Rademacher coeff. C_{01j}(N) #computed by the recurrence C01:=proc(j,N): CNi(N,N-j): end: #CNif(N,i,Dig): the coeff. of z^i in 1/*(1+(1+z))*(1+(1+z)+(1+z+z^2)*...)) #in Floating-point CNif:=proc(N,i,Dig) local r: option remember: Digits:=Dig: if i<0 then RETURN(0): fi: if i=0 then RETURN(evalf((-1)^N/N!)): fi: if N=0 then RETURN(0): fi: (-CNif(N-1,i,Dig)-add(binomial(N,r+1)*CNif(N,i-r,Dig),r=1..i))/N: end: #C01f(j,N,Dig): The Rademacher coeff. C_{01j}(N) #computed in floating-point C01f:=proc(j,N,Dig): CNif(N,N-j,Dig): end: #C01Seq(j,N0):The first N0 terms of the Rademacher sequence #C_{01j}(N). Try: #C01Seq(1,100); C01Seq:=proc(j,N0) local N: [seq(C01(j,N),N=1..N0)]: end: #C01Seqf(j,N0):The first N0 terms of the Rademacher sequence #C_{01j}(N). Try: #C01Seqf(1,100); C01Seqf:=proc(j,N0) local N: [seq(C01f(j,N),N=1..N0)]: end: #PNzold(N,z): The rational generating function Sum(CNi(N,i)*z^i,i=0..infinity) #For example, try: #PNzold(2,z); PNzold:=proc(N,z) local r: option remember: if N=0 then RETURN(1): else -PNzold(N-1,z)/(add(binomial(N,r+1)*z^r,r=0..N)): fi: end: #PNz(N,z): The rational generating function Sum(CNi(N,i)*z^i,i=0..infinity) #For example, try: #PNz(2,z); PNz:=proc(N,z) option remember: if N=0 then RETURN(1): else -PNz(N-1,z)/normal (((1+z)^N-1)/z) : fi: end: PaV1:=proc(L) local i: if nops(L)=1 then RETURN([[1,L[1]]]): fi: if L[2]>=L[1] then for i from 3 to nops(L) while L[i]>=L[i-1] do od: RETURN([i-1,L[i-1]]): else for i from 3 to nops(L) while L[i]1 then RETURN(FAIL): fi: n1:=nops(resh)-1: s1:=sn(resh,n1): s2:=sn(resh,n1-1): theta:=evalf(2*(s1+s2)/(s1-s2)^2): mu:=evalf(sqrt(op(n1+1,resh)/op(n1-1,resh))*exp(-(s1+s2)/((s1-s2)*s1))): [theta,mu]: end: #QuasiZinn(L,p): given a list of numbers L, and a pos. integer p #tries to guess, for each n mod p, the appx. asymptotics #in the form C*mu^n*n^theta, the output is a list of triples #[C,theta,mu] for each residue class. #For example, try: #read HANS500: QuasiZinn(lu[1],32); QuasiZinn:=proc(L,p) local lu,mu,i,N,T: N:=nops(L): if N-4*p-1<1 then print(`List too small`): RETURN(FAIL): fi: for i from 0 to p-1 do mu:=evalf([L[N-3*p-i],L[N-2*p-i],L[N-p-i], L[N-i]]): lu:=Zinn(mu): T[(N-i) mod p]:=[L[N-i]/lu[2]^((N-i)/p)/(N-i)^lu[1], lu[1],lu[2]^(1/p)]: od: [seq(T[i],i=1..p-1),T[0]]: end: #ChklNdirect(h,k,l,N): the coeff. of 1/(x-exp(2*Pi*h/k))^l in #1/(1-x)/(1-x^2)/.../(1-x^N). For example, try: #ChklNdirect(0,1,1,5); ChklNdirect:=proc(h,k,l,N) local gu,i,w,h1,x,gu1,r,co: co:=0: w:=exp(2*Pi*I/k): gu:=1: for i from 1 to N do if i mod k<>0 then gu:=gu/(x^i-1): else gu1:=1: co:=co+1: for r from 0 to i-1 do if r/i<>h/k then gu1:=gu1*(x-exp(2*Pi*I*r/i)): fi: od: gu:=gu/gu1: fi: od: gu:=subs(x=x+w^h,gu): (-1)^N*coeff(taylor(gu,x=0,co-l+3),x,co-l): end: #C1hklN(h,k,l,N): the coeff. of 1/(x-exp(2*Pi*h/k))^(trunc(N/k)-l) in #1/(1-x)/(1-x^2)/.../(1-x^N). Using the recurrence #For example, try: #C1hklN(0,1,1,5); C1hklN:=proc(h,k,l,N) local om,a,lu,r: option remember: om:=exp(2*Pi*I*h/k): if N=0 then if l=0 then RETURN(1): else RETURN(0): fi: fi: if l<0 then RETURN(0): fi: if N mod k=0 then lu:=om/N*C1hklN(h,k,l,N-1)- add(binomial(N,a+1)/N/om^a*C1hklN(h,k,l-a,N),a=1..l): lu:=evalc(lu): RETURN(lu): else r:=N mod k: lu:=C1hklN(h,k,l,N-1) - add(binomial(N,a)*om^(r-a)*C1hklN(h,k,l-a,N), a=1..l): lu:=lu/(om^r-1): lu:=evalc(lu): RETURN(lu): fi: end: #ChklN(h,k,l,N): the coeff. of 1/(x-exp(2*Pi*h/k))^l in #1/(1-x)/(1-x^2)/.../(1-x^N). Using the recurrence #For example, try: #ChklN(1,2,1,5); ChklN:=proc(h,k,l,N): (-1)^N*C1hklN(h,k,trunc(N/k)-l,N): end: #C1hklNfloat(h,k,l,N,Dig): the coeff. of 1/(x-exp(2*Pi*h/k))^(trunc(N/k)-l) in #1/(1-x)/(1-x^2)/.../(1-x^N). Using the recurrence with floating point #with Dig Digits #For example, try: #C1hklNfloat(0,1,1,5,20); C1hklNfloat:=proc(h,k,l,N,Dig) local om,a,lu,r: option remember: Digits:=Dig: om:=evalf(exp(2*Pi*I*h/k)): if N=0 then if l=0 then RETURN(1.): else RETURN(0): fi: fi: if l<0 then RETURN(0): fi: if N mod k=0 then lu:=om/N*C1hklNfloat(h,k,l,N-1,Dig)- add(binomial(N,a+1)/N/om^a*C1hklNfloat(h,k,l-a,N,Dig),a=1..l): lu:=evalf(lu): RETURN(lu): else r:=N mod k: lu:=C1hklNfloat(h,k,l,N-1,Dig) - add(binomial(N,a)*om^(r-a)*C1hklNfloat(h,k,l-a,N,Dig), a=1..l): lu:=lu/(om^r-1): lu:=evalf(lu): RETURN(lu): fi: end: #ChklNfloat(h,k,l,N,Dig): the coeff. of 1/(x-exp(2*Pi*h/k))^l in #1/(1-x)/(1-x^2)/.../(1-x^N). Using the recurrence #For example, try: #ChklNfloat(1,2,1,5,30); ChklNfloat:=proc(h,k,l,N,Dig): (-1)^N*C1hklNfloat(h,k,trunc(N/k)-l,N,Dig): end: #ChklSeq(h,k,l,N0):The first N0 terms of the Rademacher sequence #C_{hkl}(0,1,1,60). Try: #ChklSeq(0,1,1,100); ChklSeq:=proc(h,k,l,N0) local N: [seq(ChklN(h,k,l,N),N=1..N0)]: end: #ChklSeqF(h,k,l,N0,Dig):The first N0 terms of the Rademacher sequence #C_{hkl}(0,1,1,60) in floating-point. Try: #ChklSeqF(0,1,1,100,20); ChklSeqF:=proc(h,k,l,N0,Dig) local N: [seq(evalf(ChklNfloat(h,k,l,N,2*Dig),Dig),N=1..N0)]: end: #shk(h,k): The Dedekind sum Eq. (5.2.6) in George Andrews's #The Theory of Partitions" (Chap. 5.2, p. 72) #For example, try: #shk(11,5); shk:=proc(h,k) local m: add((m/k-trunc(m/k)-1/2)*(h*m/k-trunc(h*m/k)-1/2),m=1..k-1): end: #exp(Pi*I*shk(h,k)): The Dedekind sum Eq. (5.2.6) in George Andrews's whkF:=proc(h,k): exp(Pi*I*shk(h,k)): end: #whk(h,k): The 24-th root of unity of Rademacher whk:=proc(h,k) local h1: for h1 from 1 to k-1 while h*h1+1 mod k<>0 do od: if h mod 2=1 then legendre(-k,h)* exp( -Pi*I* ( 1/4*(2-h*k-h) +1/12*(k-1/k)*(2*h-h1+h^2*h1) ) ): elif k mod 2=1 then legendre(-h,k)* exp( -Pi*I* ( 1/4*(k-1) +1/12*(k-1/k)*(2*h-h1+h^2*h1) ) ): fi: end: #L32(y1): L_{3/2} in Rademacher's formula on p.302 #of his book "Topics in Analytic Number Theory" #(Springer, 1973, Grundlehren series Band 169) #For example, try #L32(-4); L32:=proc(y1) local y: if evalf(y1)>=0 then print(`argument of L32 must he negative`): RETURN(FAIL): fi: y:=sqrt(-y1): -(2*cos(2*y)-sin(2*y)/y)/2/sqrt(Pi)/y^2: end: #ChklInfinity(h,k,l): Rademacher's definition of C[h,k,l](infinity) #Eq. (130.6) in p. 302 #of his book "Topics in Analytic Number Theory" #(Springer, 1973, Grundlehren series Band 169) #For example, try #For example, try #ChklInfinity(0,1,1); ChklInfinity:=proc(h,k,l) local gu,a,mu,alpha: gu:=-2*Pi*(Pi/12)^(3/2)*whkF(h,k)*exp(2*Pi*I*h*l/k)/k^(5/2): alpha:=1/24: mu:=add(binomial(l-1,a)*(-1)^(l-1-a)*L32(-Pi^2/6/k^2*(alpha+1+a)), a=0..l-1): simplify(evalc(gu*mu)): end: #ClosestN(h,k,l,K): The integer N between 1 and K that comes #closest to ChklInfinity(h,k,l), together with the difference #For example, try: #ClosestN(0,1,1,100); ClosestN:=proc(h,k,l,K) local lu,i,halev,aluf,si,si1: lu:=ChklInfinity(h,k,l): aluf:=1: si:=evalf(abs(ChklNfloat(h,k,l,1)-lu)): if lu<>0 then si1:=evalf(abs(ChklNfloat(h,k,l,1)/lu)): else si1:=FAIL: fi: for i from 2 to K do halev:=evalf(abs(ChklNfloat(h,k,l,i)-lu)): if halev0 then si1:=evalf(abs(ChklNfloat(h,k,l,i)/lu)): else si1:=FAIL: fi: fi: od: aluf,si,si1: end: #CloseEncounters(K1,L1,N0): The list of the closest #approach to Rademacher's C_{hkl}(infinity) #together with the deviations for all 0<=h<=k<=K1 #1<=l<=L1 for N<=N0. For example try: #CloseEncounters(3,4,50); CloseEncounters:=proc(K1,L1,N0) local N,i,x,h,k,l,lu,t0: t0:=time(): print(`Close Encounters of the Rademacher kind`): print(``): print(`By Shalosh B. Ekhad `): print(``): print(`The great analytic number theorist, Hans Rademacher`): print(`in his posthumously published book`): print(`"Topics in Analytic Number Theory", Springer, 1973 `): print(`(pp. 299-302) conjectured that C_{hkl}(N),`): print(`the coefficient of (x-exp(2*Pi*I*h/k)^(-l) `): print(`in the partial fraction decomposition of`): print(1/Prod(1-x^i,i=1..N)): print(`converges to a certain number defined in his book and`): print(` implemented in procedure ChklInfinity(h,k,l) in`): print(`The Maple package HANS, available, free of charge from`): print(`http://www.math.rutgers.edu/~zeilberg/tokhniot/HANS .`): print(`It turns out (empirically), that these sequences diverge very`): print(`badly, and as N goes to infinity, sooner or later they oscillate`): print(`widely with peaks that go to infinity (exponentially fast!) and`): print(`valleys that go to negative infinity (also exponentially fast!)`): print(``): print(`But it seems that one can tweak Hans's conjecture`): print(`if instead of going all the way to infinity, you only go so far.`): print(`In this article we will find, for h and k, 1<=h1 then RETURN(`ERROR: h and k must be relatively prime`) fi: if r>= k then RETURN(`ERROR: r must be in the range 0,1,2,...,k-1`) fi: rcf:=[ ]; for j from 1 to R do eq:= RecEqn(h,k,r,n,j,C); adjust:={seq( C[h,k,n-y+1,k*n+r]=op(y,rcf),y=1..j-1)}; adjust:={op(adjust), C[h,k,n-j,k*n-k+r]=f(n-1), C[h,k,n-j+1,k*n+r]=f(n)}; eq:= subs(adjust,eq); # init:= f(j) = coeff(Chkgf(k*j+r,x,h,k),x-RootOf1(h,k),j-1); init:= f(j) = simplify(ChklN(h,k,1,k*j+r)): rcf:=[op(rcf), rsolve({eq,init}, f(n)) ] od; rcf end proc: RootOf1:=proc(h,k); cos(2*Pi*h/k) + I*sin(2*Pi*h/k) end proc: ## Automatically generates a paper which states and proves ## formulas for Rademacher coefficients C[h,k, n-j](k*n+r) ## for j from 0 to R-1 HansTopDownAutoPaper:=proc(h,k,r,n,C,R) local t0,i,results,j, G, N, Theorem, l, x; t0:=time(); print(` TOP DOWN FORMULAS FOR SOME RADEMACHER COEFFICIENTS `): print(` by Tapuah I. Sefer, Georgia Southern University `): print(` Dedicated to my mentor, Professor Shalosh B. Ekhad `): print(` `): print(`Abstract: In this paper, I will prove formulas for several`): print(`Rademacher coefficients.`): print(` `): print(`Definition: Let`,C[h,k,l](N) ,` be the numerator of the term `): print(`whose denominator is, ` (x-exp(2*Pi*h*I/k))^l, `in the partial fraction`): print(`decomposition of`, Product( 1/(1-x^j), j=1..N)): results:=ChkFormulaList(h,k,r,n,R); print(`Definition: Let`, G[N]=(x-RootOf1(h,k))^(floor(N/k))/Product(1-x^j, j=1..N)); print(`Thus G satsifies the following recurrence (R):`); print( mul(1-x^(j+r), j=1..k) * G[k*n+r] = (x-RootOf1(h,k)) * G[k*n-k+r] ): print(` Expand both sides of (R) as a Taylor series about x=`, RootOf1(h,k) ): for i from 1 to R do print(`-------------------------------------------------------------------`): print(`Theorem`, i): print(C[h,k,n+1-i](k*n+r)= op(i,results)): print(`Proof :`): print(`Comparing coefficients of`, (x-RootOf1(h,k))^i, `on both sides of (R)` ): print(`we obtain the recurrence (R1) :`): print(RecEqn(h,k,r,n,i,C) ): print(`Solve this recurrence (R1) with the initial condition`, C[h,k,i-i+1](k*i+r) = ChklN(h,k,1,k*i+r) ): ####formerly, coeff(Chkgf(k*i+r,x,h,k),x-RootOf1(h,k),i-1)): if i>1 then print(`taking into account the results of Theorems`, seq(j, j=1..i-1) ) fi: print(`and the result follows.`) od: print(`This paper was produced in`, time()-t0, `seconds of CPU time.`): end proc: ################C01Table #OSC01(v,N): applies one step in the recurrence with vectors for #CNi OSC01:=proc(v,N) local r,a,kha,v1: v1:=[]: for r from 0 to nops(v)-1 do kha:=-v[r+1]/N-add(binomial(N,a+1)/N*v1[r-a+1],a=1..r): v1:=[op(v1),kha]: od: v1: end: #C01Table(K,L): A table of C_{01l}(N) for l from 1 to L #and N from 1 to K, it is given as a list of lists #such that T[l][N] is C[0,1,l](N). For example, try: #C01Table(50,10); C01Table:=proc(K,L) local v,T,N,j: v:=[1,0$K]: T[1]:=[]: for j from 2 to L do T[j]:=[]: od: for N from 1 to K do v:=OSC01(v,N ): T[1]:=[op(T[1]),v[N]]: for j from 2 to L do if N-j<=-1 then T[j]:=[op(T[j]),0]: else T[j]:=[op(T[j]),v[N-j+1]]: fi: od: od: [seq(T[j],j=1..L)]: end: #OSC01f(v,N,Dig): applies one step in the recurrence with vectors for #CNi OSC01f:=proc(v,N,Dig) local r,a,kha,v1: v1:=[]: Digits:=Dig: for r from 0 to nops(v)-1 do kha:=-v[r+1]/N-add(binomial(N,a+1)/N*v1[r-a+1],a=1..r): kha:=evalf(kha): v1:=[op(v1),kha]: od: v1: end: #C01TableF(K,L,Dig): A table of C_{01l}(N) for l from 1 to L #and N from 1 to K, it is given as a list of lists #such that T[l][N] is C[0,1,l](N) C01TableF:=proc(K,L,Dig) local v,T,N,j: v:=[1.,0$K]: T[1]:=[]: for j from 2 to L do T[j]:=[]: od: for N from 1 to K do v:=OSC01f(v,N,Dig ): T[1]:=[op(T[1]),v[N]]: for j from 2 to L do if N-j<=-1 then T[j]:=[op(T[j]),0]: else T[j]:=[op(T[j]),v[N-j+1]]: fi: od: od: [seq(T[j],j=1..L)]: end: ################End C01Table #AsyB:=proc(L,n): Given a list of positive numbers, tries #to find an expression for them as C*mu^n*n^theta #it returns the approx. conjenctured asymptotics #in the form C*mu^n*n^theta. #For example, try: #AsyB([seq(5*2^i*i^2,i=1..10)],n): AsyB:=proc(L,n) local L1,s,eq,var,C1,mu1,t1,n1,i: if nops(L)<6 then print(`List too short`): RETURN(FAIL): fi: if nops({seq(sign(L[i]),i=nops(L)-5..nops(L))})<>1 then RETURN(FAIL): fi: if L[1]<0 then L1:=[seq(-L[i],i=1..nops(L))]: s:=-1: else L1:=L: s:=1: fi: eq:=evalf({seq(log(L1[n1])-C1-n1*mu1-t1*log(n1),n1=nops(L1)-2..nops(L1))}): var:={C1,mu1,t1}: var:=solve(eq,var): C1:=subs(var,C1): C1:=exp(C1): if abs(coeff(C1,I,1))>10^(-8) then RETURN(FAIL): fi: C1:=evalf(coeff(C1,I,0),10): mu1:=evalf(subs(var,mu1),10): t1:=evalf(subs(var,t1),10): s*C1*(exp(mu1))^n*n^t1: end: #Ana1(L): analyzes the pattern of differences between peaks and valleys #of the sequence L, for example, try #Ana1(HANS700[1]); Ana1:=proc(L) local lu,i: lu:=PaV(L): lu:=[seq(lu[i][1],i=1..nops(lu))]: [seq(lu[i+1]-lu[i],i=1..nops(lu)-1)]: end: #Ana(L,n): inputs a sequence of real numbers L, and outputs #conjenctured asymptotics in terms of a list of length #P, where P is the "period" and the ith-term of the #output list gives an estimate for the asymptotic behavior of #L[P*n+i], for i from 1 to P-1, and the P-th term an #for L[P*n]. #It also returns the list [max,min] where #max is the integer mod P where the local maxima occur #and min is the integer mod P where the local minima occur #For example, try: #read HANS800: Ana(HANS800[1],n); Ana:=proc(L,n) local lu,i,P,L1,n1,vu,gu,lu1,ru1,ru2,ru: L1:=evalf(L): lu:=PaV(L1): lu1:=[seq(lu[i][1],i=1..nops(lu))]: lu1:=[seq(lu1[i+1]-lu1[i],i=1..nops(lu)-1)]: if nops(lu1)<6 then print(`Not long enough`): RETURN(FAIL): fi: P:={op(nops(lu1)-4..nops(lu1),lu1)}: if nops(P)<>1 then RETURN(FAIL): fi: ru1:=lu[nops(lu)-1][1]: ru2:=lu[nops(lu)-2][1]: P:=2*P[1]: if lu[nops(lu)-1][2]>0 then ru:=[ru1 mod P, ru2 mod P]: else ru:=[ru2 mod P, ru1 mod P]: fi: gu:=[]: for i from 1 to P-1 do lu:=[seq(L1[P*n1+i],n1=1..trunc((nops(L1)-i)/P))]: vu:=AsyB(lu,n): gu:=[op(gu),vu]: od: lu:=[seq(L1[P*n1],n1=1..trunc(nops(L1)/P))]: vu:=AsyB(lu,n): gu:=[op(gu),vu]: gu,ru: end: #C01TableMinus(K,L): A table of C_{01l}(N) for l from 0 to -L #and N from 1 to K, it is given as a list of lists #such that T[l+1][N] is C[0,1,-l](N). For example, try: #C01TableMinus(50,10); C01TableMinus:=proc(K,L) local v,T,N,j: v:=[1,0$(K+L)]: for j from 0 to L do T[-j]:=[]: od: for N from 1 to K do v:=OSC01(v,N ): T[0]:=[op(T[0]),v[N+1]]: for j from 1 to L do T[-j]:=[op(T[-j]),v[N+j]]: od: od: [seq(T[-j],j=0..L)]: end: #SipurAsy01(n,N): The guessed asymptotics of C01l(32n+r) for l from 1 to 15 and r from 0 to 31. #Type: #SipurAsy01(n,N): SipurAsy01:=proc(n) local gu,l,r,ma,mi,P: read HANSC01f: print(`Conjectured Asymptotic Formulas for the Rademacher Numbers C[0,1,l](N) for l from 1 to 15`): print(``): print(`By Shalosh B. Ekhad`): print(``): for l from 1 to 15 do print(`-------------------------------------------------------`): print(`Doing l=`, l): gu:=Ana(C01f[l],n): if gu=FAIL then next: fi: P:=nops(gu[1]): ma:=gu[2][1]: mi:=gu[2][2]: gu:=gu[1]: print(`The local maxima of`, C[0,1,l](N) , `seem to be when N is congruent to`, ma, `modulo `, P): if ma<>0 then print(`The asymptotic behavior of `, C[0,1,l](P*n+ma), `is approximately`): print(gu[ma]): else print(`The asymptotic behavior of `, C[0,1,l](P*n+ma), `is approximately`): print(gu[P]): print(`and in Maple input format`): lprint(gu[P]): fi: print(`The local minima of`, C[0,1,l](N) , `seem to be when N is congruent to`, mi, `modulo `, P): if mi<>0 then print(`The asymptotic behavior of `, C[0,1,l](P*n+mi), `is approximately`): print(gu[mi]): else print(`The asymptotic behavior of `, C[0,1,l](P*n+mi), `is approximately`): print(gu[P]): print(`and in Maple input format`): lprint(gu[P]): fi: od : print(`We will now list conjectured approximate asymptotic formulas for all residue classes mod`, P): for l from 1 to 10 do print(`-------------------------------------------------------`): print(`Right now l=`, l): gu:=Ana(C01f[l],n): if gu=FAIL then next: fi: P:=nops(gu[1]): gu:=gu[1]: if gu[P]<>FAIL then print(`The asymptotic behavior of `, C[0,1,l](P*n), `is approximately`): print(gu[P]): print(`and in Maple input format`): lprint(gu[P]): fi: for r from 1 to P-1 do if gu[r]<>FAIL then print(`The asymptotic behavior of `, C[0,1,l](P*n+r), `is approximately`): print(gu[r]): print(`and in Maple input format`): lprint(gu[r]): fi: od : od: end: