[motivation]
## Motivation
The method of [matched asymptotic expansions](/pages/1162) handles problems where the solution has different structure in different *spatial regions* — a [boundary](/page/Boundary) layer here, a smooth outer solution there. But there is a second, fundamentally different way that a regular perturbation expansion can fail: the small perturbation may have a **cumulative effect over long times**, even though the solution looks well-behaved at any fixed instant.
The prototypical setting is a **weakly perturbed oscillator**:
\begin{align*}
\ddot{y} + y + \epsilon\,f(y, \dot{y}) = 0, \qquad \epsilon \ll 1.
\end{align*}
For $\epsilon = 0$, the solution oscillates forever with fixed amplitude and frequency. For small $\epsilon$, the perturbation is weak at any given instant — but over $O(1/\epsilon)$ oscillation periods, its cumulative effect becomes $O(1)$. A regular expansion $y \sim y_0 + \epsilon y_1 + \cdots$ produces terms that grow without bound in time (**secular terms**), causing the expansion to break down at $t = O(1/\epsilon)$.
The resolution is that the solution depends on $t$ and $\epsilon$ through the combination $\epsilon t$ (as well as through $t$ itself) — the secular term is a Taylor expansion of a function of $\epsilon t$ that has been artificially broken apart. The **method of multiple scales** addresses this by introducing $\tau = t$ and $T = \epsilon t$ as independent variables, constructing a Poincaré expansion of a two-variable extension $Y(\tau, T;\epsilon)$ in which $\epsilon t$ is never broken apart, and then evaluating the result on the diagonal $(\tau, T) = (t, \epsilon t)$.
This page develops the method from first principles: the failure of regular perturbation, the two-time-scale extension, the removal of secular terms, and two fully worked examples.
[/motivation]
## The Problem: Weakly Perturbed Oscillators
We consider equations of the form
\begin{align*}
\ddot{y} + y + \epsilon\,f(y, \dot{y}) = 0, \qquad \epsilon \searrow 0,
\end{align*}
where dots denote [derivatives](/page/Derivative) with respect to $t$. We have normalised so that the unperturbed oscillator ($\epsilon = 0$) has angular frequency $\omega_0 = 1$ (period $2\pi$).
A standard Poincaré expansion $y(t;\epsilon) \sim y_0(t) + \epsilon y_1(t) + \cdots$ (as $\epsilon \to 0$ with $t$ fixed) is valid for any fixed $t$, but produces secular terms — terms that grow without bound in $t$ — making the expansion useless for $t = O(1/\epsilon)$. As we will show, these secular terms arise because the true solution depends on $\epsilon t$ as a unit, and expanding in $\epsilon$ with $t$ fixed breaks this combination apart.
The method of multiple scales resolves this as follows:
1. Define a two-variable extension $Y(\tau, T;\epsilon)$ (where $\tau = t$, $T = \epsilon t$) such that $y(t;\epsilon) = Y(t, \epsilon t;\epsilon)$.
2. Find a Poincaré expansion $Y(\tau, T;\epsilon) \sim Y_0(\tau, T) + \epsilon Y_1(\tau, T) + \cdots$ as $\epsilon \to 0$ with $\tau$ and $T$ **both fixed**, subject to the condition that each $Y_k$ is bounded in $\tau$.
3. Evaluate the result on the **physical diagonal** $(\tau, T) = (t, \epsilon t)$ to obtain an approximation to $y(t;\epsilon)$.
The output of step 3 is *not* a Poincaré expansion of $y(t;\epsilon)$ — it is the restriction to a diagonal of a Poincaré expansion of a different function $Y$. The boundedness condition in step 2 is what prevents secular terms from appearing; it determines the previously free functions in the leading-order solution.
## Failure of Regular Perturbation: Secular Terms
[example: Secular Terms in the Duffing Equation]
Consider the Duffing equation
\begin{align*}
\ddot{y} + y + \epsilon y^3 = 0, \qquad \epsilon \searrow 0.
\end{align*}
We attempt a regular expansion $y(t;\epsilon) \sim y_0(t) + \epsilon\,y_1(t)$. Substituting into the ODE:
\begin{align*}
(\ddot{y}_0 + \epsilon\,\ddot{y}_1 + \cdots) + (y_0 + \epsilon\,y_1 + \cdots) + \epsilon(y_0 + \cdots)^3 = 0.
\end{align*}
Expanding $(y_0 + \epsilon y_1 + \cdots)^3 = y_0^3 + O(\epsilon)$ and grouping by powers of $\epsilon$:
\begin{align*}
\underbrace{(\ddot{y}_0 + y_0)}_{O(1)} + \epsilon\underbrace{(\ddot{y}_1 + y_1 + y_0^3)}_{O(\epsilon)} + O(\epsilon^2) = 0.
\end{align*}
*At $O(1)$:* $\ddot{y}_0 + y_0 = 0$, with general solution $y_0 = a\cos(t + \phi_0)$ for constants $a > 0$ and $\phi_0$.
*At $O(\epsilon)$:*
\begin{align*}
\ddot{y}_1 + y_1 = -y_0^3.
\end{align*}
Computing the right-hand side using the identity $\cos^3\theta = \frac{3}{4}\cos\theta + \frac{1}{4}\cos 3\theta$:
\begin{align*}
-y_0^3 = -a^3\cos^3(t + \phi_0) = -\frac{3a^3}{4}\cos(t + \phi_0) - \frac{a^3}{4}\cos 3(t + \phi_0).
\end{align*}
The equation for $y_1$ is now:
\begin{align*}
\ddot{y}_1 + y_1 = -\frac{3a^3}{4}\cos(t + \phi_0) - \frac{a^3}{4}\cos 3(t + \phi_0).
\end{align*}
The $\cos 3(t+\phi_0)$ term on the right is unproblematic — it drives a bounded particular solution proportional to $\cos 3(t+\phi_0)$. But the $\cos(t+\phi_0)$ term is a **resonant forcing**: it solves the homogeneous equation $\ddot{y}_1 + y_1 = 0$. By the theory of forced oscillators, the particular solution for this term is:
\begin{align*}
y_1^{(\text{resonant})} = -\frac{3a^3}{8}\,t\sin(t + \phi_0).
\end{align*}
This grows linearly in $t$. When $t = O(1/\epsilon)$, the "correction" $\epsilon y_1 = O(1)$ — it is no longer small compared to $y_0$. The regular expansion has broken down.
These unbounded terms are called **secular terms** (from the Latin *saeculum*, meaning "age" — they become important over long time scales).
[/example]
[remark: When Do Secular Terms Arise?]
Secular terms appear whenever the right-hand side of the $O(\epsilon)$ equation contains terms that solve the corresponding homogeneous equation. For oscillator problems with natural frequency $\omega_0 = 1$, this means terms proportional to $e^{\pm it}$ (in complex notation) or $\cos t$ and $\sin t$ (in real notation). In Fourier-[series](/page/Series) language: secular terms arise from the **first-harmonic** (fundamental frequency) components of the forcing. Higher harmonics produce bounded responses.
[/remark]
## The Two-Time-Scale Extension
### Motivation
The secular term tells us something physical: the weak perturbation slowly modulates the oscillation over $O(1/\epsilon)$ time periods. The solution involves two processes operating on different scales:
- Fast oscillation on the scale $\tau = t$ (period $2\pi$)
- Slow evolution of amplitude and phase on the scale $T = \epsilon t$
The idea of the method of multiple scales is to treat these as **independent variables** and construct a two-variable function $Y(\tau, T;\epsilon)$ such that
\begin{align*}
y(t;\epsilon) = Y(t,\, \epsilon t;\, \epsilon)
\end{align*}
on the **physical diagonal** $(\tau, T) = (t, \epsilon t)$. This is called a **two-variable extension** of $y$. Crucially, this extension is **not unique**: given $y(t;\epsilon)$, there are many two-variable [functions](/page/Function) $Y(\tau, T;\epsilon)$ that agree with $y$ when restricted to the diagonal $T = \epsilon\tau$. For instance, adding any function that vanishes on the diagonal (such as $(\tau - T/\epsilon) \cdot g(\tau, T)$) produces a different $Y$ with the same $y$. This non-uniqueness is not a defect — it is the **freedom** that the method exploits. The no-secular-terms condition (Step 5 below) will select a specific $Y$ from all the possibilities by requiring that the expansion coefficients $Y_0, Y_1, \ldots$ remain bounded in $\tau$.
### The chain rule
If $\tau = t$ and $T = \epsilon t$, then by the chain rule:
\begin{align*}
\frac{dy}{dt} = \frac{\partial Y}{\partial \tau}\frac{d\tau}{dt} + \frac{\partial Y}{\partial T}\frac{dT}{dt} = \frac{\partial Y}{\partial \tau} + \epsilon\frac{\partial Y}{\partial T}.
\end{align*}
For the second derivative, apply the operator $\frac{d}{dt} = \frac{\partial}{\partial \tau} + \epsilon\frac{\partial}{\partial T}$ twice:
\begin{align*}
\frac{d^2 y}{dt^2} &= \left(\frac{\partial}{\partial \tau} + \epsilon\frac{\partial}{\partial T}\right)\!\left(\frac{\partial Y}{\partial \tau} + \epsilon\frac{\partial Y}{\partial T}\right) \\[4pt]
&= \frac{\partial^2 Y}{\partial \tau^2} + \epsilon\frac{\partial^2 Y}{\partial T\partial\tau} + \epsilon\frac{\partial^2 Y}{\partial \tau\partial T} + \epsilon^2\frac{\partial^2 Y}{\partial T^2} \\[4pt]
&= \frac{\partial^2 Y}{\partial \tau^2} + 2\epsilon\frac{\partial^2 Y}{\partial \tau\partial T} + \epsilon^2\frac{\partial^2 Y}{\partial T^2}.
\end{align*}
These relations hold exactly on the physical diagonal. The key step is to **use them as the governing equations off the diagonal** — that is, to replace $d/dt$ by $\partial/\partial\tau + \epsilon\,\partial/\partial T$ in the original ODE and treat $\tau$ and $T$ as independent.
### The expansion
We then posit the Poincaré expansion
\begin{align*}
Y(\tau, T;\epsilon) \sim Y_0(\tau, T) + \epsilon\,Y_1(\tau, T) + \cdots \quad \text{as } \epsilon \to 0,
\end{align*}
and require that $Y_0$ and $Y_1$ are **bounded** as $\tau \to +\infty$ for each fixed $T$. This is the condition that eliminates secular terms. Since $T = \epsilon t$ is held fixed during the expansion (rather than being broken into $\epsilon \cdot t$), the combination $\epsilon t$ is never separated, and the resulting approximation $y(t;\epsilon) \approx Y_0(t, \epsilon t) + \epsilon\,Y_1(t, \epsilon t)$ formally avoids the secular growth that plagued the regular expansion.
[remark: On Uniform Validity]
It is commonly stated that the multiple-scales approximation is "uniformly valid for $t = O(1/\epsilon)$." Formally, this means $|y(t;\epsilon) - Y_0(t,\epsilon t) - \epsilon Y_1(t,\epsilon t)| \leq C\epsilon$ for $t \in [0, K/\epsilon]$. However, **at the level of this course, this is not proven** — it is a formal expectation based on the absence of secular terms. The method produces an approximation; the rigorous justification of its accuracy over long times is a separate (and harder) problem.
[/remark]
[remark: Periodicity as an Alternative to Boundedness]
An equivalent formulation restricts $\tau$ to one period of the unperturbed oscillator, $\tau \in [0, 2\pi)$, with the physical diagonal $(\tau, T) = (\text{mod}(t, 2\pi), \epsilon t)$. The secular-term condition is then replaced by requiring $Y_0, Y_1, \ldots$ to be **$2\pi$-periodic in $\tau$**. This leads to the same amplitude equations and is sometimes more natural when connecting to Fourier analysis or solvability conditions.
[/remark]
## The General Procedure
For a weakly perturbed oscillator $\ddot{y} + y + \epsilon\,f(y,\dot{y}) = 0$:
**Step 1.** Replace $y(t;\epsilon)$ by $Y(\tau, T;\epsilon)$ using $\dot{y} \to Y_\tau + \epsilon Y_T$ and $\ddot{y} \to Y_{\tau\tau} + 2\epsilon Y_{\tau T} + \epsilon^2 Y_{TT}$.
**Step 2.** Substitute $Y \sim Y_0 + \epsilon Y_1 + \cdots$ and collect by powers of $\epsilon$.
**Step 3.** At $O(1)$: solve $Y_{0\tau\tau} + Y_0 = 0$. The general solution is
\begin{align*}
Y_0(\tau, T) = A(T)\,e^{i\tau} + \overline{A(T)}\,e^{-i\tau},
\end{align*}
where $A(T)$ is an unknown complex amplitude that varies on the slow time scale. (Equivalently, $Y_0 = a(T)\cos(\tau + \theta(T))$ with real amplitude $a$ and phase $\theta$.)
**Step 4.** At $O(\epsilon)$: obtain an equation of the form $Y_{1\tau\tau} + Y_1 = h(\tau, T)$, where $h$ depends on $Y_0$ and hence on $A(T)$.
**Step 5.** Demand that $Y_1$ is bounded (no secular terms). This means: the **resonant part** of $h$ — the components proportional to $e^{i\tau}$ and $e^{-i\tau}$ — must **both** vanish. Since $Y_0$ and the original equation are real-valued, the right-hand side $h(\tau, T)$ is real, and its Fourier decomposition has the form
\begin{align*}
h = \bigl[\text{coefficient}\bigr]\,e^{i\tau} + \overline{\bigl[\text{coefficient}\bigr]}\,e^{-i\tau} + \text{non-resonant terms}.
\end{align*}
That is, the coefficient of $e^{-i\tau}$ is the **complex conjugate** of the coefficient of $e^{i\tau}$. Setting one to zero automatically forces the other to zero. So although *two* resonant terms must vanish, this yields only **one** complex equation (equivalently, two real equations: one for amplitude, one for phase). This is the **amplitude equation** (or **solvability condition**) for $A(T)$.
**Step 6.** Solve the amplitude equation for $A(T)$, then evaluate on the physical diagonal $(\tau, T) = (t, \epsilon t)$ to obtain $y(t;\epsilon)$.
[remark: Phasor Notation]
The complex representation $Y_0 = A(T)e^{i\tau} + \text{c.c.}$ (where c.c. denotes the complex conjugate) is typically the most efficient for calculations. To extract the amplitude and phase, write $A(T) = R(T)\,e^{i\phi(T)}$ with $R > 0$ real. Then $Y_0 = 2R\cos(\tau + \phi)$, so the oscillation amplitude is $2R$ and the slow phase is $\phi(T)$.
[/remark]
## Worked Example I: The Duffing Oscillator
[example: Duffing Oscillator via Multiple Scales]
Consider
\begin{align*}
\ddot{y} + y + \epsilon y^3 = 0, \qquad \epsilon \searrow 0.
\end{align*}
**Step 1: Two-variable extension.**
Replace $y(t;\epsilon)$ by $Y(\tau, T;\epsilon)$. Substituting $\ddot{y} \to Y_{\tau\tau} + 2\epsilon Y_{\tau T} + \epsilon^2 Y_{TT}$ and $y^3 \to Y^3$:
\begin{align*}
Y_{\tau\tau} + 2\epsilon\,Y_{\tau T} + \epsilon^2 Y_{TT} + Y + \epsilon Y^3 = 0.
\end{align*}
Rearranging:
\begin{align*}
Y_{\tau\tau} + Y + \epsilon(2Y_{\tau T} + Y^3) + \epsilon^2 Y_{TT} = 0.
\end{align*}
**Step 2: Expand and collect.**
Substitute $Y \sim Y_0 + \epsilon Y_1 + \cdots$. At $O(1)$, the $\epsilon$ and $\epsilon^2$ terms drop:
\begin{align*}
Y_{0\tau\tau} + Y_0 = 0.
\end{align*}
At $O(\epsilon)$: the $2Y_{\tau T}$ term acts on $Y_0$, and $Y^3 = (Y_0 + \epsilon Y_1 + \cdots)^3 = Y_0^3 + O(\epsilon)$. So:
\begin{align*}
Y_{1\tau\tau} + Y_1 = -2Y_{0\tau T} - Y_0^3.
\end{align*}
**Step 3: Solve at $O(1)$.**
\begin{align*}
Y_0(\tau, T) = A(T)\,e^{i\tau} + \overline{A}(T)\,e^{-i\tau},
\end{align*}
where $A(T)$ is an unknown complex amplitude.
**Step 4: Compute the right-hand side at $O(\epsilon)$.**
We need to compute $-2Y_{0\tau T}$ and $-Y_0^3$ separately.
*The $-2Y_{0\tau T}$ term.* From $Y_0 = Ae^{i\tau} + \bar{A}e^{-i\tau}$:
\begin{align*}
Y_{0\tau} = iA\,e^{i\tau} - i\bar{A}\,e^{-i\tau}, \qquad Y_{0\tau T} = i\frac{dA}{dT}\,e^{i\tau} - i\frac{d\bar{A}}{dT}\,e^{-i\tau}.
\end{align*}
So:
\begin{align*}
-2Y_{0\tau T} = -2i\frac{dA}{dT}\,e^{i\tau} + 2i\frac{d\bar{A}}{dT}\,e^{-i\tau} = -2i\frac{dA}{dT}\,e^{i\tau} + \text{c.c.}
\end{align*}
*The $-Y_0^3$ term.* We need to cube $Y_0 = Ae^{i\tau} + \bar{A}e^{-i\tau}$. Writing $z = Ae^{i\tau}$, so $Y_0 = z + \bar{z}$, and using the identity $(z + \bar{z})^3 = z^3 + 3z^2\bar{z} + 3z\bar{z}^2 + \bar{z}^3 = z^3 + 3|z|^2 z + \text{c.c.}$:
\begin{align*}
Y_0^3 = A^3 e^{3i\tau} + 3A|A|^2 e^{i\tau} + \text{c.c.}
\end{align*}
(Here $|z|^2 = z\bar{z} = |A|^2$ since $|e^{i\tau}|^2 = 1$.) Therefore:
\begin{align*}
-Y_0^3 = -A^3 e^{3i\tau} - 3A|A|^2 e^{i\tau} + \text{c.c.}
\end{align*}
*Combining:*
\begin{align*}
Y_{1\tau\tau} + Y_1 = \left(-2i\frac{dA}{dT} - 3A|A|^2\right)e^{i\tau} - A^3 e^{3i\tau} + \text{c.c.}
\end{align*}
**Step 5: Remove secular terms.**
The operator $\partial_{\tau\tau} + 1$ has homogeneous solutions $e^{\pm i\tau}$. On the right-hand side, the $e^{3i\tau}$ term is non-resonant (it produces a bounded particular solution proportional to $e^{3i\tau}/(1 - 9) = -e^{3i\tau}/8$). But the $e^{i\tau}$ term is **resonant** — it would produce a secular term $\propto \tau\,e^{i\tau}$.
For $Y_1$ to remain bounded, the coefficient of $e^{i\tau}$ must vanish:
\begin{align*}
-2i\frac{dA}{dT} - 3A|A|^2 = 0.
\end{align*}
This is the **amplitude equation**:
\begin{align*}
\frac{dA}{dT} = \frac{3i}{2}A|A|^2.
\end{align*}
**Step 6: Solve the amplitude equation.**
Write $A(T) = R(T)\,e^{i\phi(T)}$ with $R > 0$ real. Then $|A|^2 = R^2$ and:
\begin{align*}
\frac{dA}{dT} = \left(\frac{dR}{dT} + iR\frac{d\phi}{dT}\right)e^{i\phi}.
\end{align*}
Substituting into $dA/dT = \frac{3i}{2}R^2 \cdot Re^{i\phi} = \frac{3i}{2}R^3 e^{i\phi}$ and dividing both sides by $e^{i\phi}$:
\begin{align*}
\frac{dR}{dT} + iR\frac{d\phi}{dT} = \frac{3i}{2}R^3.
\end{align*}
Separating real and imaginary parts:
\begin{align*}
\text{Real:} \quad \frac{dR}{dT} = 0, \qquad \text{Imaginary:} \quad R\frac{d\phi}{dT} = \frac{3}{2}R^3.
\end{align*}
From the real part: $R(T) = \text{const} =: a/2$ (where the factor of $2$ is chosen so that $a$ is the oscillation amplitude, since $Y_0 = 2R\cos(\tau + \phi)$).
From the imaginary part: $d\phi/dT = \frac{3}{2}R^2 = \frac{3a^2}{8}$. Integrating:
\begin{align*}
\phi(T) = \frac{3a^2}{8}T + \phi_0.
\end{align*}
**Step 7: Return to the physical diagonal.**
Set $(\tau, T) = (t, \epsilon t)$:
\begin{align*}
y(t;\epsilon) &\sim Y_0(t, \epsilon t) = 2R\cos(\tau + \phi(T)) = a\cos\!\left(t + \frac{3a^2}{8}\epsilon t + \phi_0\right) \\[4pt]
&= a\cos\!\left[\left(1 + \frac{3\epsilon a^2}{8}\right)t + \phi_0\right].
\end{align*}
**Interpretation.** The amplitude $a$ is constant — the Duffing oscillator does not grow or decay (it is conservative). The frequency, however, is shifted from $1$ to $1 + \frac{3\epsilon a^2}{8}$: the weak cubic nonlinearity produces an **amplitude-dependent frequency shift** of order $\epsilon$. Larger oscillations have a slightly higher frequency.
Note that this result is *not* a Poincaré expansion of $y(t;\epsilon)$ in powers of $\epsilon$ — if you expanded $\cos\!\left[\left(1 + \frac{3\epsilon a^2}{8}\right)t + \phi_0\right]$ for small $\epsilon$ with $t$ fixed, you would recover the secular terms we started with. It is the evaluation on the diagonal of a Poincaré expansion of the two-variable function $Y(\tau, T;\epsilon)$.
[/example]
## Worked Example II: The Van der Pol Oscillator
[example: Van der Pol Oscillator via Multiple Scales]
Consider the Van der Pol equation
\begin{align*}
\ddot{y} + \epsilon\,\dot{y}(y^2 - 1) + y = 0, \qquad \epsilon \searrow 0.
\end{align*}
The perturbation represents **nonlinear damping**: for $|y| > 1$, the $\dot{y}(y^2 - 1)$ term has the same sign as $\dot{y}$, so it acts as positive damping (dissipation). For $|y| < 1$, it has the opposite sign — negative damping (amplification). We therefore expect small oscillations to grow and large oscillations to decay, with a stable **[limit](/page/Limit) cycle** in between.
**Step 1: Two-variable extension.**
Replace $\dot{y} \to Y_\tau + \epsilon Y_T$ and $\ddot{y} \to Y_{\tau\tau} + 2\epsilon Y_{\tau T} + \epsilon^2 Y_{TT}$. The term $\epsilon\dot{y}(y^2 - 1)$ becomes:
\begin{align*}
\epsilon(Y_\tau + \epsilon Y_T)(Y^2 - 1) = \epsilon Y_\tau(Y^2 - 1) + \epsilon^2 Y_T(Y^2 - 1).
\end{align*}
The full equation is:
\begin{align*}
Y_{\tau\tau} + Y + \epsilon\bigl[2Y_{\tau T} + Y_\tau(Y^2 - 1)\bigr] + \epsilon^2\bigl[Y_{TT} + Y_T(Y^2 - 1)\bigr] = 0.
\end{align*}
**Step 2: Expand and collect.**
Substitute $Y \sim Y_0 + \epsilon Y_1$. At $O(1)$:
\begin{align*}
Y_{0\tau\tau} + Y_0 = 0.
\end{align*}
At $O(\epsilon)$: we need $2Y_{0\tau T} + Y_{0\tau}(Y_0^2 - 1)$. The $Y^2$ in the perturbation is $(Y_0 + \epsilon Y_1 + \cdots)^2 = Y_0^2 + O(\epsilon)$:
\begin{align*}
Y_{1\tau\tau} + Y_1 = -2Y_{0\tau T} - Y_{0\tau}(Y_0^2 - 1).
\end{align*}
**Step 3: Solve at $O(1)$.**
For variety, we use the real-valued form:
\begin{align*}
Y_0(\tau, T) = a(T)\cos(\tau + \theta(T)),
\end{align*}
where $a(T) > 0$ is a real amplitude and $\theta(T)$ is a slow phase.
**Step 4: Compute the right-hand side at $O(\epsilon)$.**
*The $-2Y_{0\tau T}$ term.* From $Y_0 = a\cos(\tau + \theta)$:
\begin{align*}
Y_{0\tau} = -a\sin(\tau + \theta).
\end{align*}
Since $a$ and $\theta$ depend on $T$ (not $\tau$):
\begin{align*}
Y_{0\tau T} = \frac{\partial}{\partial T}\bigl[-a\sin(\tau + \theta)\bigr] = -\frac{da}{dT}\sin(\tau + \theta) - a\frac{d\theta}{dT}\cos(\tau + \theta).
\end{align*}
So:
\begin{align*}
-2Y_{0\tau T} = 2\frac{da}{dT}\sin(\tau + \theta) + 2a\frac{d\theta}{dT}\cos(\tau + \theta).
\end{align*}
*The $-Y_{0\tau}(Y_0^2 - 1)$ term.* We have $Y_{0\tau} = -a\sin(\tau + \theta)$ and $Y_0^2 = a^2\cos^2(\tau + \theta)$. So:
\begin{align*}
-Y_{0\tau}(Y_0^2 - 1) = a\sin(\tau + \theta)\bigl[a^2\cos^2(\tau + \theta) - 1\bigr].
\end{align*}
Expanding:
\begin{align*}
= a^3\sin(\tau + \theta)\cos^2(\tau + \theta) - a\sin(\tau + \theta).
\end{align*}
Now use the identity $\sin\alpha\cos^2\alpha = \frac{1}{4}\sin\alpha + \frac{1}{4}\sin 3\alpha$ (which follows from $\cos^2\alpha = \frac{1}{2} + \frac{1}{2}\cos 2\alpha$ and $\sin\alpha\cos 2\alpha = \frac{1}{2}\sin 3\alpha - \frac{1}{2}\sin\alpha$):
\begin{align*}
a^3\sin(\tau+\theta)\cos^2(\tau+\theta) = \frac{a^3}{4}\sin(\tau+\theta) + \frac{a^3}{4}\sin 3(\tau+\theta).
\end{align*}
Therefore:
\begin{align*}
-Y_{0\tau}(Y_0^2 - 1) = \frac{a^3}{4}\sin(\tau+\theta) + \frac{a^3}{4}\sin 3(\tau+\theta) - a\sin(\tau+\theta).
\end{align*}
Collecting the $\sin(\tau+\theta)$ terms:
\begin{align*}
-Y_{0\tau}(Y_0^2 - 1) = \left(\frac{a^3}{4} - a\right)\sin(\tau + \theta) + \frac{a^3}{4}\sin 3(\tau + \theta).
\end{align*}
*Combining both contributions.* The full right-hand side of the $O(\epsilon)$ equation is:
\begin{align*}
\text{RHS} = \left(2\frac{da}{dT} + \frac{a^3}{4} - a\right)\sin(\tau + \theta) + 2a\frac{d\theta}{dT}\cos(\tau + \theta) + \frac{a^3}{4}\sin 3(\tau + \theta).
\end{align*}
**Step 5: Remove secular terms.**
The resonant terms are those proportional to $\sin(\tau + \theta)$ and $\cos(\tau + \theta)$ (these solve the homogeneous equation $Y_{1\tau\tau} + Y_1 = 0$). Setting their coefficients to zero:
\begin{align*}
\text{Coefficient of } \sin(\tau + \theta): \quad 2\frac{da}{dT} + \frac{a^3}{4} - a &= 0, \\[6pt]
\text{Coefficient of } \cos(\tau + \theta): \quad 2a\frac{d\theta}{dT} &= 0.
\end{align*}
From the second equation (with $a > 0$): $d\theta/dT = 0$, so $\theta(T) = \theta_0$ (constant). From the first:
\begin{align*}
\frac{da}{dT} = \frac{1}{2}a\!\left(1 - \frac{a^2}{4}\right).
\end{align*}
These constitute the **amplitude equations** for the Van der Pol oscillator.
**Step 6: Analyse the amplitude equation.**
The equation $da/dT = \frac{1}{2}a(1 - a^2/4)$ has two fixed points:
- $a = 0$: **unstable** (for $a \ll 1$, $da/dT \approx a/2 > 0$, so small perturbations grow).
- $a = 2$: **stable** (linearising $a = 2 + u$ with $|u| \ll 1$: $du/dT = \frac{1}{2}(2+u)(1 - (2+u)^2/4) = \frac{1}{2}(2+u)(1 - 1 - u - u^2/4) \approx \frac{1}{2}(2)(-u) = -u$, so $u \to 0$ exponentially).
This confirms the physical expectation: small oscillations grow, large ones decay, and all solutions approach the **limit cycle** at amplitude $a = 2$.
To solve explicitly, we separate variables:
\begin{align*}
\frac{da}{a(1 - a^2/4)} = \frac{1}{2}\,dT.
\end{align*}
Using partial fractions (writing $\frac{1}{a(1 - a^2/4)} = \frac{1}{a} + \frac{a/4}{1 - a^2/4}$) and integrating:
\begin{align*}
\ln a - \frac{1}{2}\ln\!\left|1 - \frac{a^2}{4}\right| = \frac{T}{2} + C.
\end{align*}
Exponentiating and solving for $a$:
\begin{align*}
a(T) = \frac{2}{\sqrt{1 + b\,e^{-T}}},
\end{align*}
where $b > 0$ is an [integration](/page/Integral) constant determined by the initial conditions. As $T \to \infty$: $a \to 2$ (the limit cycle).
**Step 7: Return to the physical diagonal.**
Setting $T = \epsilon t$:
\begin{align*}
y(t;\epsilon) \sim \frac{2}{\sqrt{1 + b\,e^{-\epsilon t}}}\cos(t + \theta_0).
\end{align*}
**Interpretation.** The phase $\theta_0$ and constant $b$ depend on initial conditions, but at long times ($t \gg 1/\epsilon$) both become irrelevant: $b\,e^{-\epsilon t} \to 0$, and $y \to 2\cos(t + \theta_0)$. The solution is **entrained** to the limit cycle — a periodic oscillation of amplitude $2$ determined by the equation itself, not by initial conditions.
As with the Duffing example, this approximation is not a Poincaré expansion of $y(t;\epsilon)$ — it is the evaluation on the diagonal of a Poincaré expansion of $Y(\tau, T;\epsilon)$.
[/example]
## Comparison with [Matched Asymptotic Expansions](/page/Matched%20Asymptotic%20Expansions)
The method of multiple scales and the method of [matched asymptotic expansions](/pages/1162) address two fundamentally different types of singular perturbation:
**Matched asymptotics** handles problems where the solution has different structure in different *regions of the domain* — an outer region where the regular expansion works, and a thin inner region (boundary layer) where it doesn't. Two separate Poincaré expansions are constructed and stitched together by matching.
**Multiple scales** handles problems where the solution has multiple scales operating *simultaneously at every point* — fast oscillations modulated by slow evolution. A single expansion is constructed using two (or more) independent variables, and secular terms are removed by a solvability condition on the slow scale.
In principle, each method could be applied to the other's domain. Boundary-layer problems can be treated by multiple scales (seeking uniform approximations directly), and oscillator problems can be treated by matched asymptotics (matching an early-time expansion to a late-time one). In practice, matched asymptotics is more natural for boundary layers (because Poincaré expansions are simpler than uniform approximations), while multiple scales is more natural for oscillators (because the two scales coexist rather than being spatially separated).
## References
1. C. M. Bender and S. A. Orszag, *Advanced Mathematical Methods for Scientists and Engineers*, Springer, 1999, Chapter 11.
2. E. J. Hinch, *Perturbation Methods*, Cambridge University Press, 1991, Chapter 6.
3. S. H. Strogatz, *Nonlinear Dynamics and Chaos*, Westview Press, 2015, Chapter 7.