Generating Functions

Only yesterday in a post titled Linear Recurrence Revisited did I look at this recurrence in relation to my diurnal age of 28072:$$ \begin{align} &\text{a}(n+3) = 2\text{a}(n+2) - 3\text{a}(n+1) + 2\text{a}(n)\\ & \text{a}(0) = 1, \text{a}(1) = 1, \text{a}(2) = 0 \end{align}$$Today, a day later, I find that 28073 can also be linked to a linear recurrence via OEIS A157651:$$ \begin{align} &\text{a}(n) = 3\text{a}(n-1) -3\text{a}(n-2) +\text{a}(n-3) \\ &\text{a}(0) =57, \text{a}(1) =308, \text{a}(3)=759 \end{align}$$The initial terms are:

57, 308, 759, 1410, 2261, 3312, 4563, 6014, 7665, 9516, 11567, 13818, 16269, 18920, 21771, 24822, 28073, 31524, 35175, 39026, 43077, 47328, 51779, 56430, 61281, 66332, 71583, 77034, 82685, 88536, 94587, 100838, 107289, 113940, 120791, 127842

However, this time I want to look at ordinary generating functions and exponential generating functions in relation to this linear recurrence. The following was created by Gemini.

$$
\begin{array}{@{}l}



\textbf{Derivation of the Ordinary Generating Function } G(x) \\[12pt]


\text{Let } G(x) = \displaystyle\sum_{n=0}^{\infty} a_n x^n = a_0 + a_1 x + a_2 x^2 + \dots \\[12pt]
\text{The recurrence relation is: } a_n - 3a_{n-1} + 3a_{n-2} - a_{n-3} = 0 \text{ for } n \ge 3. \\[12pt]


\text{We multiply } G(x) \text{ by the polynomial } (1 - 3x + 3x^2 - x^3): \\[6pt]
(1 - 3x + 3x^2 - x^3)G(x) = \\[6pt]
\quad \phantom{-} (a_0 + a_1 x + a_2 x^2 + a_3 x^3 + \dots) \\
\quad - \displaystyle 3x (a_0 + a_1 x + a_2 x^2 + \dots) \\
\quad + \displaystyle 3x^2 (a_0 + a_1 x + \dots) \\
\quad - \displaystyle x^3 (a_0 + \dots) \\[18pt]


\text{Grouping coefficients by powers of } x: \\[6pt]
\text{Constant: } a_0 = 57 \\[6pt]
x^1: a_1 - 3a_0 = 308 - 3(57) = 308 - 171 = 137 \\[6pt]
x^2: a_2 - 3a_1 + 3a_0 = 759 - 3(308) + 3(57) = 759 - 924 + 171 = 6 \\[6pt]
x^n \text{ (for } n \ge 3\text{): } a_n - 3a_{n-1} + 3a_{n-2} - a_{n-3} = 0 \\[18pt]


\text{Thus:} \\[6pt]
(1 - x)^3 G(x) = 57 + 137x + 6x^2 \\[12pt]


\text{The Ordinary Generating Function is:} \\[6pt]
G(x) = \displaystyle\frac{57 + 137x + 6x^2}{(1 - x)^3} \\[24pt]




\textbf{Derivation of the Exponential Generating Function } E(x) \\[12pt]


\text{Let } E(x) = \displaystyle\sum_{n=0}^{\infty} a_n \frac{x^n}{n!}. \\[12pt]
\text{The recurrence } a_n - 3a_{n-1} + 3a_{n-2} - a_{n-3} = 0 \\ \text{ becomes the differential equation:} \\[6pt]
E'''(x) - 3E''(x) + 3E'(x) - E(x) = 0 \\[12pt]


\text{The characteristic equation is } r^3 - 3r^2 + 3r - 1 = (r-1)^3 = 0. \\[6pt]
\text{This has a triple root at } r = 1. \\[12pt]


\text{The general solution for a triple root is:} \\[6pt]
E(x) = (C_1 + C_2 x + C_3 x^2)e^x \\[12pt]


\text{We find constants } C_1, C_2, C_3 \text{ using initial conditions } a_0=57, a_1=308, a_2=759. \\[6pt]
\text{Note that } a_n \text{ is the coefficient of } \frac{x^n}{n!}. \\[6pt]
\text{Expanding } E(x) = (C_1 + C_2 x + C_3 x^2) \displaystyle \sum \frac{x^n}{n!} \\ \text{ allows us to match coefficients directly.} \\[18pt]


\text{Solving for the constants (via derivatives at } x=0 \text{):} \\[6pt]
1. \quad E(0) = a_0 \implies C_1 = 57 \\[6pt]
2. \quad E'(0) = a_1 \implies C_1 + C_2 = 308 \implies 57 + C_2 = 308 \implies C_2 = 251 \\[6pt]
3. \quad E''(0) = a_2 \implies C_1 + 2C_2 + 2C_3 = 759 \\[6pt]
\phantom{3. \quad} 57 + 2(251) + 2C_3 = 759 \\[6pt]
\phantom{3. \quad} 559 + 2C_3 = 759 \implies 2C_3 = 200 \implies C_3 = 100 \\[18pt]


\text{The Exponential Generating Function is:} \\[6pt]
E(x) = (57 + 251x + 100x^2)e^x \\[24pt]




\textbf{Expansion of the Ordinary Generating Function} \\[12pt]


\text{The correct OGF is: } G(x) = \displaystyle\frac{57 + 137x + 6x^2}{(1-x)^3} \\[12pt]


\text{We use the binomial series expansion for } (1-x)^{-3}: \\[6pt]
(1-x)^{-3} = \displaystyle\sum_{n=0}^{\infty} \binom{n+2}{2} x^n = \displaystyle\sum_{n=0}^{\infty} \frac{(n+1)(n+2)}{2} x^n \\[12pt]


\text{Multiply this series by the numerator } (57 + 137x + 6x^2): \\[6pt]
a_n = 57 \cdot \binom{n+2}{2} + 137 \cdot \binom{n+1}{2} + 6 \cdot \binom{n}{2} \\[12pt]


\text{Let's verify the first few terms again with the correct formula:} \\[12pt]


\textbf{For } n=0: \\[6pt]
a_0 = 57 \cdot \frac{(1)(2)}{2} = 57(1) = \mathbf{57} \\[12pt]


\textbf{For } n=1: \\[6pt]
\text{The term } \binom{n}{2} \text{ is 0 for } n < 2. \\[6pt]
a_1 = 57 \cdot \frac{(2)(3)}{2} + 137 \cdot \frac{(1)(2)}{2} \\[6pt]
a_1 = 57(3) + 137(1) = 171 + 137 = \mathbf{308} \\[12pt]


\textbf{For } n=2: \\[6pt]
a_2 = 57 \cdot \frac{(3)(4)}{2} + 137 \cdot \frac{(2)(3)}{2} + 6 \cdot \frac{(2)(1)}{2} \\[6pt]
a_2 = 57(6) + 137(3) + 6(1) \\[6pt]
a_2 = 342 + 411 + 6 = \mathbf{759} \\[12pt]


\textbf{For } n=3: \\[6pt]
a_3 = 57 \cdot \frac{(4)(5)}{2} + 137 \cdot \frac{(3)(4)}{2} + 6 \cdot \frac{(3)(2)}{2} \\[6pt]
a_3 = 57(10) + 137(6) + 6(3) \\[6pt]
a_3 = 570 + 822 + 18 = \mathbf{1410} \\[24pt]




\textbf{Expansion of the Exponential Generating Function} \\[12pt]


\text{We expand } e^x = \displaystyle\sum_{n=0}^{\infty} \dfrac{x^n}{n!} = 1 + x + \dfrac{x^2}{2!} + \dfrac{x^3}{3!} + \dots \\[12pt]


\text{Multiply by the polynomial } (57 + 251x + 100x^2): \\[6pt]
E(x) = 57 \displaystyle\sum \dfrac{x^n}{n!} + 251x \displaystyle\sum \dfrac{x^n}{n!} + 100x^2 \displaystyle\sum \dfrac{x^n}{n!} \\[12pt]


\text{To find } a_n, \text{ we look for the coefficient of } \frac{x^n}{n!} \text{ in the combined sum.} \\[6pt]
\text{Note that } x \cdot \dfrac{x^{n-1}}{(n-1)!} = n \dfrac{x^n}{n!} \text{ and } x^2 \cdot \dfrac{x^{n-2}}{(n-2)!} = n(n-1) \dfrac{x^n}{n!}. \\[12pt]


\text{This gives us the direct formula:} \\[6pt]
a_n = 57(1) + 251(n) + 100(n)(n-1) \\[6pt]
a_n = 100n^2 + 151n + 57 \\[12pt]


\text{Let's calculate the terms using this formula:} \\[12pt]


\textbf{For } n=0: \\[6pt]
a_0 = 0 + 0 + 57 = \mathbf{57} \\[12pt]


\textbf{For } n=1: \\[6pt]
a_1 = 100(1) + 151(1) + 57 = 251 + 57 = \mathbf{308} \\[12pt]


\textbf{For } n=2: \\[6pt]
a_2 = 100(4) + 151(2) + 57 = 400 + 302 + 57 = \mathbf{759} \\[12pt]


\textbf{For } n=3: \\[6pt]
a_3 = 100(9) + 151(3) + 57 = 900 + 453 + 57 = \mathbf{1410}
\end{array}
$$

I prefer to work with the ordinary generating function in which case the following simple code will generate the terms: