#C26.txt, Dr. Z.'s Experimental Mathematics Class (Math 640), Fall 2016, Dec. 5, 2016
Help:=proc(): print(`Da(f,x,n,a), ToMv(f,x,n) , DaM(a,n), JackM(a,n) , Dqt(f,x,n,q,t) `): 
print(` DqtM(q,t,n), MacM(q,t,n) `):
end:

with(linalg):
#MacM(q,t,n): the p(n) by p(n) matrix whose i-j entry is the coeff. of M[P(n)[j]] of
 #Macdonald_{q,t}[P(n)[i]]:
 MacM:=proc(q,t,n) local S,i:
  S:=[eigenvectors(DqtM(q,t,n))]:
  S:={seq(S[i][3][1],i=1..nops(S))}:
 end:

#ToMv(f,n,M):inputs a linear combination of M[p] and outputs its 
ToMv:=proc(f,x,n) local M,i,f1,Pa:
 f1:=ToM(expand(f),x,n,M):
 Pa:=P(n):
 [seq(coeff(f1,M[Pa[i]],1),i=1..nops(Pa))]:
end:

#Dqt(f,x,n,q,t): the differential operator used to define Macdonald polynomials
Dqt:=proc(f,x,n,q,t) local i,j:
normal(
 add(mul((t*x[i]-x[j])/(x[i]-x[j]),j=1..i-1)*
mul((t*x[i]-x[j])/(x[i]-x[j]),j=i+1..n)*
subs(x[i]=q*x[i],f),i=1..n)):
 end:

#DqtM(q,t,n): the matrix version of Dqt, a p(n) by p(n) matrix whose i-th row
#is the vector of coefficients of P(n) of Da(mL(P(n)[i]):
DqtM:=proc(q,t,n) local i,Pa:
 Pa:=P(n):
  transpose([seq(ToMv(Dqt(mL(x,Pa[i],n),x,n,q,t),x,n), i=1..nops(Pa))]):
end:

#Da(f,x,n,a): the differential operator used to define Jack polynomials
Da:=proc(f,x,n,a) local i,j:
normal(
a/2*add(x[i]^2*diff(f,x[i]$2),i=1..n)+
add(add(x[i]^2/(x[i]-x[j])*diff(f,x[i]),j=1..i-1)+
add(x[i]^2/(x[i]-x[j])*diff(f,x[i]),j=i+1..n),i=1..n)):
end:

#DaM(a,n): the matrix version of Da, a p(n) by p(n) matrix whose i-th row
#is the vector of coefficients of P(n) of Da(mL(P(n)[i]):
DaM:=proc(a,n) local i,Pa:
 Pa:=P(n):
  transpose([seq(ToMv(Da(mL(x,Pa[i],n),x,n,a),x,n), i=1..nops(Pa))]):
end:

 #JackM(a,n): the p(n) by p(n) matrix whose i-j entry is the coeff. of M[P(n)[j]] of
 #J[P(n)[i]]:
 JackM:=proc(a,n) local S,i:
  S:=[eigenvectors(DaM(a,n))]:
  S:={seq(S[i][3][1],i=1..nops(S))}:
 end:

##old stuff
#C24.txt: Schur Functions and beyond, Nov. 28, 2016
Help24:=proc() print(` MonToV(Mo,H,n), PolToV(Pol,H,n), sToHvector(L), sToMvector(L) , ProveOrt(n)`): end:

#ProveOrt(n): inputs a pos. integer n and outputs true if  (4.8) on p. 35 of Macdonald (1st)
#is true for n 
ProveOrt:=proc(n) local Pa,i,j,v1,v2,k:
Pa:=P(n):
 for i from 1 to nops(Pa) do
   for j from 1 to nops(Pa) do
   v1:=sToHvector(Pa[i]):
    v2:=sToMvector(Pa[j]):
     if i=j then
    if add(v1[k]*v2[k],k=1..nops(v1))<>1 then
     RETURN(false):
    fi:
     elif
   add(v1[k]*v2[k],k=1..nops(v1))<>0 then
       RETURN(false):
    fi:
  od:
  od:
true:
end:


  

#sToMvector(L):inputs a partition L and outouts a list of numbers giving the
#coeff. of M[P(n)[i]]:
sToMvector:=proc(L) local n,Pa,f,M,i:
 n:=convert(L,`+`):
 Pa:=P(n):
 f:=sLm(L,M):
 [seq(coeff(f,M[Pa[i]],1),i=1..nops(Pa))]:
end:

#sToHvector(L):inputs a partition L and outouts a list of numbers giving the
#coeff. of H[P(n)[i]]:
sToHvector:=proc(L) local n,Pa,f,H,i:
 n:=convert(L,`+`):
 Pa:=P(n):
 f:=sLhF(L,H):
 [seq(coeff(f,H[Pa[i]],1),i=1..nops(Pa))]:
end:

#MonToV(Mo,H,n): inputs a monomial Mo in H[1], ..., H[n] and outputs c*H[L] where
#c is the coefficient and L is the decreasing degree sequence
MonToV:=proc(Mo,H,n) local L,c:
L:=[seq(degree(Mo,H[i]),i=1..n)]:
c:=normal(Mo/mul(H[i]^L[i],i=1..n)):
 if not type(c,numeric) then
   print(`Not a monomial`):
   RETURN(FAIL):
 fi:
c*H[[seq((n-i+1)$L[n-i+1],i=1..n)]]:
end:

#PolToV(Pol,H,n): inputs a polynomial Pol in H[1], ..., H[n] outputs it as
#linear combination of H[L]
PolToV:=proc(Pol,H,n)  local Pol1:
  Pol1:=expand(Pol):
 if not type(Pol1,`+`) then
    RETURN(MonToV(Pol1,H,n)):
 fi:

  add(MonToV(op(i,Pol1),H,n),i=1..nops(Pol1)):
  end:

#Old stuff from  previous classes
HelpOld:=proc(): print(`  sLm(L,M), sLe(L,E),sLeF(L,E), sLh(L,H), sLhF(L,H),poL(x,L,n) `): 
end:
with(linalg):

#sLhF(L,H): sLhF(L,H) via the Jacobi-Trudy formula
sLhF:=proc(L,H) local i,j,J,n:
 n:=convert(L,`+`):
J:=expand(det([seq([seq(H[L[i]-i+j],i=1..nops(L))],j=1..nops(L))])):
 J:=subs({seq(H[-i]=0,i=1..n), H[0]=1},J):
 PolToV(J,H,n):
end:


#sLeF(L,E): sLe(L,E) via the Jacobi-Trudy formula
sLeF:=proc(L,E) local i,j,L1,J,n:
L1:=Conj(L):
 n:=convert(L1,`+`):
 J:=expand(det([seq([seq(E[L1[i]-i+j],i=1..nops(L1))],j=1..nops(L1))])):
 J:=subs({seq(E[-i]=0,i=1..n), E[0]=1},J):
 PolToV(J,E,n):
end:

#inv1(pi): the number of inversions of the permutation pi
inv1:=proc(pi) local i,j,n,co:
 n:=nops(pi):
 co:=0:

   for i from 1 to n do
     for j from i+1 to n do
       if pi[i]>pi[j] then
         co:=co+1:
       fi:
     od:
od:

 co:
end:


with(combinat):
#AS(f,x,n): the anti-symmstrizer of the function f of x[1], ..., x[n]
AS:=proc(f,x,n) local S:

S:=permute(n):

add((-1)^inv1(pi)* ApplyPi(f,x,n,pi),pi in S):


end:

   #sL(L,x): Inputs an integer partition L (of n, say) and a letter x
    #outputs the Schur symmetric polynomial s_L(x1, ..., xn) 
   sL:=proc(L,x) local n,L1, rich,i:
    n:=convert(L,`+`):
    L1:=[op(L),0$(n-nops(L))]:
    rich:=[seq(n-i,i=1..n)]:
    L1:=L1+rich:
    
    normal(AS(mul(x[i]^L1[i],i=1..n),x,n)/mul(mul(x[i]-x[j],j=i+1..n),i=1..n-1)):
    end:
  
 #sLm(L,M): s_L(x) in terms of the Monomial Symmetric Functions
 sLm:=proc(L,M) local n,x:
   n:=convert(L,`+`):
   ToM(sL(L,x),x,n,M):
  end:

 #sLe(L,E): s_L(x) in terms of the e-basis
 sLe:=proc(L,E) local n,x:
   n:=convert(L,`+`):
   ToE(sL(L,x),x,n,E):
  end:
  

 #sLh(L,H): s_L(x) in terms of the h-basis
 sLh:=proc(L,E) local n,x:
   n:=convert(L,`+`):
   ToH(sL(L,x),x,n,E):
  end:
   
#stuff from previous classes
#C18.txt, Nov. 7, 2016, Dr. Z.'s Experimental Mathematics Class, Bases of Symmetric Polynomials

Help18:=proc(): print(`  Toe(f,x,n,e),  ToE(f,x,n,E), Andrew(F,x,n), mToe(n) `): end:



#Stuff from C15.txt ##################

   #C15.txt: ExpMath, Rutgers, Oct. 27, 2016, Bases of the alg. of Symmetric Polynomials
   Help15:=proc(): print(`P(n), ApplyPi(f,x,n,pi), Sym(f,x,n), mL(x,L,n)`):
   print(`IsSym(f,x,n), ToM(f,x,n,M), ei(x,i,n), eL(x,L,n), eTom(n)  `):
   end:
    with(combinat):
   #Rev(L): reverse of L
   Rev:=proc(L) local n,i: n:=nops(L): [seq(L[n-i+1],i=1..n)]: end:
 
    #P(n): the list of (integer) partitons of n in non-increasing order
    P:=proc(n) local S,i:
    S:=partition(n):
    Rev([seq(Rev(S[i]),i=1..nops(S))]):
    end:

    #ApplyPi(f,x,n,pi): inputs f, x, n, pi, where f is an expression in x[1],...,x[n]
    #pos. integer n, letter x, and a permutation pi of {1, ...,n}, applies x[i]->x[pi[i]]
    #to f
    ApplyPi:=proc(f,x,n,pi) local i:
    if nops(pi)<>n then
      RETURN(FAIL):
    fi:

     subs({seq(x[i]=x[pi[i]],i=1..n)},f):
     end:
   
      #Sym(f,x,n): the symmetrizer of f
      Sym:=proc(f,x,n) local pi,S:
      S:=permute(n):
      
       convert({seq(ApplyPi(f,x,n,pi), pi in S)},`+`):


      end:

        #mL(x,L,n): The monomial symmetric function m_L(x1,..., xn) where n>=nops(L)
        mL:=proc(x,L,n) local i:
        if n<nops(L) then
         RETURN(FAIL):
        fi:

         Sym(mul(x[i]^L[i],i=1..nops(L)), x,n):
        
        end:

         #IsSym(f,x,n): inputs f and x (a lettter), f is an expression in x[1],..., x[n], n pos.int.
         #outputs true (or false) iff f is symmetric
         IsSym:=proc(f,x,n) local i:
          for i from 1 to n-1 do
           if subs({x[i]=x[i+1],x[i+1]=x[i]},f)-f<>0 then
             RETURN(false):
           fi:
          od:
           true:
          end:
           
             KickZ:=proc(L) local i:
              for i from 1 to nops(L) while L[i]<>0 do  od:
               [op(1..i-1,L)]:
              end:

                            
             
           #ToM(f,x,n,M): inputs a symmetric polynomial f in x[1]..., x[n]
            #outputs it as a linear combination of M[L]
           ToM:=proc(f,x,n,M) local F,f1,mono,L,c,i:
             
             if not IsSym(f,x,n) then
                RETURN(FAIL):
             fi:

 

              if not type(f,`+`) then
               c:=normal(f/mul(x[i],i=1..n)):
               if not type(c,integer) then
                    RETURN(FAIL):
                else
                RETURN(c*M[[1$n]]):
              fi:
             fi:
             
            f1:=f:
             F:=0:
             while f1<>0 do
               mono:=op(1,f1):
               L:=[seq(degree(mono,x[i]),i=1..n)]:
                c:=normal(mono/mul(x[i]^L[i],i=1..nops(L))):

                L:=KickZ(sort(L,`>`)):

                   F:=F+c*M[L]:
                  f1:=expand(f1-c*mL(x,L,n)):
           
              
               
              od:

             F:

            end:
       #ei(x,i,n): the elementary symmetric function e_i(x[1], ..., x[n])
      ei:=proc(x,i,n) :
       if i>n then
        0:
       else
      mL(x,[1$i],n):
      fi:

         end:
       
 poi:=proc(x,i,n) local j: add(x[j]^i,j=1..n):end:
 

       #poL(x,L,n): the mul(p_L[i], i=1..n)
       poL:=proc(x,L,n) local i:
       expand(mul(poi(x,L[i] ,n),i=1..nops(L))):
        end:

       #eL(x,L,n): the mul(e_L[i], i=1..n)
       eL:=proc(x,L,n) local i:
       mul(ei(x,L[i] ,n),i=1..nops(L)):
        end:

    




#eTom(n): code by JinYoung Park
eTom:=proc(n) local Pn,k,A,B,M,i,M1,j,x:

Pn:=P(n):

k:=nops(Pn):

A:=[seq(M[Pn[i]],i=1..k)]:

B:=[seq(ToM(expand(eL(x,Pn[i],n)),x,n,M),i=1..k)]:

M:=[]:

for i from 1 to nops(A) do
 M1:=[seq(coeff(B[i],A[j]),j=1..nops(A))]:
 M:=[op(M),M1]:
od:
M:
end:


#poTom(n): code by JinYoung Park
poTom:=proc(n) local Pn,k,A,B,M,i,M1,j,x:

Pn:=P(n):

k:=nops(Pn):

A:=[seq(M[Pn[i]],i=1..k)]:

B:=[seq(ToM(expand(poL(x,Pn[i],n)),x,n,M),i=1..k)]:

M:=[]:

for i from 1 to nops(A) do
 M1:=[seq(coeff(B[i],A[j]),j=1..nops(A))]:
 M:=[op(M),M1]:
od:
M:
end:




       #ei(x,i,n): the elementary symmetric function e_i(x[1], ..., x[n])
      ei:=proc(x,i,n) :
       if i>n then
        0:
       else
      mL(x,[1$i],n):
      fi:

         end:
       
#hi(x,i,n): the complete homog. symmetric function h_i(x[1], ..., x[n])
hi:=proc(x,i,n) local lu,lu1:
lu:=P(i):
add( mL(x,lu1,n),lu1 in lu):
end:

       #hL(x,L,n): the mul(h_L[i], i=1..n)
       hL:=proc(x,L,n) local i:
       mul(hi(x,L[i] ,n),i=1..nops(L)):
        end:

######Code from previou classes

#End Stuff from C15.txt ##################

#Start Stuff from C18.txt
#mToe(n): the transition matrix from the m-basis to the e-basis
mToe:=proc(n)
inverse(eTom(n)):
end:


#mTopo(n): the transition matrix from the m-basis to the e-basis
mTopo:=proc(n)
inverse(poTom(n)):
end:


#mToh(n): the transition matrix from the m-basis to the e-basis
mToh:=proc(n)
inverse(hTom(n)):
end:


#Andrew(F,x,n): inputs a symmetric polynomial F, in x[1],..., x[n]  abd ouputs the vector of coefficients
Andrew:=proc(F,x,n) local A,Pars,M,i:

A:=ToM(expand(F),x,n,M):
Pars:=P(n):

 [seq(coeff(A,M[Pars[i]],1),i=1..nops(Pars))]:

end:




#ToH(f,x,n,E): inputs a symmetric function f, homog. of total degree n, of x[1]. ..., x[n]
#and a positive integer n, and a symbol E, expresses it as a linear combination of
#the H-basis, where H[L] means e_L[1]*...e_L[nops(L)]
ToH:=proc(f,x,n,E) local v,Pars,i:

v:=multiply(Andrew(f,x,n),mToh(n)):
Pars:=P(n):

add(v[i]*E[Pars[i]],i=1..nops(Pars)):


end:

#ToE(f,x,n,E): inputs a symmetric function f, homog. of total degree n, of x[1]. ..., x[n]
#and a positive integer n, and a symbol E, expresses it as a linear combination of
#the E-basis, where E[L] means e_L[1]*...e_L[nops(L)]
ToE:=proc(f,x,n,E) local v,Pars,i:

v:=multiply(Andrew(f,x,n),mToe(n)):
Pars:=P(n):

add(v[i]*E[Pars[i]],i=1..nops(Pars)):


end:

#Toe(f,x,n,e): inputs a symmetric function f, homog. of total degree n, of x[1]. ..., x[n]
#and a positive integer n, and a symbol e, expresses it as polynomial expression in e[1],e[2],...,e[n]
#thereby confirming the fundamental theorem of Symmetric Functions that any summetric polynomial
#in x[1],..,x[n] of degree<=n can be thus expressed

Toe:=proc(f,x,n,e) local v,Pars,i,j:

v:=multiply(Andrew(f,x,n),mToe(n)):
Pars:=P(n):

add(v[i]*mul(e[Pars[i][j]],j=1..nops(Pars[i])),i=1..nops(Pars)):

end:




#hTom(n): code by JinYoung Park
hTom:=proc(n) local Pn,k,A,B,M,i,M1,j,x:

Pn:=P(n):

k:=nops(Pn):

A:=[seq(M[Pn[i]],i=1..k)]:

B:=[seq(ToM(expand(hL(x,Pn[i],n)),x,n,M),i=1..k)]:

M:=[]:

for i from 1 to nops(A) do
 M1:=[seq(coeff(B[i],A[j]),j=1..nops(A))]:
 M:=[op(M),M1]:
od:
M:
end:

###end stuff from C18.txt

#Conj(L): the conjugate partition of L  (due to Yonah Biers-Ariel)
Conj:=proc(L) local i:add([1$L[i],0$(L[1]-L[i])],i=1..nops(L)):end: