# # jacks - accelerated computation of Jack symmetric functions # # Calling Sequence: # jacks(n); # jacks(n,verbose); # # Parameters: # n - a positive integer # # Given a positive integer n, jacks(n) first uses add_basis (if necessary) # to insert Jack symmetric functions into the current SF session, using 'J' # as the basis name, and 'a' as the Jack parameter. It then computes the # e-function expansion of J[mu] for each partition of n, and caches the # results in a remember table. One can then process symmetric function # expressions involving the J[mu]'s in the same way as with any other # symmetric functions defined by SF. # # The algorithm is *much* faster and more space efficient than what would # be achieved by using add_basis alone. In particular, it exploits # (1) the duality between J[mu] and J[mu'], where ' denotes conjugation. # (2) the fact that the coefficients of J[mu] are polynomials in 'a' of # degree at most n-length(mu). # (3) the fact that there is a known formula for the norm of J[mu]. # # In the second (verbose) form, diagnostic information is printf'ed each # time a new function J[mu] has been decomposed into the e-basis. # # Reference: # I. Macdonald, "Symmetric Functions and Hall Polynomials", Sec. VI.10. # # Example: # with(SF); # jacks(6,verbose); # tos(J[3,2,1]); jacks:=proc(n) local st,st0,j,vars,mu,nu; interface(quiet=true); gc(); st0:=time(); if not member(J[],`SF/Bases`) then SF['add_basis'](J,mu->SF['zee'](mu,a),mu->SF['hooks'](mu,a)) fi; `SF/added`(J,[]); # avoid clobbering the remember table vars:=[seq(cat('e',j),j=1..n)]; nu:=subs(0=NULL,[n]); while nu<>NULL do st:=time(); mu:=SF['conjugate'](nu); if `jacks/lex`(mu,nu) then `SF/added/gs`(J,nu):=[0,`jacks/dual`(mu,nu,n,vars)] else `SF/added/gs`(J,nu):=[`jacks/norm`(mu,nu),`jacks/gs`(mu,nu,n)] fi; if nargs>1 then printf(`Time %7.3f Space %8d %a\n`, time()-st,`SF/bytes`(),mu) fi; nu:=SF['nextPar'](nu); od; interface(quiet=false); printf(`Done with n=%-2d Time %9.3f Space %8d\n`, n,time()-st0,`SF/bytes`()); end; # Return true if mu < nu in the lex ordering used by 'Par'. `jacks/lex`:=proc(mu,nu) local i,d; d:=zip((x,y)->x-y,mu,nu,0); for i in d do if i<>0 then RETURN(evalb(i>0)) fi od; false; end: # The squared norm of J[mu] `jacks/norm`:=proc(mu,nu) local i,j; convert([seq(seq(a*(mu[i]-j)+nu[j]-i+1,j=1..mu[i]),i=1..nops(mu)), seq(seq(nu[i]-j+a*(mu[j]-i+1),j=1..nu[i]),i=1..nops(nu))],`*`) end: # Compute J[mu] via the Gram-Schmidt algorithm, but use taylor series # up to degree n-nops(mu) in place of "normal" field arithmetic. # Might want to tinker with the length trigger (25000) for *large* n. `jacks/gs`:=proc(mu,nu,n) local g,lc,f,f0,sc,tau,c,N,m; lc:=`SF/lcoeff`[J](mu); f:=lc*e[op(nu)]; g:=`SF/added/etable`(nu)[2]; sc:=[]; N:=n-nops(mu)+1; tau:=subs(0=NULL,[n]); while tau<>nu do if SF['dominate'](tau,nu) then c:=`SF/added/etable`(tau)[2]; c:=SF['scalar'](g,c,0,0,`SF/iprod`[J]); m:=a^ldegree(c,a); c:=m*expand(numer(c)/m)/denom(c); sc:=[op(sc),e[op(tau)]=c]; f0:=`SF/added/gs`(J,tau); c:=subs(sc,f0[2]); m:=ldegree(f0[1],a); if type(c,`+`) then c:=a^m*map((x,y)->x/y,c,a^m) fi; c:=convert(taylor(lc*c/f0[1],a,N),'polynom'); f:=f-c*f0[2]; if length(f)>25000 then f:=convert(taylor(f,a,N),'polynom') fi; fi; tau:=SF['nextPar'](tau) od; f:=convert(taylor(f,a,N),'polynom'); collect(f,indets(f,'indexed')); end: # Compute the e-expansion of J[mu] by twisting the e-expansion # of J[nu]=J[mu']. Again use the taylor series trick. `jacks/dual`:=proc(mu,nu,n,vars) local f,j; f:=`SF/added`(J,mu); f:=subs({seq(vars[j]=`jacks/twist`(j),j=1..n)},f); f:=taylor((-1)^n*f,a,n+1); f:=[seq(expand(coeff(f,a,j))*a^(n-j),j=nops(mu)..n)]; collect(convert(f,`+`),vars,'distributed'); end: # Remember the e-expansion of theta(e.n,-a). `jacks/twist`:=proc(n) local j,f; option remember; f:=SF['top'](cat('e',n)); f:=subs({seq(cat('p',j)=-a*cat('p',j),j=1..n)},f); SF['toe'](f,'p'); end: