# Finding recursive relations in Mathematica

1. Aug 13, 2011

### mnb96

Hello,
I have a sequence of polynomials defined in the following way:
$$P_k(x) = \frac{\partial^k}{\partial x^k}e^{s(x)}\vert_{x=0}$$

Essentially the polynomial Pk is the k-th derivative of $\exp(s(x))$ evaluated at x=0. The function s(x) is a polynomial of 2nd degree in x.

In mathematica I define the polynomials with the following code:

D[Exp[s[x]], {x, k}]

For k=1..n one obtains a list of n polynomials P1(x),...,Pn(x).
My question is: is it possible to ask Mathematica to find a recurrent relation that expresses any Pk as a function of the previous polynomials Pk-1, Pk-2... ? (assuming it exists).

A similar well-known problem exists for http://en.wikipedia.org/wiki/Hermite_polynomials#Definition"

Thanks.

Last edited by a moderator: Apr 26, 2017
2. Aug 13, 2011

### Simon_Tyler

Try

Code (Text):
In[1]:= s[x_] := a x^2+b x+c

In[2]:= poly[k_] = Simplify[SeriesCoefficient[Exp[s[x]], {x, 0, k}], k>0]

Out[2]= DifferenceRoot[Function[{y,n},
{-2 a y[n] - b y[1+n] + (2+n) y[2+n]==0, y[0]==E^c,  y[1]==b E^c} ]]

The above shows the difference equation (recurrence relation) and the initial conditions that the polynomials must satisfy for any given quadratic s[x]. We can now use the recurrence relation to calculate the sequence. Below we compare calculating using poly[n] (the DifferenceRoot object) and Recurrence table

Code (Text):

In[11]:= dr = Table[poly[n], {n, 0, 20}]; // Timing
drExp = Expand[dr]; // Timing

Out[11]= {0.21, Null}
Out[12]= {1.57, Null}

In[13]:= rt = RecurrenceTable[{-2 a y[n] - b y[1 + n] + (2 + n) y[2 + n] == 0,
y[0] == E^c, y[1] == b E^c}, y, {n, 0, 20}]; // Timing
rtExp = Expand[rt]; // Timing

Out[13]= {0.01, Null}
Out[14]= {1.66, Null}

In[15]:= dr == rt
Out[15]= True

In[16]:= drExp == rtExp
Out[16]= True

The problem is that both evaluating the DifferenceRoot poly[k] and using RecurrenceTable, whilst great for purely numeric functions, are horrible for symbolic functions such as the one above. the reason for this is that they don't simplify y[n] after every step and end up with a terrible mess, e.g.

Code (Text):

In[17]:= dr[[6]]
drExp[[6]]
Out[17]= 1/5 (2/3 a (2 a b E^c+1/2 b (2 a E^c+b^2 E^c))+1/4 b (a (2 a E^c+b^2 E^c)+1/3 b (2 a b E^c+1/2 b (2 a E^c+b^2 E^c))))
Out[18]= 1/2 a^2 b E^c+1/6 a b^3 E^c+(b^5 E^c)/120
Finally, compare calculating by taking derivatives with calculating using http://reference.wolfram.com/mathematica/tutorial/FunctionsThatRememberValuesTheyHaveFound.html" [Broken] (in the latter, it's the Expand that makes the difference compared to the DifferenceRoot and RecurrenceTable options given above)

Code (Text):

In[25]:= deriv=Table[1/k! D[Exp[s[x]], {x,k}] /. x->0, {k, 0, 100}]//Expand;//Timing
Out[25]= {0.36, Null}

In[26]:= Clear[y]
y[n_Integer?Positive] := y[n] = 1/n (2 a y[-2+n] + b y[-1+n]) //Expand
y[0] = E^c;
y[1] = b E^c;

In[30]:= mem = Table[y[n], {n, 0, 100}]//Expand;//Timing
Out[30]= {0.3, Null}

In[31]:= deriv==mem
Out[31]= True

Last edited by a moderator: May 5, 2017
3. Aug 13, 2011

### Simon_Tyler

Also, if you want the recurrence relation for the polynomials not evaluated at x=0, try

Code (Text):

s[x_] := a x^2 + b x + c
poly[k_] = Simplify[SeriesCoefficient[Exp[s[x]], {x, x, k}], k > 0]
which returns

Code (Text):
DifferenceRoot[Function[{y,n},
{-2 a y[n] - (b + 2 a x) y[1+n] + (2+n) y[2+n]==0,
y[0]==E^(c+b x+a x^2),  y[1]==(b+2 a x) E^(c+b x+a x^2)}
]][k]

4. Aug 14, 2011

### Simon_Tyler

I must have had my brain turned off this morning, since I didn't mention that what you supplied was an exponential generating function. Just google "generating function recursion relation" to find stuff on how to approach these problems more generally. In particular the free book "generatingfunctionology" is worth reading (not that I can remember much of it)

5. Aug 16, 2011

### mnb96

Thanks a lot for the answers!