Power-Series Solutions

Every method met so far — integrating factors, separable variables, characteristic equations — is really a lookup table: it works precisely when the answer happens to be built from the functions already in the toolbox (e^x, \sin x, polynomials, logarithms). Most differential equations that arise in real science are not so obliging. Sooner or later you meet an equation whose solution simply isn't any combination of familiar functions — and you need a method that doesn't care.

Here it is: assume the answer is a power series with unknown coefficients, and let the equation tell you what they must be. Write

y = \sum_{n=0}^{\infty} a_n x^n = a_0 + a_1 x + a_2 x^2 + \cdots,

plug it in, and matching the coefficient of each power of x turns the differential equation into a recurrence for the a_n. Solve the recurrence and you have built the solution one coefficient at a time — often recovering a function you already know, but sometimes defining a function nobody had ever written down before. Some of the most important functions in physics — the Bessel functions describing waves on a drum, the Legendre functions describing the shape of a planet's gravity — were born exactly this way: not looked up, but manufactured, term by term, straight out of an equation with nothing else to offer.

Worked example: y' = y rebuilds e^x

We know y' = y has solution y = a_0 e^x. Let us pretend we don't, and discover it from the series.

Step 1 — write the series and its derivative.

y = \sum_{n=0}^{\infty} a_n x^n, \qquad y' = \sum_{n=1}^{\infty} n\,a_n x^{n-1}.

Step 2 — reindex y' so both sides carry x^n. Let m = n - 1 (so n = m+1), then rename m back to n:

y' = \sum_{n=0}^{\infty} (n+1)\,a_{n+1}\, x^n.

Step 3 — impose y' = y. Both are now series in x^n, so match coefficients power by power:

(n+1)\,a_{n+1} = a_n \qquad \text{for every } n \ge 0.

Step 4 — read off the recurrence.

a_{n+1} = \frac{a_n}{n+1}.

Step 5 — unwind it. Starting from a free constant a_0:

a_1 = \frac{a_0}{1}, \quad a_2 = \frac{a_1}{2} = \frac{a_0}{2!}, \quad a_3 = \frac{a_2}{3} = \frac{a_0}{3!}, \quad \dots, \quad a_n = \frac{a_0}{n!}.

Step 6 — sum the series.

y = a_0 \sum_{n=0}^{\infty} \frac{x^n}{n!} = a_0 e^x.

The series is the exponential. The lone free constant a_0 = y(0) is exactly the one integration constant a first-order equation should have.

The same machine on a second-order equation gives two free constants and recovers both trig functions. From y'' + y = 0, matching x^n yields the two-step recurrence

(n+2)(n+1)\,a_{n+2} = -a_n \;\Rightarrow\; a_{n+2} = -\frac{a_n}{(n+1)(n+2)}.

It couples even indices to even and odd to odd. Seeding with a_0 gives a_2 = -a_0/2!,\ a_4 = +a_0/4!, \dots — the cosine series; seeding with a_1 gives the sine series. So

y = a_0 \cos x + a_1 \sin x,

the familiar general solution, assembled purely from a recurrence. The two free seeds a_0 = y(0) and a_1 = y'(0) are the two constants a second-order equation demands.

Second worked example: Airy's equation — a genuinely new function

Not every ODE hands back an old friend. y'' = xy — Airy's equation — looks almost as tame as y' = y, yet no combination of polynomials, exponentials, sines or logarithms solves it. Run the very same power-series machine anyway, and watch it manufacture a brand-new function, coefficient by coefficient.

Step 1 — write the series and its second derivative.

y = \sum_{n=0}^{\infty} a_n x^n, \qquad y'' = \sum_{n=2}^{\infty} n(n-1)\,a_n x^{n-2}.

Step 2 — reindex y'' with m = n-2, then rename m back to n:

y'' = \sum_{n=0}^{\infty} (n+2)(n+1)\,a_{n+2}\, x^n.

Step 3 — write the right-hand side xy, which starts one power of x higher than y itself:

xy = \sum_{n=0}^{\infty} a_n x^{n+1} = \sum_{n=1}^{\infty} a_{n-1}\, x^n.

Step 4 — match coefficients of each power of x. The two sums start at different powers, so the very first one has no partner and must be handled on its own:

n = 0: \quad 2\,a_2 = 0 \;\Rightarrow\; a_2 = 0. \qquad n \ge 1: \quad (n+2)(n+1)\,a_{n+2} = a_{n-1}.

Step 5 — read off the recurrence and unwind it.

a_{n+2} = \frac{a_{n-1}}{(n+1)(n+2)}.

This links every third coefficient: a_3 back to a_0, a_4 back to a_1, and a_5 back to a_2 — but Step 4 already forced a_2 = 0, so the whole chain a_2, a_5, a_8, \dots is zero forever. Two independent series survive: one seeded by the free constant a_0, one by the free constant a_1:

y = a_0\Big(1 + \frac{x^3}{6} + \frac{x^6}{180} + \cdots\Big) + a_1\Big(x + \frac{x^4}{12} + \frac{x^7}{504} + \cdots\Big).

Compute as many terms as you like — neither bracket ever collapses into a disguised polynomial, exponential or trig function. This is the Airy function, defined by its series precisely because no simpler closed form exists. It matters far beyond a classroom curiosity: Airy's equation appears wherever a wave meets a turning point between oscillating and dying away — classically, it produces the fine bright-and-dark fringes of light just inside a rainbow, and the same mathematics governs how a quantum particle behaves near the edge of a region it classically couldn't enter.

Two mix-ups catch almost everyone the first time through this method.

Index-shifting is fiddly bookkeeping — and it's the #1 source of errors. Every sum has to be rewritten so it reads \sum(\cdots)x^n with the same power and the same starting index before you're allowed to match coefficients. Shift by the wrong amount, or forget that the two sums above start at different values of n, and the recurrence you read off is simply wrong. The safest habit: after reindexing, write out the first two or three terms of each sum longhand and check they really do line up power by power before trusting the general formula.
The recurrence never invents a_0 out of thin air. It only tells you how each later coefficient depends on earlier ones — the first one or two coefficients (a_0 alone for a first-order equation, a_0 and a_1 for a second-order one) are free choices that the method simply cannot make for you. That is exactly where an initial condition such as y(0) or y'(0) enters: it tells you which free constants to plug in before you start unwinding the recurrence.

Bessel's equation is born from exactly this kind of series bookkeeping, and its solutions — the Bessel functions — lead a bizarre double life. Friedrich Bessel first wrote them down in the 1820s while cataloguing the wobbles in planetary orbits caused by planets tugging on one another, a purely astronomical problem with not a drum in sight.

Yet the very same functions turn out to be the natural building blocks for a vibrating circular drumhead. A guitar string is straight, so its modes are plain sines; a drum is round, so separating its wave equation in circular coordinates turns the radial part into Bessel's equation, and the drum's rings-within-rings vibration patterns are literally graphs of Bessel functions. One power series, discovered while tracking planets across the sky, ended up describing every timpani, snare and kettle drum ever built.

When the method is guaranteed

Consider y'' + p(x)\,y' + q(x)\,y = 0 near a point x_0.

Ordinary point (p, q analytic at x_0): there is a power-series solution y = \sum a_n (x - x_0)^n with two free constants a_0, a_1, and it converges at least as far as the nearest singularity of p or q.
Regular singular point ((x-x_0)p and (x-x_0)^2 q analytic): the Frobenius method applies — seek y = (x-x_0)^r \sum a_n (x-x_0)^n, where the exponent r solves the indicial equation.
Matching like powers turns the equation into a recurrence for the a_n; the number of free constants equals the order.

The dividing line is whether the coefficient functions p, q stay finite. At an ordinary point a plain power series works. At a singular point — where a coefficient blows up — the naive series can fail, and Frobenius's (x-x_0)^r prefactor rescues it.

This is where mathematics' celebrated special functions are born — each the series solution of an equation that has no elementary answer:

Airy y'' = xy — the turning point between oscillation and decay (rainbows, quantum tunnelling).
Bessel x^2 y'' + x y' + (x^2 - \nu^2)y = 0 — vibrating drums and circular waveguides; a regular singular point at the origin.
Legendre (1 - x^2)y'' - 2x y' + \ell(\ell+1)y = 0 — the angular part of problems on the sphere.

None has a closed form in elementary functions, yet the recurrence delivers every coefficient. The power-series method is how these functions were first defined and how their values are still computed.

Watch the partial sums converge

Here is the series solution of y' = y with a_0 = 1 — the partial sum \sum_{n=0}^{N} x^n/n! against the true e^x (dashed). Add terms with the slider and watch the polynomial hug the exponential over a wider and wider window.

Notice that this particular series never gives up: however far you push x from the origin, adding enough terms eventually catches up with e^x, because e^x has no bad point anywhere to stop it. That is a special luxury. The Airy and Bessel series above are guaranteed by the same theorem to converge only out to the nearest place where the equation's coefficients misbehave — for Airy's equation, which has no singular point at all, that radius is infinite too, but for an equation like (1-x)y' = y the series would need to stop trusting itself at x = 1, exactly where the coefficient blows up.