%---begin plain TeX file------------ 
%%A plain TeX file of a preprint by D. Zeilberger:"The J.C.P. Miller
%% Recurrence for exponentiating a polynomial, and its q-Analog"
 
 
%begin macros
\baselineskip=14pt
\parskip=10pt
\def\ha{{ {1} \over {2} } }
\def\epsilon{\varepsilon}
\def\centerthis{\leftskip=0pt plus 1fill \rightskip=\leftskip \relax}
\def\Re{{\rm Re}\,}
\def\Im{{\rm Im}\,}
\def\halmos{\hbox{\vrule height0.15cm width0.01cm\vbox{\hrule height
 0.01cm width0.2cm \vskip0.15cm \hrule height 0.01cm width0.2cm}\vrule
 height0.15cm width 0.01cm}}
\font\eightrm=cmr8 \font\sixrm=cmr6 
\font\eighttt=cmtt8
\magnification=\magstephalf
\def\lheading#1{\bigbreak \noindent {\bf#1} 
}
\parindent=0pt
\overfullrule=0in
%end macros
\bf
The J.C.P. Miller Recurrence 
for Exponentiating a Polynomial, and its $q$-Analog
 
\it
 \quad \quad \quad
\qquad \qquad \qquad \qquad \qquad \qquad Doron ZEILBERGER \footnote{$^1$}
{\eightrm  \raggedright
Department of Mathematics, Temple University,
Philadelphia, PA 19122, USA. 
{\eighttt zeilberg@math.temple.edu .}
Supported in part by the NSF. I wish to thank Herb Wilf for many
useful comments.
}
 
\rm
In my research on constant term conjectures, I often need to
expand powers of polynomials  $P^m$, where $m$ is very large,
and $P$ is (usually) a polynomial of several variables. I was frustrated
by the slowness of all the commercial computer algebra packages. For example,
in Maple$^{TM}$, it takes several days to do 
{\tt expand (( 1+ 3*x+ 2*x**2 )**3000)}.
When I asked Herb Wilf if there was any faster way of doing it, he
referred me to  section 4.3 if [1] and 
chapter $21$ of [2], but added, ``It is pretty fast,
but not as fast as FFT (Fast Fourier Transform)''.  When I looked it up,
I realized that, for a {\it fixed} polynomial $P(x)$ of degree $L$,
the time-complexity and space-complexity are both superior to FFT!
The time-complexity of doing $P(x)^m$ using FFT is 
$O(mL \log_2 (mL))=O(m \log (m) )$,
and its space-complexity is $O(mL)=O(m)$.
On the other hand, the time- and space- complexity
of expanding $P(x)^m$, using the recurrence of J.C.P. Miller given below,
are $O(L^2 m)=O(m)$ and $O(L)=O(1)$ respectively.  
 
Although Herb Wilf claims
that [1] and [2] sold lots of copies, and nothing in them
can be called {\it obscure},
my own inquiries among the experts of computer algebra revealed that few people
are aware of the method.
The main purpose of this note is to bring this beautiful
and simple method to the attention of the users of computer algebra.
Its second purpose is to point out a $q-$analog, that should
be useful to people, like myself,
who investigate $q-$constant term identities.
 
We all know that the coefficients of $(1+x)^m$ are given by the binomial
coefficients. One way to compute them is via Pascal's triangle. Alas,
this will give us lots of useless information, namely, the coefficients
of $(1+x)^{m'}$ for {\it all} $m' < m$. A much better way is to use
the {\it closed form} formula ${{m} \choose {k}}= m!/(k!(m-k)!)$ for
the coefficient, let's call it $a(m,k)$,
of $x^k$ in $(1+x)^m$. An even better method is,
to use the equivalent statement that 
$a(m,k)={{m-k+1} \over {k}} a(m,k-1)$, starting
with the obvious initial condition $a(m,0)=1$.
 
In the J.C.P. Miller method, one can replace $(1+x)$ by {\it any}
polynomial $P(x)$, and get a recurrence for the coefficients
$a(m,k)$ of $P(x)^m$, whose order equals the degree of $P(x)$, and
for which $m$ stays fixed, so we never have to compute values of
$a(m',k')$ for $m' < m$.
 
Recall that ${{m} \choose {k}}/2^m$, the coefficient of $x^k$ in
$(\ha + \ha x)^m$ is the probability of scoring
$k$ Heads upon tossing $m$ (fair) coins.  Likewise, the probability
of scoring a total of $k$ points, upon rolling a (fair) die
$m$ times, is the coefficients of $x^k$ in $(x+x^2+x^3+x^4+x^5+x^6)^m/6^m$.
More generally, for rolling $m$ times an $L-faced$ die, with faces marked
$a_1 , \dots , a_L$, with probabilities $p_1 , \dots , p_L$, the
probability of scoring a total of $k$ points is the coefficients of $x^k$ in
$P(x)^m$, where $P(x)= p_1 x^{a_1} + \dots + p_L x^{a_L}$.
 
From now, fix a polynomial $P(x):= \sum_{i=0}^{L} p_i x^i $, and let
$P(x)^m := \sum_{k=0}^{mL} a(m,k) x^k$.
The analog of the Pascal triangle recurrence is the following.
 
{\bf Partial Recurrence:} $a(m,k)=\sum_{i=0}^{L} p_i \, a(m-1,k-i)$.
 
{\bf Proof:} Compare the coefficient of $x^k$ in $P(x)^m=P(x)^{m-1}P(x)$. 
\halmos
 
The analog of 
${{m} \choose {k}}={{m-k+1} \over {k}} {{m} \choose {k-1}}$ is the following
$L$-order recurrence, that enables one to recursively compute the coefficients
$a(m,k)$ of $P(x)^m$, starting with the obvious initial values
$a(m,0)=p_0^m$, and $a(m,k)=0$, for $k < 0$. It requires $2L^2 m$ 
additions and $3L^2 m$ multiplications.
 
{\bf The J.C.P. Miller Pure Recurrence}([1][2]): Assume that $p_0 \neq 0$, then
 
$$
a(m,k)= {{1} \over {k p_0}} \sum_{i=1}^{L}
p_i \, [(m+1)i-k] \, a(m,k-i) \quad .
$$
 
{\bf Proof:} The original proof in [2] was similar to the one above:
compare the coefficient of $x^k$ in $(P(x)^m)'=mP'(x)P(x)^{m-1}$.
Another proof, that is easier to q-analogize, is as follows.
$a(m,k)$ is the coefficient of $x^0$ in the Laurent
polynomial $P(x)^m/x^k$. For any Laurent polynomial $f(x)$, we have
that the coefficient of $x^0$ in $x {{d} \over {dx}} f(x)$ equals 0.
Hence
 
$$
0= C_{x^0} [ x {{d} \over {dx}} {{P(x)^{m+1}} \over {x^k}}]=
C_{x^0} \{
-k[ p_0 + p_1 x + \dots + p_L x^L ]  {{P(x)^{m}} \over {x^k}}
+(m+1)[ p_1 + 2 p_2 x + \dots L p_L x^{L-1} ] {{P(x)^{m}} \over {x^{k-1}}}
\}=
$$
$$
C_{x^0} \{
-k[ p_0 {{P(x)^{m}} \over {x^k}} + p_1 {{P(x)^{m}} \over {x^{k-1}}}
+ \dots + p_L {{P(x)^{m}} \over {x^{k-L}}}]  
+(m+1)[ p_1 {{P(x)^{m}} \over {x^{k-1}}}+ 2 p_2 {{P(x)^{m}} \over {x^{k-2}}}
+ \dots L p_L  {{P(x)^{m}} \over {x^{k-L}}}] 
\} \quad ,
$$
 
from which the pure recurrence follows. \halmos
 
Given a polynomial of several variables $P( x_1 , \dots , x_k )$, we
can similarly compute  the coefficient
of $x_1^0 \dots x_k^0$ in 
$ x_1 {{\partial} \over { \partial x_1}} ( P^{m+1}/(x^{n_1}  \dots x^{n_k}) )$,
to get a partial recurrence in $n_1 , \dots , n_k$,
but with {\it fixed} $m$. The
details are left to the readers.
 
Finally, here is a $q$-analog. The $q$-analog of raising a polynomial
to a power, $P(x)^m$, is 
 
$$
P(x) P(qx) P( q^2 x) \dots P( q^{m-1} x) \quad .
$$
 
{\bf The q-analog of the J.C.P. Miller Pure Recurrence:} 
Let $P(x)= \sum_{i=0}^{L} p_i x^i$ be a polynomial,  with $p_0 \neq 0$, and let
 
$$
P(x)P(qx) \dots P( x^{m-1} x) =
\sum_{i=0}^{Lm} a(m,k) x^i \quad.
$$
 
The coefficients $a(m,k)$ can be successively computed from the following
recurrence.
 
$$
a(m,k)= {{(-1) } \over { p_0 (1- q^{-k}) }} \sum_{i=1}^{L}
p_i \, [ q^{-i} - q^{mi-k} ] \, a(m,k-i) \quad .
$$
 
{\bf Proof:} Use the {\it Stembridge-Stanton trick}[3] that for any 
Laurent polynomial
$F(x)$, the coefficient of $x^0$ in $[F(x)-F(qx)]$ equals 0, with
 
$$
F(x)= 
{{P( q^{-1} x) P(x) \dots P( q^{m-1} x) } \over {x^k}} \quad ,
$$
 
and proceed as before. \halmos 
 
{\bf Appendix:} The following are Maple procedures, {\tt Pxm} and
{\tt qPxm}, implementing
the two recurrences.   {\tt Pxm} is much faster than both
{\tt expand(P**m)} and {\tt taylor(P**m,x=0,L*m+1)}.
The function call for {\tt Pxm} is is {\tt Pxm(POLY,x,m)}, where
{\tt POLY} is the polynomial $P(x)$, {\tt x} is the
variable, $x$, and {\tt m} is the power. It outputs $P(x)^m$. 
Similarly for {\tt qPxm}.
The
programs are also available via anonymous {\tt ftp}\footnote{$^2$}
{\eightrm For the sake of novices, in case they still exist,
let me recall the steps. 1) {\eighttt ftp math.temple.edu}
2)You logon as {\eighttt anonymous}, with password [your e-mail address].
3) You type: {\eighttt cd pub/zeilberger/programs}. 4) You type:
{\eighttt get Pxm.maple}. 
(If you want all the programs, type {\eighttt mget *.*})
5) You type {\eighttt quit}.
}
to {\tt math.temple.edu}
in directory {\tt pub/zeilberger/programs},  by {\tt get}ing  file
{\tt Pxm.maple}.
This directory contains many other of the author's Maple programs.
 
{\eighttt
Pxm:=proc(POLY,x,m) local n,g,i,L,p,P: 
P:=expand(POLY/x**ldegree(POLY,x)):L:=degree(P,x): 
for i from 0 to L do p[i]:=coeff(P,x,i): od:
g:=p[0]**m: for n from 1 to L*m do 
g:=g+ sum(((m+1)*'i'-n)*p['i']
* coeff(g,x,n-'i'),'i'=1..L)
/n/p[0]*x**n: od:
g:= expand( g*x**(ldegree(POLY,x)*m)):g: end:}
 
{\eighttt
qPxm:=proc(POLY,x,m) local n,g,i,L,p,P: 
P:=expand(POLY/x**ldegree(POLY,x)):L:=degree(P,x): 
for i from 0 to L do p[i]:=coeff(P,x,i): od:
g:=p[0]**m: for n from 1 to L*m do 
g:=g- normal(sum( (q**(-'i')-q**(m*'i'-n))/(1-q**(-n))*p['i']
* coeff(g,x,n-'i'),'i'=1..L))
/p[0]*x**n: od:
g:= expand( g*x** ( ldegree (POLY,x)*m)):
taylor(g, x=0, degree(g,x)+1): end:}
 
{\bf References}
 
[1] D. E. Knuth, {\it ``Seminumerical Algorithms''}, (``The art
of Computer Programming v. {\bf 2} ''), second edition, Addison Wesley,
Reading, 1981.
 
[2] A. Nijenhuis and H. Wilf, {\it Combinatorial Algorithms}
($2^{nd}$ edition), Academic Press, New York, 1978.
 
[3] D. Zeilberger,
{\it A Stembridge-Stanton style proof of the Habsieger-Kadell
q-Morris identity}, Discrete Math., {\bf 79}, 313-322 (1989/90).
 
{\bf April 18, 1994; 8 Iyar, 5754.}
 
\bye
 
 
 
 
 
 
 
 
