# Integral of GeneratingFunction

I know that GeneratingFunction can be used to compute the generating function $\sum_{n=0}^\infty a_n x^n$ of a sequence $(a_n)_n$ via

GeneratingFunction[a[n],n,x]


I also know that polynomial (in $n$) coefficients of the general term $a_n$ of the sequence result in an expression involving derivatives of the generating function, e.g.

GeneratingFunction[n a[n],n,x]
(* x (GeneratingFunction^(0,0,1))[a[n],n,x] *)


This does not seem to happen with rational (in $n$) coefficients, which would lead to indefinite integrals of the generating function. For example I would expect the input

GeneratingFunction[a[n]/(n + 1), n, x]


to result in

(*1/x Integrate[GeneratingFunction[a[n],n,y],{y,0,x}]*)


or

(*1/x GeneratingFunction^(0,0,-1))[a[n],n,x]*)


but this does not seem to be the case.

I thus have the following question: Is it possible to get Mathematica to compute the generating function of a linear combination of lagged values of $a_n$ with rational (in $n$) coefficients in terms of derivatives and integrals of the generating function of the sequence $(a_n)_n$?

-
GeneratingFunction[a[n]/(n + 1), n, x] === 1/x Integrate[GeneratingFunction[a[n], n, y], {y, 0, x}] is False. Are you sure that it lands there ? – Rorschach Sep 2 '13 at 18:18
Thanks for your comment, @Blackbird. That's what I obtain through fairly easy manual rearrangement of the series defining the generating function. The first example that came to my mind, $a_n=(n+1)^{-k}$ supports the claim. The question is, how do I get Mathematica to recognise the equality and, ideally, to simplify the LHS of your equation accordingly into the RHS. – Eckhard Sep 2 '13 at 18:28
GeneratingFunction^(0,0,1) is not recognized by M, what do you mean by this ? – Rorschach Sep 2 '13 at 19:11
I mean the derivative with respect to the third argument, $x$ in this case. – Eckhard Sep 2 '13 at 19:57
For me GeneratingFunction[n a[n],n,x] evaluates to x GeneratingFunction[(1+n)a[1+n], n, x]. Which version of Mma are you using? (I'm using 9.) – Andrew MacFie Sep 18 '13 at 0:08

In this answer, I handle coefficients of the form

$$\frac{p(n)}{q(n)},$$

where $p(n)$ and $q(n)$ are polynomials in the index symbol $n$, $\deg(p(n)) < \deg(q(n))$, and all of the roots of $q(n)$ are negative integers. If $\deg(p(n)) \geq \deg(q(n))$, we may have to differentiate a power series, and as I mentioned in my comment, Mathematica 9 doesn't handle that conveniently. If $q(n)$ has positive integer roots, we have to handle undefined expressions and things are more complicated.

First, some setup code

(* http://mathematica.stackexchange.com/questions/32531/how-to-unprotect-generatingfunction *)
GeneratingFunction;
protected = Unprotect[GeneratingFunction];

(* what can usually be used as a variable, according to
ref/message/General/ivar *)
variablePattern = Except[_String | _?NumberQ | _Plus | _Times |
_Sum | _Product | _^_Integer];

arg2 = Sequence[n:variablePattern, x:variablePattern];


Now, we add a downvalue for GeneratingFunction which expands $p(n)/q(n)$ into partial fractions.

GeneratingFunction[expr_, arg2,
opts:OptionsPattern[]] /; expr =!= Apart[expr, n] :=
GeneratingFunction[Apart[expr, n], n, x, opts];


I will assume that the coefficients are now in the form

$$\frac{c}{(n + k)^j}.$$

To go from here, in the next rule we add a downvalue saying

$$\sum_{n \geq 0} \frac{1}{n+k} f(n) x^n = x^{-k} \int_0^x t^{k-1} \left( \sum_{n \geq 0} f(n) t^n \right) dt,$$

which gets used repeatedly until $j=0$.

GeneratingFunction[expr_ * (n_ + k_?Positive)^(j_?Negative), arg2, opts:OptionsPattern[]] := With[
{t = Unique[]},
x^(-k) * Integrate[t^(k-1) * GeneratingFunction[expr * (n + k)^(j + 1),
n, t, opts], {t, 0, x}]
];


And we're done.

Protect[Evaluate[protected]];

-