Assume, one is given a linear recursion with polynomial coefficients for a sequence $(a_i)_i$, such as

a[i] == i a[i-1]

I would like to convert this recursion into a differential equation for the (formal) generating function $G(x)=\sum_{i}a_i x^i$.

In my example above, direct manipulation leads me to the equation $$ G(x)=2x^2 a_1+x^2 G'(x)+xG(x) $$ which Mathematica can solve in terms of special functions. Interestingly, the power series with coefficients $i!$ (the solution to $a_i=i a_{i-1}$) doesn't converge for any $x$. Since I'm only interested in formal power series, however, I still believe that my question makes sense.

I'm particular interested in a method that also works for higher order recursions with higher degree polynomial coefficients.


I would like to add that I am familiar with the functions GeneratingFunction[] and RSolve[] but neither seems to be of much help for this problem.

Edit 18 April 2014: In version 9 the behaviour of GeneratingFunction appears to have changed in such a way that it doesn't directly solve this problem any longer. I therefore reactivate this question and ask for a way to transform a linear recursion with polynomial coefficients into a differential equation for the generating function.

  • 1
    $\begingroup$ GeneratingFunction[] does what I want. Apparently I wasn't sufficiently familiar with this function after all. $\endgroup$
    – Eckhard
    Aug 26, 2013 at 17:40

2 Answers 2


The usual approach for problems like these is to combine GeneratingFunction[] (as already noted by the OP) with DifferenceRoot[], as in the following:

GeneratingFunction[DifferenceRoot[Function[{a, k}, {a[k] == a[k - 1]/k, a[0] == 1}]][k],
                   k, x]
   DifferentialRoot[Function[{y, x}, {-y[x] + y'[x] == 0, y[0] == 1}]][x]

after which one can use FunctionExpand[]:


However, for the OP's specific case involving the exponential integral, GeneratingFunction[] refuses to budge at all, even after setting VerifyConvergence -> False. One possible recourse is to directly construct the regularized sum:

Sum[DifferenceRoot[Function[{a, k}, {a[k] == k a[k - 1], a[0] == 1}]][k] x^k,
    {k, 0, ∞}, Regularization -> "Borel"]
   -((E^(-1/x) Gamma[0, -(1/x)])/x)

and then you can use DifferentialRootReduce[] to see the associated ODE:

DifferentialRootReduce[%, x]
   DifferentialRoot[Function[{y, x}, {y[x] + (-1 + 3 x) y'[x] + x^2 y''[x] == 0,
                                      y[1] == -(Gamma[0, -1]/E), y'[1] == 1}]][x]

Note that the initial conditions were not placed at $x=0$ even if the function is actually well-defined there.


Beware of "formal power series" of this type. I have an example of a very competent mathematics professional that led him into ?? land. Concerning your question let's look at the condition $a_0=0\cdot a_{-1} $ which is off the end but doesn't matter because it sets $a_0 = 0$

But by recursion this sets all

$a_x = 0 $

So the relation can only be non-trivially true for i>0. BTW: direct calculation doesn't seem to need $2x^2a_1$ In other words.

$G(x)-x^2 G'(x)-xG(x)$ seems to be zero. Gathering the $x^i$ terms $a_i x^i - x^2 x^{i-2} a_{i-1}(i-1) - x^1 x^{i-1} a_{i-1} = a_i x^i - i a_{i-1}x^i $

Which is the requirement. Adding another term makes in non-homogenious. In addition evaluating the left and right hand sides for $x^2$ gives

$a_2 x^2 = 2a_1 x^2 +a_1 x^2 + a_1 x^2$ $a_2 = 4 a_1 $ together with $a_2=2a_1 $ Doesn't work! Unless you start the recursion at i>2.

Sorry for being a little "long winded"


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.