[motivation] ## Motivation Many differential equations arising in physics and engineering contain a small dimensionless parameter $\epsilon$ — the ratio of two disparate scales. It might be the ratio of viscous to inertial forces in a fluid, the ratio of a reaction rate to a diffusion rate, or the ratio of a membrane thickness to its lateral extent. When $\epsilon$ is small, we seek an *asymptotic approximation*: an expression that becomes increasingly accurate as $\epsilon \to 0$, even though we cannot solve the equation exactly. The most natural first attempt is to expand the solution in a [power series](/page/Power%20Series) in $\epsilon$ and substitute into the equation. This is called a **regular perturbation expansion**, and it works when the character of the equation does not change as $\epsilon \to 0$. But it can fail catastrophically when $\epsilon$ multiplies the highest-order derivative. Setting $\epsilon = 0$ then *reduces the order* of the equation, and the resulting lower-order equation has fewer free constants than [boundary](/page/Boundary) or initial conditions to satisfy. No single expansion in $\epsilon$ is valid across the entire domain. The resolution is that the solution develops structure on two different scales: a "slow" outer scale where the regular expansion works, and a "fast" inner scale — a thin region of width $\delta(\epsilon) \to 0$ — where the solution varies rapidly. The method of **matched asymptotic expansions** constructs separate Poincaré expansions on each scale and stitches them together by demanding consistency in the overlap region where both are valid. This page develops the method from first principles: the formal setup, the rescaling procedure, dominant balance, matching, and composite approximations. [/motivation] ## The Problem We consider a differential equation of the form \begin{align*} F\!\left(x, y, y', y'', \ldots; \epsilon\right) = 0, \qquad x \in \Omega \subseteq \mathbb{R}, \end{align*} subject to boundary or initial conditions, in the [limit](/page/Limit) $\epsilon \searrow 0$. Here $y = y(x;\epsilon)$ depends on the independent variable $x$ and parametrically on $\epsilon$. Our goal is to find a **Poincaré asymptotic expansion** of $y$ — that is, an expression of the form \begin{align*} y(x;\epsilon) \sim \sum_{k=0}^{N} \mu_k(\epsilon) \, y_k(x) \quad \text{as } \epsilon \to 0, \end{align*} where $\{\mu_k(\epsilon)\}$ is an asymptotic [sequence](/page/Sequence) (i.e. $\mu_{k+1}(\epsilon) = o(\mu_k(\epsilon))$ as $\epsilon \to 0$) called the **gauge [functions](/page/Function)**, and each $y_k(x)$ is independent of $\epsilon$. The gauge functions are not known in advance — they are determined by the structure of the equation itself. In a singular perturbation problem, no single such expansion is valid for all $x \in \Omega$. Instead, we will construct *two* Poincaré expansions — one valid in an "outer" region and one in an "inner" region — and match them in their common domain of validity. ## Failure of Regular Perturbation Theory [example: Failure of Regular Perturbation — A Singularly Perturbed IVP] Consider the initial-value problem \begin{align*} \epsilon \frac{d^2 y}{dt^2} + \frac{dy}{dt} + y = 0, \qquad y(0) = 0, \quad y'(0) = 1, \end{align*} in the limit $\epsilon \searrow 0$. We attempt a regular perturbation expansion \begin{align*} y(t;\epsilon) \sim y_0(t) + \epsilon \, y_1(t) + \cdots \quad \text{as } \epsilon \to 0, \end{align*} with $t$ fixed. Substituting $y \sim y_0 + \epsilon\,y_1 + \cdots$ into $\epsilon y'' + y' + y = 0$: \begin{align*} \epsilon(y_0'' + \epsilon\,y_1'' + \cdots) + (y_0' + \epsilon\,y_1' + \cdots) + (y_0 + \epsilon\,y_1 + \cdots) = 0. \end{align*} Expanding and grouping by powers of $\epsilon$: \begin{align*} \underbrace{(y_0' + y_0)}_{O(1)} + \epsilon\underbrace{(y_0'' + y_1' + y_1)}_{O(\epsilon)} + \cdots = 0. \end{align*} At $O(1)$, the $\epsilon y_0''$ term is absorbed into the $O(\epsilon)$ group, leaving: \begin{align*} \frac{dy_0}{dt} + y_0 = 0, \end{align*} a first-order equation with general solution $y_0 = A e^{-t}$. There is one free constant and two initial conditions — $y(0) = 0$ forces $A = 0$, while $y'(0) = 1$ forces $y_0'(0) = -A = 1$, so $A = -1$. Both cannot hold simultaneously. The cause is clear: $\epsilon$ multiplies the highest derivative $y''$. Setting $\epsilon = 0$ reduces the equation from second order (two-parameter family of solutions) to first order (one-parameter family). The expansion cannot hold uniformly on any interval containing $t = 0$. [/example] [remark: Signals of Singular Perturbation] The most common signal of a singular perturbation is that the small parameter multiplies the highest-order derivative. But more generally, any situation where setting $\epsilon = 0$ changes the *qualitative character* of the equation — its order, its type (elliptic to parabolic), or the nature of its singularities — should be suspected as singular. The failure of the regular expansion is the diagnostic, not the cause. [/remark] ## The Inner Region: Rescaling and Dominant Balance The failure of the regular expansion tells us that the solution must vary rapidly in some narrow region — a **boundary layer** or **initial layer** — where the highest-order derivative cannot be neglected. The goal of this section is to determine the width of this layer rigorously. ### The rescaling We make the following **structural assumption**: there exists a point $\bar{x} \in \overline{\Omega}$ and a function $\delta(\epsilon) \to 0$ as $\epsilon \to 0$, such that the solution varies on the scale $\delta$ near $\bar{x}$. We introduce the rescaled variable and function \begin{align*} \xi := \frac{x - \bar{x}}{\delta(\epsilon)}, \qquad Y(\xi;\epsilon) := y\!\left(\bar{x} + \delta(\epsilon)\,\xi;\,\epsilon\right). \end{align*} This is an exact change of variables. By the chain rule: \begin{align*} \frac{dy}{dx} = \frac{1}{\delta}\frac{dY}{d\xi}, \qquad \frac{d^2 y}{dx^2} = \frac{1}{\delta^2}\frac{d^2 Y}{d\xi^2}. \end{align*} We then assume that $Y(\xi;\epsilon)$ admits a Poincaré expansion \begin{align*} Y(\xi;\epsilon) \sim \mu_0(\epsilon)\,Y_0(\xi) + \mu_1(\epsilon)\,Y_1(\xi) + \cdots \quad \text{as } \epsilon \to 0, \; \xi \text{ fixed}, \end{align*} where each $Y_k(\xi)$ is independent of $\epsilon$. In particular, $Y_0(\xi)$, $Y_0'(\xi)$, $Y_0''(\xi)$ are all $O(1)$ as $\epsilon \to 0$. This is the precise content of "the rescaling has captured the correct scale of variation." ### Dominant balance After substituting $y(x;\epsilon) = Y(\xi;\epsilon)$ into the ODE using the chain rule, each term in the equation acquires a coefficient involving $\epsilon$ and $\delta$. Since $Y$ and its $\xi$-[derivatives](/page/Derivative) are $O(1)$ by assumption, the relative sizes of the terms are determined entirely by these coefficients. The function $\delta(\epsilon)$ is determined by **dominant balance**: we require that the leading-order equation for $Y_0$ is *nontrivial* — that is, at least two terms survive at leading order. This means the two largest coefficients must be of the same asymptotic order. The procedure is: 1. List all possible pairings of two or more terms. 2. For each pairing, solve for $\delta(\epsilon)$ by equating the corresponding coefficients. 3. Check consistency: every term *not* in the pairing must have a coefficient that is asymptotically smaller than the balanced terms. 4. Discard any pairing that leads to a contradiction (i.e. an "unbalanced" term that is larger than the supposed dominant pair). Typically only one pairing is consistent; this determines $\delta(\epsilon)$. [example: Dominant Balance for the Singularly Perturbed IVP] Returning to $\epsilon y'' + y' + y = 0$ with $\xi = (t - \bar{t})/\delta$ and $Y(\xi;\epsilon) = y(\bar{t} + \delta\xi;\epsilon)$, the transformed equation is \begin{align*} \frac{\epsilon}{\delta^2}\frac{d^2 Y}{d\xi^2} + \frac{1}{\delta}\frac{dY}{d\xi} + Y = 0. \end{align*} The three terms have coefficients $\epsilon/\delta^2$, $1/\delta$, and $1$, respectively. **Balance $\epsilon/\delta^2 \sim 1/\delta$:** gives $\delta \sim \epsilon$. Then both coefficients are $O(1/\epsilon)$, while the third term is $O(1) \ll O(1/\epsilon)$. ✓ Consistent. **Balance $\epsilon/\delta^2 \sim 1$:** gives $\delta \sim \epsilon^{1/2}$. Then $1/\delta \sim \epsilon^{-1/2}$, which is larger than both of the supposedly balanced terms. ✗ The $Y'$ term dominates with nothing to balance it. **Balance $1/\delta \sim 1$:** gives $\delta \sim 1$. Then $\epsilon/\delta^2 \sim \epsilon \ll 1$, so the $Y''$ term is negligible — we recover the regular expansion. ✗ Contradicts the assumption that $Y''$ participates. Only the first balance is consistent: $\delta(\epsilon) \sim \epsilon$. Setting $\delta = \epsilon$ in the transformed equation: \begin{align*} \frac{\epsilon}{\epsilon^2}\frac{d^2 Y}{d\xi^2} + \frac{1}{\epsilon}\frac{dY}{d\xi} + Y = 0 \quad \Longrightarrow \quad \frac{1}{\epsilon}\frac{d^2 Y}{d\xi^2} + \frac{1}{\epsilon}\frac{dY}{d\xi} + Y = 0. \end{align*} Multiplying through by $\epsilon$: \begin{align*} \frac{d^2 Y}{d\xi^2} + \frac{dY}{d\xi} + \epsilon Y = 0, \end{align*} where $\epsilon$ has migrated from the highest derivative to the lowest-order term — confirming that the rescaling has "unfolded" the layer. [/example] ## Determining the Layer Location In an initial-value problem where all conditions are imposed at one point, the layer location is immediate. In a boundary-value problem with conditions at both endpoints, the layer could in principle be at either end (or at an interior point). Its location must be determined by a consistency requirement that comes from **matching** — the procedure developed in the next section. We preview the key consequence here, since it is needed before the full matching machinery is assembled. The matching procedure (whether by intermediate variables or Van Dyke's principle) requires us to compare the inner and outer expansions in an overlap region. In the intermediate variable framework, this means introducing $\eta = (x - \bar{x})/\delta(\epsilon)^\alpha$ with $0 < \alpha < 1$ and $\eta$ fixed, and then **composing** the inner solution $Y_0$ with the coordinate change $\xi = \delta(\epsilon)^{\alpha - 1}\eta$ — note that $\xi$ and $\eta$ are both independent coordinates related by a rescaling that depends on $\epsilon$. The matching step then requires expanding the resulting expression $Y_0(\delta(\epsilon)^{\alpha-1}\eta)$ as $\epsilon \to 0$ with $\eta$ fixed. Since $\alpha < 1$ and $\delta(\epsilon) \to 0$, the prefactor $\delta(\epsilon)^{\alpha-1} \to +\infty$ as $\epsilon \to 0$. Therefore, even though $\eta$ is held fixed, the **argument** of $Y_0$ grows without bound during this $\epsilon \to 0$ limit. Crucially, we do **not** need this limit to be finite for all $\eta$ — only for $\eta$ in the **overlap region**, which has a definite sign determined by the geometry. If the outer region lies to the right of $\bar{x}$, then $\eta > 0$ in the overlap, and matching requires $Y_0(s) \to \text{finite}$ as $s \to +\infty$. If the outer region lies to the left, then $\eta < 0$ in the overlap, and matching requires $Y_0(s) \to \text{finite}$ as $s \to -\infty$. The behaviour of $Y_0$ for the opposite sign of $\eta$ — which corresponds to the boundary side, not the outer region — is irrelevant to matching. This one-sided finiteness condition selects the layer location $\bar{x}$. For a given inner solution $Y_0$, only those $\bar{x}$ survive for which $Y_0(s)$ remains bounded as $s \to +\infty$ or $s \to -\infty$ in the direction of the outer region. [remark: Boundedness Must Hold at Every Order] The boundedness condition applies not only to $Y_0$ but to **every** $Y_n$ in the inner expansion. The reason is the integrity of the Poincaré expansion itself. The inner expansion \begin{align*} Y(\xi;\epsilon) \sim Y_0(\xi) + \epsilon\,Y_1(\xi) + \epsilon^2\,Y_2(\xi) + \cdots \end{align*} is a Poincaré expansion only if each successive term is asymptotically smaller than the previous one. In the matching step, we compose with the outer coordinate by evaluating at $\xi = (x - \bar{x})/\delta$, where $\delta \to 0$, and re-expand as $\epsilon \to 0$ with $x$ fixed. The $n$-th term becomes $\epsilon^n Y_n(x/\delta)$, and for this to remain $O(\epsilon^n)$, the function $Y_n$ must not blow up at the relevant end (toward the outer region). For example, if $Y_2(\xi) \sim c\,e^{\xi}$ as $\xi \to +\infty$ and the outer region lies to the right, then $\epsilon^2 Y_2(x/\epsilon) \sim \epsilon^2 e^{x/\epsilon}$, which is **beyond all algebraic orders** — it overwhelms every term in the outer expansion, and matching is impossible. By contrast, polynomial growth of $Y_n(\xi)$ is harmless: if $Y_1(\xi) \sim c\xi$ as $\xi \to \infty$, then $\epsilon Y_1(x/\epsilon) \sim cx = O(1)$, which simply migrates from $O(\epsilon)$ into $O(1)$ during the re-expansion step and gets absorbed into the matching at that order. In practice, since the inner equations at each order are typically the same homogeneous operator (with different inhomogeneous terms), the growing and decaying exponentials are the same at every order. The exponential that was discarded at $O(1)$ must be discarded at every subsequent order for the same reason. [/remark] [example:Layer Location Depends On The Equation] For the equation $\epsilon y'' + y' = x$ (positive coefficient on $y'$), the inner equation $Y_0'' + Y_0' = 0$ has solutions $1$ and $e^{-\xi}$. The exponential $e^{-\xi}$ is bounded as $\xi \to +\infty$ but unbounded as $\xi \to -\infty$. So matching works when the outer region is to the right ($\eta > 0$), forcing the layer to the **left** endpoint. For the equation $\epsilon y'' - y' = x$ (negative coefficient on $y'$), the inner equation $Y_0'' - Y_0' = 0$ has solutions $1$ and $e^{+\xi}$. The exponential $e^{+\xi}$ is bounded as $\xi \to -\infty$ but unbounded as $\xi \to +\infty$. So matching works when the outer region is to the left ($\eta < 0$), forcing the layer to the **right** endpoint. The sign of the first-order coefficient determines which exponential appears in the inner solution, which determines the direction in which the composition stays bounded, which determines the layer location. This is the precise mechanism behind the physical intuition that the layer forms "where the flow pushes against the boundary." [/example] [remark: No Coordinate Is Sent To Infinity] A common source of confusion: the condition "the inner solution must be bounded as $\xi \to +\infty$" does **not** mean we are sending the inner coordinate $\xi$ to infinity inside the inner expansion. The inner expansion $Y(\xi;\epsilon) \sim Y_0(\xi) + \cdots$ is valid for $\xi$ fixed as $\epsilon \to 0$, and $\xi$ is simply an independent coordinate. The large-argument behaviour of $Y_0$ enters at a different point: when we **compose** $Y_0$ with the coordinate change $\xi = \delta^{\alpha-1}\eta$ and take $\epsilon \to 0$. It is the $\epsilon \to 0$ limit (with $\eta$ fixed) that forces the argument of $Y_0$ to grow, because $\delta^{\alpha-1} \to +\infty$. The shorthand "$Y_0$ must decay as $\xi \to +\infty$" is a statement about the function $Y_0$ evaluated at large arguments — which is what the matching limit requires — not a statement about the inner coordinate itself becoming infinite. [/remark] [example: Layer Location for a Boundary-Value Problem] Consider $\epsilon y'' + y' = x$, $y(0) = 0$, $y(1) = 1$. The rescaled leading-order inner equation is $Y_0'' + Y_0' = 0$, with general solution \begin{align*} Y_0(\xi) = c_1 + c_2 \, e^{-\xi}. \end{align*} Now test each candidate layer location by checking whether matching is consistent. **Case $\bar{x} = 0$.** Here $\xi = x/\epsilon$, and the outer region lies to the right ($x > 0$). In the matching step, we introduce $\eta = x/\epsilon^\alpha$ with $0 < \alpha < 1$ and evaluate $Y_0$ at $\xi = \epsilon^{\alpha-1}\eta$: \begin{align*} Y_0(\epsilon^{\alpha-1}\eta) = c_1 + c_2 \, e^{-\epsilon^{\alpha-1}\eta}. \end{align*} For $\eta > 0$ (the outer side), we take $\epsilon \to 0$. Since $\alpha - 1 < 0$, the exponent $\epsilon^{\alpha-1} \to +\infty$, so $-\epsilon^{\alpha-1}\eta \to -\infty$, and $e^{-\epsilon^{\alpha-1}\eta} \to 0$. Therefore: \begin{align*} Y_0(\epsilon^{\alpha-1}\eta) \to c_1 \quad \text{as } \epsilon \to 0 \; (\eta > 0 \text{ fixed}). \end{align*} This is finite — the composition has a well-defined $\epsilon \to 0$ limit, so matching can proceed. The outer solution will determine $c_1$. **Case $\bar{x} = 1$.** Here $\xi = (x-1)/\epsilon$, and the outer region lies to the left ($x < 1$). In the intermediate variable, $\eta = (x-1)/\epsilon^\alpha < 0$ for $x < 1$. We evaluate $Y_0$ at $\xi = \epsilon^{\alpha-1}\eta$ with $\eta < 0$: \begin{align*} Y_0(\epsilon^{\alpha-1}\eta) = c_1 + c_2 \, e^{-\epsilon^{\alpha-1}\eta}. \end{align*} For $\eta < 0$, the exponent $-\epsilon^{\alpha-1}\eta = \epsilon^{\alpha-1}|\eta| \to +\infty$ as $\epsilon \to 0$. Therefore: \begin{align*} Y_0(\epsilon^{\alpha-1}\eta) = c_1 + c_2 \, e^{\epsilon^{\alpha-1}|\eta|} \to \begin{cases} +\infty & \text{if } c_2 \neq 0, \\ c_1 & \text{if } c_2 = 0. \end{cases} \end{align*} If $c_2 \neq 0$, the composition diverges and matching fails — the inner expansion cannot agree with the bounded outer solution. If $c_2 = 0$, the inner solution is the constant $c_1$, which has no boundary-layer structure and cannot satisfy the boundary condition at $x = 1$ nontrivially. Therefore $\bar{x} = 0$: the boundary layer is at the left endpoint. [/example] [remark: Boundary Layers at Infinity] The discussion above assumes the domain is a finite interval with the layer at one endpoint. But boundary layers can also occur **at infinity**: the outer solution (obtained by setting $\epsilon = 0$) satisfies the condition at the finite boundary but fails to meet the far-field condition as $x \to \infty$. There is no new machinery required. The observation is that **a boundary layer at infinity becomes a standard boundary layer at the origin after inversion**. If $y(x)$ is defined on $[1, \infty)$ with a layer as $x \to \infty$, define $t = 1/x$ and $H(t) = y(1/t)$. The domain maps to $t \in (0, 1]$, and the far-field condition at $x \to \infty$ becomes a condition at $t = 0$. The problem for $H(t)$ is a standard singularly perturbed problem with a boundary layer at $t = 0$, and the entire procedure — dominant balance, inner expansion, matching — applies verbatim. If dominant balance in the inverted problem yields a layer width $\delta(\epsilon)$ at $t = 0$, the inner variable is $\tau = t/\delta$. Translating back to the original coordinate: $\tau = 1/(\delta x)$, so the inner region corresponds to $x = O(1/\delta)$. The natural inner variable in the original coordinate is therefore $X = \delta x$, not $x/\delta$ — one **multiplies** by the small scale rather than dividing, because one is "zooming out" to capture a distant region rather than "zooming in" to resolve a thin one. For instance, if $\delta = \epsilon$, then $X = \epsilon x$. When $x$ is $O(1/\epsilon)$, $X$ is $O(1)$, and the rescaled equation typically promotes a previously negligible $O(\epsilon^2)$ term to leading order — this term encodes the slow decay mechanism (exponential, algebraic, or oscillatory) that the outer expansion cannot capture. [/remark] ## The Outer and Inner Expansions We are now in a position to define the two expansions precisely. [definition: Outer Expansion] The **outer expansion** is the Poincaré expansion of $y(x;\epsilon)$ as $\epsilon \to 0$ with $x$ held fixed (away from $\bar{x}$): \begin{align*} y(x;\epsilon) \sim \sum_{k=0}^{N} \mu_k^{\mathrm{out}}(\epsilon)\,y_k(x) \quad \text{as } \epsilon \to 0, \; x \text{ fixed}, \end{align*} where $\{\mu_k^{\mathrm{out}}(\epsilon)\}$ is an asymptotic sequence and each $y_k$ is $\epsilon$-independent. The gauge functions and the $y_k$ are determined by substituting the expansion into the original ODE, collecting terms at each asymptotic order, and solving the resulting hierarchy of equations. The outer expansion inherits the boundary condition(s) that lie *outside* the layer. [/definition] [definition: Inner Expansion] Let $\delta(\epsilon) \to 0$ be the layer width determined by dominant balance, $\bar{x}$ the layer location, $\xi = (x - \bar{x})/\delta(\epsilon)$, and $Y(\xi;\epsilon) = y(\bar{x} + \delta\xi;\epsilon)$ the rescaled function. The **inner expansion** is the Poincaré expansion of $Y(\xi;\epsilon)$ as $\epsilon \to 0$ with $\xi$ held fixed: \begin{align*} Y(\xi;\epsilon) \sim \sum_{k=0}^{N} \mu_k^{\mathrm{in}}(\epsilon)\,Y_k(\xi) \quad \text{as } \epsilon \to 0, \; \xi \text{ fixed}, \end{align*} where $\{\mu_k^{\mathrm{in}}(\epsilon)\}$ is an asymptotic sequence and each $Y_k$ is $\epsilon$-independent. The gauge functions are determined by substituting into the transformed ODE and collecting terms. The inner expansion inherits the boundary or initial condition(s) at $\bar{x}$. [/definition] [remark: Gauge Functions Are Problem-Dependent] The gauge functions $\mu_k^{\mathrm{out}}$ and $\mu_k^{\mathrm{in}}$ are not known a priori. They emerge from the structure of the equation and its boundary data. In many problems the gauge functions are integer powers of $\epsilon$ (i.e. $1, \epsilon, \epsilon^2, \ldots$), but they need not be. For example, initial conditions of size $\epsilon$ might force the inner expansion to begin at $O(\epsilon)$ rather than $O(1)$. In nonlinear problems, half-integer powers $\epsilon^{1/2}$ can appear (see the nonlinear example below). In problems with degenerate operators, [logarithmic gauge functions](/pages/1138) such as $1/\ln\epsilon$ arise. [/remark] [remark: [Distribution](/page/Distribution) of Boundary Conditions] In a typical singularly perturbed two-point boundary-value problem, one boundary condition is imposed on the inner expansion (the condition at $\bar{x}$) and the other on the outer expansion (the condition at the far boundary). Each expansion is thereby partially determined, but free constants remain. These are fixed by matching. For initial-value problems where the layer is at the initial point, the inner expansion inherits all initial conditions. The outer expansion has no imposed conditions at all — its free constants are determined entirely by matching. [/remark] ## Matching Neither expansion is valid everywhere. The outer expansion breaks down near $\bar{x}$, and the inner expansion breaks down far from $\bar{x}$. But if the assumed asymptotic structure is correct, there must be an **overlap region** — a range of $x$ values where both approximations are valid and agree. The requirement that they agree in the overlap determines the remaining free constants. This procedure is called **asymptotic matching**. ### The intermediate variable method The most rigorous matching approach introduces a third variable that lives *between* the inner and outer scales. The idea is: if both the inner and outer expansions are valid approximations to $y$, then they must both give the same approximation to $y$ when evaluated at any intermediate scale. **Setup.** Let $\alpha \in (0,1)$ and define the intermediate variable \begin{align*} \eta := \frac{x - \bar{x}}{\delta(\epsilon)^\alpha}, \end{align*} and the function $Z(\eta;\epsilon) := y(\bar{x} + \delta^\alpha \eta;\epsilon)$. The point is that $\delta^\alpha$ is an intermediate scale: since $0 < \alpha < 1$, we have $\delta \ll \delta^\alpha \ll 1$, so the intermediate region is wider than the inner region but narrower than the outer region. In terms of the inner variable, $\xi = \delta^{\alpha-1}\eta \to +\infty$ as $\epsilon \to 0$ (since $\alpha - 1 < 0$), so we are far from the centre of the layer. In terms of the outer variable, $x = \bar{x} + \delta^\alpha \eta \to \bar{x}$ as $\epsilon \to 0$, so we are close to the layer. This is the overlap zone. **The procedure.** We construct two approximations to $Z(\eta;\epsilon)$: 1. Define $Z_{\text{from outer}}(\eta;\epsilon) := y_{\text{out}}^{(N)}(\bar{x} + \delta^\alpha\eta;\epsilon)$ — the outer expansion evaluated at $x = \bar{x} + \delta^\alpha\eta$. Then expand $Z_{\text{from outer}}$ as $\epsilon \to 0$ with $\eta$ fixed. 2. Define $Z_{\text{from inner}}(\eta;\epsilon) := Y_{\text{in}}^{(N)}(\delta^{\alpha-1}\eta;\epsilon)$ — the inner expansion evaluated at $\xi = \delta^{\alpha-1}\eta$. Then expand $Z_{\text{from inner}}$ as $\epsilon \to 0$ with $\eta$ fixed. If both expansions are valid in the overlap, these two functions must agree to the required order. Setting $Z_{\text{from outer}} = Z_{\text{from inner}}$ (order by order in $\epsilon$) determines the free constants. The range of $\alpha$ for which the comparison is consistent defines the overlap region. [example: Intermediate Variable Matching for the IVP] For the problem $\epsilon y'' + y' + y = 0$, $y(0) = 0$, $y'(0) = 1$, the two-term expansions are: \begin{align*} y^{(2)}_{\mathrm{out}}(t;\epsilon) &= \epsilon \, e^{-t} + \epsilon^2(-t e^{-t} + b e^{-t}), \\[4pt] Y^{(2)}_{\mathrm{in}}(\tau;\epsilon) &= \epsilon(1 - e^{-\tau}) + \epsilon^2(2 - \tau - (2+\tau)e^{-\tau}), \end{align*} where $\tau = t/\epsilon$ and $b$ is an unknown eigensolution constant. We define the intermediate variable $\eta = t/\epsilon^\alpha$ with $0 < \alpha < 1$, and the function $Z(\eta;\epsilon) = y(\epsilon^\alpha\eta;\epsilon)$. **Approximation 1:** Define $Z_{\text{from outer}}(\eta;\epsilon) := y_{\text{out}}^{(2)}(\epsilon^\alpha\eta;\epsilon)$ — evaluate the outer expansion at $t = \epsilon^\alpha\eta$: \begin{align*} Z_{\text{from outer}}(\eta;\epsilon) = \epsilon\,e^{-\epsilon^\alpha\eta} + \epsilon^2\bigl(-\epsilon^\alpha\eta\,e^{-\epsilon^\alpha\eta} + b\,e^{-\epsilon^\alpha\eta}\bigr). \end{align*} Now expand for small $\epsilon$ with $\eta$ fixed. Since $\epsilon^\alpha \to 0$, use $e^{-\epsilon^\alpha\eta} = 1 - \epsilon^\alpha\eta + \frac{1}{2}\epsilon^{2\alpha}\eta^2 - \cdots$: *First term:* \begin{align*} \epsilon\,e^{-\epsilon^\alpha\eta} = \epsilon(1 - \epsilon^\alpha\eta + O(\epsilon^{2\alpha})) = \epsilon - \epsilon^{1+\alpha}\eta + O(\epsilon^{1+2\alpha}). \end{align*} *Second term:* \begin{align*} -\epsilon^2 \cdot \epsilon^\alpha\eta\,e^{-\epsilon^\alpha\eta} = -\epsilon^{2+\alpha}\eta(1 + O(\epsilon^\alpha)) = -\epsilon^{2+\alpha}\eta + O(\epsilon^{2+2\alpha}). \end{align*} *Third term:* \begin{align*} b\epsilon^2\,e^{-\epsilon^\alpha\eta} = b\epsilon^2(1 - \epsilon^\alpha\eta + \cdots) = b\epsilon^2 + O(\epsilon^{2+\alpha}). \end{align*} Combining all three: \begin{align*} Z_{\text{from outer}}(\eta;\epsilon) = \epsilon - \epsilon^{1+\alpha}\eta + b\epsilon^2 + O(\epsilon^{\min(1+2\alpha,\,2+\alpha)}). \end{align*} **Approximation 2:** Define $Z_{\text{from inner}}(\eta;\epsilon) := Y_{\text{in}}^{(2)}(\epsilon^{\alpha-1}\eta;\epsilon)$ — evaluate the inner expansion at $\tau = \epsilon^{\alpha-1}\eta$. (Note: $\tau = t/\epsilon = \epsilon^\alpha\eta/\epsilon = \epsilon^{\alpha-1}\eta$.) Since $\alpha < 1$, we have $\alpha - 1 < 0$, so $\tau = \epsilon^{\alpha-1}\eta \to +\infty$ as $\epsilon \to 0$. Every $e^{-\tau}$ term vanishes beyond all algebraic orders (exponentially small, denoted e.s.t.): \begin{align*} Z_{\text{from inner}}(\eta;\epsilon) &= \epsilon\bigl(1 - \underbrace{e^{-\epsilon^{\alpha-1}\eta}}_{\text{e.s.t.}}\bigr) + \epsilon^2\bigl(2 - \epsilon^{\alpha-1}\eta - (2 + \epsilon^{\alpha-1}\eta)\underbrace{e^{-\epsilon^{\alpha-1}\eta}}_{\text{e.s.t.}}\bigr). \end{align*} Dropping all e.s.t. terms: \begin{align*} &= \epsilon + \epsilon^2(2 - \epsilon^{\alpha-1}\eta). \end{align*} Expanding the last term: $\epsilon^2 \cdot \epsilon^{\alpha-1}\eta = \epsilon^{2+\alpha-1}\eta = \epsilon^{1+\alpha}\eta$. So: \begin{align*} Z_{\text{from inner}}(\eta;\epsilon) = \epsilon + 2\epsilon^2 - \epsilon^{1+\alpha}\eta + \text{e.s.t.} \end{align*} **Comparison.** Both $Z_{\text{from outer}}$ and $Z_{\text{from inner}}$ are now functions of $\eta$ and $\epsilon$. Setting them equal: \begin{align*} \underbrace{\epsilon - \epsilon^{1+\alpha}\eta + b\epsilon^2}_{Z_{\text{from outer}}} + \cdots = \underbrace{\epsilon - \epsilon^{1+\alpha}\eta + 2\epsilon^2}_{Z_{\text{from inner}}} + \cdots \end{align*} The $\epsilon$ terms match. The $\epsilon^{1+\alpha}\eta$ terms match. What remains is: \begin{align*} b\epsilon^2 + \cdots = 2\epsilon^2 + \cdots \end{align*} But this comparison is only valid if $\epsilon^2$ is actually the next term in the hierarchy — i.e. if $\epsilon^{1+\alpha}$ is *larger* than $\epsilon^2$ (so it has already been accounted for, and $\epsilon^2$ is genuinely the next order). This requires $1 + \alpha < 2$, i.e. $\alpha < 1$ (which we already assumed). We also need $\epsilon^{1+\alpha}$ to be larger than $\epsilon^2$, which requires $1 + \alpha < 2$, always true. However, the error terms in $Z_{\text{from outer}}$ are $O(\epsilon^{\min(1+2\alpha,\,2+\alpha)})$. For these to be *smaller* than $\epsilon^2$ (so they don't contaminate our comparison), we need $\min(1+2\alpha, 2+\alpha) > 2$. The binding constraint is $1 + 2\alpha > 2$, i.e. $\alpha > 1/2$. So the comparison is valid for $1/2 < \alpha < 1$, and matching the $\epsilon^2$ coefficients gives $b = 2$. The restriction $\alpha > 1/2$ corresponds to the overlap region $\epsilon \ll t \ll \epsilon^{1/2}$ (since $t = \epsilon^\alpha\eta$ with $\eta = O(1)$, and $\alpha \in (1/2, 1)$ means $t$ ranges between $O(\epsilon)$ and $O(\epsilon^{1/2})$). Note that this overlap region is narrower than the leading-order overlap — it is typical for the overlap domain to shrink as the matching order increases. [/example] ### Van Dyke's matching principle The intermediate variable method, while rigorous, requires choosing $\alpha$ and tracking the overlap region. **Van Dyke's matching principle** is a simpler procedure that produces the same result in most cases. The principle states: > *The inner expansion to order $\mu_i$, re-expanded in the outer limit to order $\mu_o$, equals the outer expansion to order $\mu_o$, re-expanded in the inner limit to order $\mu_i$.* **The procedure.** Van Dyke matching at orders $\mu_o : \mu_i$ consists of three steps. **Step 1 (Re-expand outer in inner coordinates).** Take the outer expansion $y_{\text{out}}^{(\mu_o)}(x;\epsilon)$ and compose it with the inner coordinate by setting $x = \bar{x} + \delta\xi$: \begin{align*} W_1(\xi;\epsilon) := y_{\text{out}}^{(\mu_o)}(\bar{x} + \delta\xi;\epsilon). \end{align*} Expand $W_1$ for small $\epsilon$ with $\xi$ **fixed**, keeping terms to order $\mu_i$. This strips the outer expansion down to what it looks like from inside the layer — high-order $x$-structure simplifies because $x = \bar{x} + \delta\xi$ is close to $\bar{x}$. **Step 2 (Re-expand inner in outer coordinates).** Take the inner expansion $Y_{\text{in}}^{(\mu_i)}(\xi;\epsilon)$ and compose it with the outer coordinate by setting $\xi = (x - \bar{x})/\delta$: \begin{align*} W_2(x;\epsilon) := Y_{\text{in}}^{(\mu_i)}\!\left(\frac{x - \bar{x}}{\delta};\epsilon\right). \end{align*} Expand $W_2$ for small $\epsilon$ with $x$ **fixed**, keeping terms to order $\mu_o$. This strips the inner expansion down to what it looks like from outside the layer — since $\xi = (x - \bar{x})/\delta \to \infty$, exponentials like $e^{-\xi}$ become e.s.t. and vanish, leaving only the algebraic skeleton. **Step 3 (Compare).** After steps 1 and 2, both limits have been taken: $W_1$ is an algebraic expression in $(\xi, \epsilon)$ and $W_2$ is an algebraic expression in $(x, \epsilon)$. At this stage, convert one into the other's coordinate using $x = \bar{x} + \delta\xi$ and compare term by term. Setting the expressions equal determines the free constants. [remark: Why This Works] Steps 1 and 2 do the real work: each limit kills the parts of the respective expansion that do not belong in the overlap region. Step 1 kills the outer expansion's large-scale $x$-structure (reducing it to its near-$\bar{x}$ form), and step 2 kills the inner expansion's boundary-layer structure (exponentials that decay away from $\bar{x}$). What remains on both sides is the overlap form of $y$ — the part that both expansions share. Step 3 is then just algebra: the two overlap forms, expressed in different coordinates, must agree. [/remark] The orders $\mu_o$ and $\mu_i$ can be chosen freely, subject to the important restriction discussed in the following remark. [remark: Do Not Cut Between Logarithms — and Why] **The order-jumping phenomenon.** The root cause of logarithmic difficulty is that the same term can change its asymptotic order under rescaling. Consider $\epsilon \ln x$. With $x$ fixed, this is simply $O(\epsilon)$. But substituting $x = \epsilon\xi$ (i.e. evaluating in the inner variable): \begin{align*} \epsilon \ln(\epsilon\xi) = \epsilon\ln\epsilon + \epsilon\ln\xi. \end{align*} The first term is $O(\epsilon\ln\epsilon)$, which is *larger* than $O(\epsilon)$ as $\epsilon \to 0$ — the logarithm has promoted the term to a higher asymptotic order. No power of $\epsilon$ has this property: $\epsilon \cdot \epsilon^a = \epsilon^{1+a}$, which is always *smaller* than $\epsilon$. Logarithms are the unique source of order-jumping under rescaling. **The asymptotic sequence.** When logarithms appear, the clean power-[series](/page/Series) sequence $1, \epsilon, \epsilon^2, \ldots$ is replaced by a finer sequence that includes logarithmic sub-orders within each power of $\epsilon$: \begin{align*} \underbrace{\ln\epsilon, \; 1}_{\epsilon^0 \text{ group}}; \qquad \underbrace{\epsilon\ln\epsilon, \; \epsilon}_{\epsilon^1 \text{ group}}; \qquad \underbrace{\epsilon^2\ln^2\epsilon, \; \epsilon^2\ln\epsilon, \; \epsilon^2}_{\epsilon^2 \text{ group}}; \qquad \ldots \end{align*} Within each $\epsilon^n$ group, successive terms differ by a factor of $1/\ln\epsilon$ — a gap that shrinks very slowly. Between [groups](/page/Group), the gap is a full factor of $\epsilon$ — a genuine asymptotic separation. The rule "do not cut between logarithms" means: **you can safely truncate between $\epsilon^n$ groups** (at $1$, $\epsilon$, $\epsilon^2$, ...) **but not within a group** (not at $\epsilon\ln\epsilon$, for instance). **Why the gap matters: the ordering-matching circularity.** The warning "do not cut between logarithms" is often stated as a rule to memorise. The usual justification — that the gap between successive logarithmic orders like $1/\ln\epsilon$ and $1/\ln^2\epsilon$ is too small for clean truncation — is correct but incomplete. The deeper issue is a **circularity** in the matching procedure itself. Consider a typical situation: the outer expansion contains a term $a_0 \ln x$, where $a_0$ is a constant to be determined by matching. The matching procedure requires evaluating this outer expansion at $x = \bar{x} + \delta^\alpha\eta$ (or equivalently at $x = \delta\xi$ if $\bar{x} = 0$), which produces a term $a_0 \ln(\delta^\alpha\eta) = a_0 \alpha \ln\delta + a_0\ln\eta$. The first piece, $a_0 \alpha\ln\delta$, depends on $\epsilon$ (since $\delta = \delta(\epsilon)$), and its asymptotic order depends on the **size of $a_0$**: - If $a_0 = O(1)$, then $a_0\ln\delta = O(\ln\epsilon)$, which **blows up** as $\epsilon \to 0$ — it is larger than any $O(1)$ term. - If $a_0 = O(1/\ln\epsilon)$, then $a_0\ln\delta = O(1)$ — it is a bounded correction. This creates a circular dependency: **you cannot determine the asymptotic ordering of terms until you know $a_0$, but you cannot find $a_0$ until you have fixed the ordering and performed the matching.** In a standard power-series problem, this circularity does not arise: all constants are $O(1)$, and multiplying an $O(1)$ constant by a power of $\epsilon$ does not change its order. But when $\ln\epsilon$ appears, constants can absorb or cancel logarithmic factors, shifting terms between asymptotic orders. The switchback example below (Steps 3–6) demonstrates this concretely. Naive matching assumes $a_0 = O(1)$, which forces $a_0 \ln\epsilon$ into the ordering as a term that is larger than $O(1)$ — and the resulting matching equations are inconsistent ($A_0 = 0$ and $1 = 0$). The resolution is that $a_0 = -1/\ln\epsilon$, which is *not* $O(1)$ — it sits at an intermediate logarithmic order that the naive truncation did not anticipate. The fix is to **group all terms by integer powers of $\epsilon$** — that is, collect everything multiplying $\epsilon^0$, everything multiplying $\epsilon^1$, everything multiplying $\epsilon^2$, and so on — while **allowing the coefficients to depend on $\ln\epsilon$**. Instead of treating $1$ and $1/\ln\epsilon$ as separate asymptotic orders, write: \begin{align*} y \sim y_0(x;\epsilon) + \epsilon\,y_1(x;\epsilon) + \epsilon^2\,y_2(x;\epsilon) + \cdots, \end{align*} where each $y_k(x;\epsilon)$ may depend on $\ln\epsilon$ but not on powers of $\epsilon$. For instance, $y_0$ might be $x + 1 - (\ln x)/\ln\epsilon$. You then match at each power of $\epsilon$ ($\epsilon^0$, $\epsilon^1$, $\epsilon^2$, ...) where genuine asymptotic gaps exist, and never at intermediate logarithmic orders like $\epsilon\ln\epsilon$. The underlying principle is not specific to logarithms. *Any* function of $\epsilon$ that decays more slowly than every positive power of $\epsilon$ — such as $1/\ln\epsilon$, $1/\ln\ln\epsilon$, or $1/\ln^k\epsilon$ — would cause the same issue if it appeared in an expansion. Logarithms are simply the most common example in practice, arising from problems with degenerate operators or cylindrical/spherical geometries. See the page on [matched asymptotic expansions with logarithmic terms](/pages/1138) for detailed worked examples. **Intermediate variable matching and the truncation issue.** The restriction on truncation order applies specifically to Van Dyke matching: it is about which orders $\mu_o : \mu_i$ you choose before performing the three-step procedure. Intermediate variable matching does not require a truncation decision — it simply expands two compositions as $\epsilon \to 0$ and compares them — so it avoids this particular pitfall. However, **both methods** suffer from the more fundamental **separation problem** when logarithms are present. In either method, the comparison step requires you to match terms on the left and right sides. If you naively separate "coefficient of $\ln\epsilon$" from "remaining constants," you get contradictions — as the switchback example (Step 3) demonstrates using intermediate variables, not Van Dyke. The root cause is the same in both cases: $a_0\ln\epsilon$ looks like it should be matched separately, but $a_0 = O(1/\ln\epsilon)$ makes $a_0\ln\epsilon = O(1)$, so it belongs with the constants. The resolution is also the same for both methods: use **generalised expansions** — group terms by algebraic power of $\epsilon$, allow coefficients to depend on $\ln\epsilon$, and match only at algebraic orders ($1$, $\epsilon$, $\epsilon^2$, ...), never at intermediate logarithmic sub-orders. Once you adopt this approach, either matching method works without modification. The worked examples below demonstrate both. **Alternative: separating logarithmic orders explicitly.** Instead of collecting terms by algebraic order, one can explicitly include logarithmic gauge functions in the expansions. For instance, if the inner solution has a prefactor $A_0 = O(1/\ln\epsilon)$, one could pose \begin{align*} \text{outer: } y \sim y_0(x) + \frac{1}{\ln\epsilon}\,y_{0,1}(x) + \cdots, \qquad \text{inner: } Y \sim \frac{1}{\ln\epsilon}\,Y_{0,1}(\xi) + \cdots \end{align*} and match at each logarithmic sub-order separately. The functions $y_{0,1}$ and $Y_{0,1}$ satisfy the *same* differential equations as $y_0$ and $Y_0$ (since the ODE is linear and the $1/\ln\epsilon$ factor divides out), typically with different inhomogeneous terms or boundary conditions. This approach is completely rigorous and can be matched using intermediate variables without any special precautions. However, it is less economical than the generalised-expansion approach: one must anticipate which logarithmic orders appear and solve a separate problem for each. For most purposes, the generalised approach (grouping by algebraic order) is preferred. [/remark] [example: Van Dyke Matching for the IVP] For the problem $\epsilon y'' + y' + y = 0$, $y(0) = 0$, $y'(0) = 1$, we perform $\epsilon^2 : \epsilon^2$ Van Dyke matching. **Direction 1:** Define $W_1(\tau;\epsilon) := y_{\text{out}}^{(\epsilon^2)}(\epsilon\tau;\epsilon)$ — the outer expansion evaluated at $t = \epsilon\tau$: \begin{align*} W_1(\tau;\epsilon) = \epsilon \, e^{-\epsilon\tau} + \epsilon^2(-\epsilon\tau \, e^{-\epsilon\tau} + b\,e^{-\epsilon\tau}). \end{align*} Expand $W_1$ for small $\epsilon$ with $\tau$ fixed. Using $e^{-\epsilon\tau} = 1 - \epsilon\tau + \frac{1}{2}\epsilon^2\tau^2 - \cdots$: \begin{align*} W_1(\tau;\epsilon) &= \epsilon(1 - \epsilon\tau + O(\epsilon^2)) + \epsilon^2\bigl(-\epsilon\tau(1 - \epsilon\tau + \cdots) + b(1 - \epsilon\tau + \cdots)\bigr) \\ &= \epsilon - \epsilon^2\tau + O(\epsilon^3) + \epsilon^2 b - \epsilon^3 b\tau - \epsilon^3\tau + O(\epsilon^4) \\ &= \epsilon - \epsilon^2\tau + b\epsilon^2 + O(\epsilon^3). \end{align*} **Direction 2:** Define $W_2(t;\epsilon) := Y_{\text{in}}^{(\epsilon^2)}(t/\epsilon;\epsilon)$ — the inner expansion evaluated at $\tau = t/\epsilon$: \begin{align*} W_2(t;\epsilon) = \epsilon(1 - e^{-t/\epsilon}) + \epsilon^2\bigl(2 - t/\epsilon - (2 + t/\epsilon)e^{-t/\epsilon}\bigr). \end{align*} Expand $W_2$ for small $\epsilon$ with $t > 0$ fixed. Since $t/\epsilon \to +\infty$, every $e^{-t/\epsilon}$ term is exponentially small (e.s.t.): \begin{align*} W_2(t;\epsilon) &= \epsilon(1 - \text{e.s.t.}) + \epsilon^2(2 - t/\epsilon - \text{e.s.t.}) \\ &= \epsilon + 2\epsilon^2 - \epsilon^2 \cdot \frac{t}{\epsilon} + \text{e.s.t.} \\ &= \epsilon + 2\epsilon^2 - \epsilon t + \text{e.s.t.} \end{align*} **Set equal.** Van Dyke's principle requires $W_1 = W_2$ (to the required order). To compare, note that $W_1$ is a function of $\tau$ while $W_2$ is a function of $t$, related by $t = \epsilon\tau$. Substituting $\epsilon t = \epsilon^2\tau$ in $W_2$: \begin{align*} \epsilon - \epsilon^2\tau + b\epsilon^2 = \epsilon - \epsilon^2\tau + 2\epsilon^2, \end{align*} giving $b = 2$. [/example] ## Composite Approximations Having determined both the inner and outer expansions, it is often useful to form a single approximation valid across the entire domain. The idea is straightforward: add the inner and outer expansions, and subtract their common part (the "overlap" value) so as not to double-count in the region where both are valid. At leading order, if $y_{\text{out}}$ and $Y_{\text{in}}$ share a common limiting value $y_{\text{overlap}}$ in the overlap region, the **composite approximation** is: \begin{align*} y_{\text{comp}}(x;\epsilon) = y_{\text{out}}(x;\epsilon) + Y_{\text{in}}\!\left(\frac{x - \bar{x}}{\delta};\epsilon\right) - y_{\text{overlap}}(\epsilon). \end{align*} [example: Composite Approximation for the IVP] For $\epsilon y'' + y' + y = 0$, $y(0)=0$, $y'(0)=1$, the leading-order expansions are: \begin{align*} y_{\text{out}} = \epsilon \, e^{-t}, \qquad Y_{\text{in}} = \epsilon(1 - e^{-t/\epsilon}). \end{align*} In the overlap region ($\epsilon \ll t \ll 1$): $y_{\text{out}} \approx \epsilon$ and $Y_{\text{in}} \approx \epsilon$. So $y_{\text{overlap}} = \epsilon$. The composite is: \begin{align*} y_{\text{comp}} = \epsilon \, e^{-t} + \epsilon(1 - e^{-t/\epsilon}) - \epsilon = \epsilon(e^{-t} - e^{-t/\epsilon}). \end{align*} This reduces to the inner approximation when $t = O(\epsilon)$ (since $e^{-t} \approx 1$) and to the outer approximation when $t = O(1)$ (since $e^{-t/\epsilon}$ is exponentially small). [/example] ## Worked Example: A Boundary-Value Problem [example: Singularly Perturbed Boundary-Value Problem] Consider the boundary-value problem \begin{align*} \epsilon \frac{d^2 y}{dx^2} + \frac{dy}{dx} = x, \qquad y(0) = 0, \quad y(1) = 1, \qquad \epsilon \searrow 0. \end{align*} **Step 1: Regular expansion fails.** We pose the outer expansion $y(x;\epsilon) \sim y_0(x) + \epsilon\,y_1(x)$ with $x$ fixed. Substituting into $\epsilon y'' + y' = x$: \begin{align*} \epsilon(y_0'' + \epsilon\,y_1'' + \cdots) + (y_0' + \epsilon\,y_1' + \cdots) = x. \end{align*} Grouping by powers of $\epsilon$: \begin{align*} \underbrace{y_0'}_{O(1)} + \epsilon\underbrace{(y_0'' + y_1')}_{O(\epsilon)} + \cdots = x. \end{align*} At $O(1)$:\begin{align*} \frac{dy_0}{dx} = x. \end{align*} Integrating directly: $y_0 = \frac{1}{2}x^2 + C$, where $C$ is an integration constant. Attempting to impose both boundary conditions: $y_0(0) = 0$ requires $C = 0$, while $y_0(1) = 1$ requires $\frac{1}{2} + C = 1$, i.e. $C = \frac{1}{2}$. These are incompatible — the expansion cannot hold on all of $[0,1]$. **Step 2: Rescaling and dominant balance.** We assume a layer of width $\delta(\epsilon) \to 0$ near some point $\bar{x}$. Define $\xi = (x - \bar{x})/\delta$ and $Y(\xi;\epsilon) = y(\bar{x} + \delta\xi;\epsilon)$. By the chain rule: \begin{align*} \frac{dy}{dx} = \frac{1}{\delta}\frac{dY}{d\xi}, \qquad \frac{d^2 y}{dx^2} = \frac{1}{\delta^2}\frac{d^2 Y}{d\xi^2}. \end{align*} Substituting into $\epsilon y'' + y' = x$: \begin{align*} \frac{\epsilon}{\delta^2}\frac{d^2 Y}{d\xi^2} + \frac{1}{\delta}\frac{dY}{d\xi} = \bar{x} + \delta\xi. \end{align*} We assume $Y$ has a Poincaré expansion with $Y_0, Y_0', Y_0''$ all $O(1)$. The coefficients are $\epsilon/\delta^2$, $1/\delta$, and $O(1)$ (RHS). We check all pairings: **Balance $\epsilon/\delta^2 \sim 1/\delta$:** gives $\delta \sim \epsilon$. Then both left-hand coefficients are $O(1/\epsilon)$, while the right-hand side is $O(1) \ll O(1/\epsilon)$. ✓ Consistent — the RHS enters as a lower-order correction. **Balance $\epsilon/\delta^2 \sim 1$ (RHS):** gives $\delta \sim \epsilon^{1/2}$. Then $1/\delta \sim \epsilon^{-1/2}$, which is larger than $\epsilon/\delta^2 = O(1)$. The $Y'$ term dominates with nothing to balance it. ✗ Inconsistent. **Balance $1/\delta \sim 1$ (RHS):** gives $\delta \sim 1$. Then $\epsilon/\delta^2 \sim \epsilon \ll 1$, so $Y''$ is negligible — this is the regular expansion, not a boundary layer. ✗ Contradicts the assumption that $Y''$ participates. Only the first balance is consistent: $\delta(\epsilon) = \epsilon$. Setting $\delta = \epsilon$ in the transformed equation and multiplying through by $\epsilon$: \begin{align*} \frac{d^2 Y}{d\xi^2} + \frac{dY}{d\xi} = \epsilon(\bar{x} + \epsilon\xi). \end{align*} **Step 3: Layer location.** We pose the inner Poincaré expansion $Y(\xi;\epsilon) \sim Y_0(\xi) + \epsilon\,Y_1(\xi) + \cdots$. At $O(1)$, the right-hand side drops out, giving the leading-order inner equation: \begin{align*} Y_0'' + Y_0' = 0. \end{align*} This is a second-order linear ODE with constant coefficients. The characteristic equation is $r^2 + r = r(r+1) = 0$, with roots $r = 0$ and $r = -1$. The general solution is therefore: \begin{align*} Y_0(\xi) = c_1 + c_2\,e^{-\xi}, \end{align*} or equivalently, $Y_0'(\xi) = A_0\,e^{-\xi}$ where $A_0 = -c_2$. Now we determine the layer location using the matching consistency argument from the previous section. The inner solution is $Y_0(\xi) = c_1 + c_2\,e^{-\xi}$. Matching requires composing $Y_0$ with $\xi = \delta^{\alpha-1}\eta$ and taking $\epsilon \to 0$ with $\eta$ fixed. Since $\delta = \epsilon$ and $\alpha < 1$, we have $\delta^{\alpha-1} = \epsilon^{\alpha-1} \to +\infty$. If $\bar{x} = 0$, the outer region has $x > 0$, so $\eta > 0$. Then $Y_0(\epsilon^{\alpha-1}\eta) = c_1 + c_2\,e^{-\epsilon^{\alpha-1}\eta} \to c_1$ as $\epsilon \to 0$, which is finite. Matching can proceed. If $\bar{x} = 1$, the outer region has $x < 1$, so $\eta < 0$. Then $Y_0(\epsilon^{\alpha-1}\eta) = c_1 + c_2\,e^{\epsilon^{\alpha-1}|\eta|} \to +\infty$ (for $c_2 \neq 0$). The composition diverges, so matching with the bounded outer solution is impossible. Setting $c_2 = 0$ avoids the divergence but leaves $Y_0$ constant, unable to form a boundary layer. If $\bar{x}$ were an interior point, the outer region would lie on both sides, requiring the $\epsilon \to 0$ limit of the composition to be finite for both signs of $\eta$. The exponential $e^{-\epsilon^{\alpha-1}\eta}$ diverges for $\eta < 0$ and vanishes for $\eta > 0$, so it cannot be bounded in both directions — only a trivial inner solution would survive. Therefore $\bar{x} = 0$: the boundary layer is at the left endpoint, and $\xi = x/\epsilon$. **Step 4: Distribute boundary conditions.** The inner expansion lives near $x = 0$, so it inherits $y(0) = 0$. Expanding: $Y(0;\epsilon) = 0$ gives $Y_0(0) = 0$ at $O(1)$ and $Y_1(0) = 0$ at $O(\epsilon)$. The outer expansion is valid on $(0,1]$, so it inherits $y(1) = 1$. Expanding: $y_0(1) = 1$ at $O(1)$ and $y_1(1) = 0$ at $O(\epsilon)$. **Step 5: Solve the outer expansion.** *Leading order.* We have $y_0' = x$ with $y_0(1) = 1$. Integrating: \begin{align*} y_0(x) = \frac{1}{2}x^2 + C. \end{align*} The boundary condition $y_0(1) = 1$ gives $\frac{1}{2} + C = 1$, so $C = \frac{1}{2}$: \begin{align*} y_0(x) = \frac{1}{2}x^2 + \frac{1}{2}. \end{align*} *Order $\epsilon$.* From the grouped expansion above, the $O(\epsilon)$ equation is $y_0'' + y_1' = 0$, i.e.: \begin{align*} y_1' = -y_0''. \end{align*} Computing: $y_0 = \frac{1}{2}x^2 + \frac{1}{2}$, so $y_0' = x$, and therefore $y_0'' = 1$. Thus $y_1' = -1$. Integrating: \begin{align*} y_1(x) = -x + D, \qquad y_1(1) = -1 + D = 0 \implies D = 1. \end{align*} Therefore $y_1(x) = 1 - x$. The two-term outer expansion is: \begin{align*} y(x;\epsilon) \sim \frac{1}{2}x^2 + \frac{1}{2} + \epsilon(1 - x) \quad \text{as } \epsilon \to 0, \; x > 0 \text{ fixed}. \end{align*} **Step 6: Solve the inner expansion.** *Leading order.* We have $Y_0(\xi) = c_1 + c_2\,e^{-\xi}$ from Step 3. Applying $Y_0(0) = 0$: \begin{align*} Y_0(0) = c_1 + c_2 = 0 \implies c_1 = -c_2. \end{align*} Writing $A_0 = c_1 = -c_2$: \begin{align*} Y_0(\xi) = A_0 - A_0\,e^{-\xi} = A_0(1 - e^{-\xi}). \end{align*} The constant $A_0$ is **not** determined — we have one boundary condition ($Y_0(0) = 0$) and a second-order ODE, so one constant remains free. This will be fixed by matching. *Order $\epsilon$.* The inner equation is $Y'' + Y' = \epsilon(\bar{x} + \epsilon\xi)$ with $\bar{x} = 0$, i.e. $Y'' + Y' = \epsilon^2\xi$. Substituting $Y \sim Y_0 + \epsilon Y_1 + \cdots$: \begin{align*} (Y_0'' + \epsilon\,Y_1'' + \cdots) + (Y_0' + \epsilon\,Y_1' + \cdots) = \epsilon^2\xi. \end{align*} At $O(1)$ we recover $Y_0'' + Y_0' = 0$ (already solved). At $O(\epsilon)$, the right-hand side does not contribute (it is $O(\epsilon^2)$), so: \begin{align*} Y_1'' + Y_1' = 0, \qquad Y_1(0) = 0. \end{align*} is the same ODE as before, with the same boundary condition. The solution is: \begin{align*} Y_1(\xi) = A_1(1 - e^{-\xi}), \end{align*} with $A_1$ a new free constant. The two-term inner expansion is: \begin{align*} Y(\xi;\epsilon) \sim A_0(1 - e^{-\xi}) + \epsilon\,A_1(1 - e^{-\xi}) \quad \text{as } \epsilon \to 0, \; \xi \text{ fixed}. \end{align*} **Step 7: Van Dyke matching at $\epsilon : \epsilon$.** We perform Van Dyke matching with $\mu_o = \epsilon$ and $\mu_i = \epsilon$. *Direction 1 (Outer $\to$ inner).* Define $W_1(\xi;\epsilon) := y_{\text{out}}^{(\epsilon)}(\epsilon\xi;\epsilon)$ — the outer expansion evaluated at $x = \epsilon\xi$: \begin{align*} W_1(\xi;\epsilon) = \frac{1}{2}\epsilon^2\xi^2 + \frac{1}{2} + \epsilon(1 - \epsilon\xi) = \frac{1}{2}\epsilon^2\xi^2 + \frac{1}{2} + \epsilon - \epsilon^2\xi. \end{align*} Expand $W_1$ for small $\epsilon$ with $\xi$ fixed, keeping terms to order $\epsilon$ (so the $\epsilon^2$ terms drop): \begin{align*} W_1(\xi;\epsilon) = \frac{1}{2} + \epsilon + O(\epsilon^2). \end{align*} *Direction 2 (Inner $\to$ outer).* Define $W_2(x;\epsilon) := Y_{\text{in}}^{(\epsilon)}(x/\epsilon;\epsilon)$ — the inner expansion evaluated at $\xi = x/\epsilon$: \begin{align*} W_2(x;\epsilon) = A_0(1 - e^{-x/\epsilon}) + \epsilon\,A_1(1 - e^{-x/\epsilon}). \end{align*} Expand $W_2$ for small $\epsilon$ with $x > 0$ fixed. Since $x > 0$ is fixed and $\epsilon \to 0$, we have $x/\epsilon \to +\infty$, so $e^{-x/\epsilon}$ is exponentially small (beyond all algebraic orders). Therefore: \begin{align*} W_2(x;\epsilon) = A_0 + \epsilon\,A_1 + \text{e.s.t.} \end{align*} *Set equal.* Van Dyke's principle requires $W_1 = W_2$ (to order $\epsilon$). Both are now expressed as functions of $\epsilon$ alone (the $\xi$-dependence in $W_1$ and the $x$-dependence in $W_2$ have been absorbed into higher-order terms): \begin{align*} \frac{1}{2} + \epsilon = A_0 + \epsilon\,A_1. \end{align*} Matching at $O(1)$: $A_0 = \frac{1}{2}$. Matching at $O(\epsilon)$: $A_1 = 1$. All constants are now determined. The final inner expansion is: \begin{align*} Y(\xi;\epsilon) \sim \frac{1}{2}(1 - e^{-\xi}) + \epsilon(1 - e^{-\xi}). \end{align*} [/example] ## Worked Example: A Nonlinear Problem with an Interior Layer [example: Nonlinear Problem with a Singularity in the Outer Solution] Consider the initial-value problem \begin{align*} (x + \epsilon y)\frac{dy}{dx} + y = 1, \qquad y(1) = 2, \qquad \epsilon \searrow 0. \end{align*} The domain of $x$ is the entire real line, and the nonlinear term $\epsilon y\,y'$ makes this problem qualitatively different from the linear examples above. **Step 1: Outer expansion.** We pose $y(x;\epsilon) \sim y_0(x) + \epsilon\,y_1(x)$ with $x$ fixed. Substituting into $(x + \epsilon y)y' + y = 1$: \begin{align*} \bigl(x + \epsilon(y_0 + \epsilon y_1 + \cdots)\bigr)(y_0' + \epsilon y_1' + \cdots) + (y_0 + \epsilon y_1 + \cdots) = 1. \end{align*} Expanding the product $(x + \epsilon y_0 + \cdots)(y_0' + \epsilon y_1' + \cdots) = x\,y_0' + \epsilon(x\,y_1' + y_0\,y_0') + O(\epsilon^2)$. Adding the remaining $y$ terms: \begin{align*} \underbrace{x\,y_0' + y_0}_{O(1)} + \epsilon\underbrace{(x\,y_1' + y_0\,y_0' + y_1)}_{O(\epsilon)} + O(\epsilon^2) = 1. \end{align*} *At $O(1)$:* \begin{align*} x\,y_0' + y_0 = 1. \end{align*} The left-hand side is the derivative of $x\,y_0$: $(x\,y_0)' = 1$. Integrating: $x\,y_0 = x + c$, so $y_0 = 1 + c/x$. Applying $y_0(1) = 2$: $2 = 1 + c$, giving $c = 1$: \begin{align*} y_0(x) = 1 + \frac{1}{x}. \end{align*} This diverges as $x \to 0$. Unlike the previous examples, the layer is not at a boundary — it is at an **interior point** where the outer solution itself becomes singular. *At $O(\epsilon)$:* from the grouped expansion above, $x\,y_1' + y_0\,y_0' + y_1 = 0$, i.e.: \begin{align*} x\,y_1' + y_1 = -y_0\,y_0'. \end{align*} Computing $y_0' = -1/x^2$ (from $y_0 = 1 + 1/x$), so: \begin{align*} y_0\,y_0' = \left(1 + \frac{1}{x}\right)\!\left(-\frac{1}{x^2}\right) = -\frac{1}{x^2} - \frac{1}{x^3}. \end{align*} The equation becomes $(x\,y_1)' = -1/x^2 - 1/x^3$ (recognising $x\,y_1' + y_1 = (x\,y_1)'$). Integrating each term: $\int -x^{-2}\,dx = x^{-1}$ and $\int -x^{-3}\,dx = \frac{1}{2}x^{-2}$. So: \begin{align*} x\,y_1 = \frac{1}{x} + \frac{1}{2x^2} + d. \end{align*} Applying $y_1(1) = 0$: $1 \cdot y_1(1) = 1 + \frac{1}{2} + d = 0$, giving $d = -\frac{3}{2}$. Therefore $x\,y_1 = 1/x + 1/(2x^2) - 3/2$, and dividing by $x$: \begin{align*} y_1(x) = \frac{1}{x^2} + \frac{1}{2x^3} - \frac{3}{2x} = \frac{3}{2x} - \frac{1}{x^2} - \frac{1}{2x^3}. \end{align*} Notice that $y_1 \sim -1/(2x^3)$ as $x \to 0$, which is **more singular** than $y_0 \sim 1/x$. The perturbation expansion is breaking down near $x = 0$. **Step 2: Determine the inner scaling.** The expansion becomes disordered when $\epsilon\,y_1 \sim y_0$. Since $y_1 \sim -1/(2x^3)$ as $x \to 0$ and $y_0 \sim 1/x$, this occurs when $\epsilon/x^3 \sim 1/x$, i.e. $x^2 \sim \epsilon$, giving $x \sim \epsilon^{1/2}$. This is also a **distinguished limit** of the ODE. To verify, substitute $x = O(\epsilon^{1/2})$ and $y = O(\epsilon^{-1/2})$ (from $y_0 \sim 1/x = 1/\epsilon^{1/2}$) into $(x + \epsilon y)y' + y = 1$. Estimating $y' \sim y/x \sim \epsilon^{-1/2}/\epsilon^{1/2} = \epsilon^{-1}$, the three terms on the left-hand side are: \begin{align*} \underbrace{x\,y'}_{\epsilon^{1/2} \cdot \epsilon^{-1} = \epsilon^{-1/2}} + \underbrace{\epsilon y\,y'}_{\epsilon \cdot \epsilon^{-1/2} \cdot \epsilon^{-1} = \epsilon^{-1/2}} + \underbrace{y}_{\epsilon^{-1/2}} = 1. \end{align*} All three are $O(\epsilon^{-1/2})$, confirming a three-way balance. We therefore define the rescaled variables: \begin{align*} \xi = \frac{x}{\epsilon^{1/2}}, \qquad Y(\xi;\epsilon) = \epsilon^{1/2}\,y(\epsilon^{1/2}\xi;\epsilon). \end{align*} That is, $y = \epsilon^{-1/2}\,Y$. Computing the derivatives by the chain rule: since $x = \epsilon^{1/2}\xi$, we have $dx = \epsilon^{1/2}\,d\xi$, so \begin{align*} \frac{dy}{dx} = \frac{d(\epsilon^{-1/2}Y)}{dx} = \epsilon^{-1/2}\frac{dY}{d\xi}\frac{d\xi}{dx} = \epsilon^{-1/2} \cdot \frac{1}{\epsilon^{1/2}}\frac{dY}{d\xi} = \epsilon^{-1}\,Y'. \end{align*} Substituting $x = \epsilon^{1/2}\xi$, $y = \epsilon^{-1/2}Y$, and $y' = \epsilon^{-1}Y'$ into $(x + \epsilon y)y' + y = 1$: \begin{align*} (\epsilon^{1/2}\xi + \epsilon \cdot \epsilon^{-1/2}Y) \cdot \epsilon^{-1}Y' + \epsilon^{-1/2}Y &= 1 \\ (\epsilon^{1/2}\xi + \epsilon^{1/2}Y) \cdot \epsilon^{-1}Y' + \epsilon^{-1/2}Y &= 1 \\ \epsilon^{-1/2}(\xi + Y)\,Y' + \epsilon^{-1/2}Y &= 1. \end{align*} Multiplying through by $\epsilon^{1/2}$: \begin{align*} (\xi + Y)\frac{dY}{d\xi} + Y = \epsilon^{1/2}. \end{align*} The gauge functions go in powers of $\epsilon^{1/2}$: $Y(\xi;\epsilon) \sim Y_0(\xi) + \epsilon^{1/2}\,Y_1(\xi) + \cdots$. **Step 3: Leading-order inner solution.** At $O(1)$, the right-hand side drops out: \begin{align*} (\xi + Y_0)\,Y_0' + Y_0 = 0. \end{align*} We observe that the left-hand side is an exact derivative: \begin{align*} \xi\,Y_0' + Y_0 + Y_0\,Y_0' = \frac{d}{d\xi}(\xi\,Y_0) + \frac{d}{d\xi}\!\left(\frac{1}{2}Y_0^2\right) = \frac{d}{d\xi}\!\left(\xi\,Y_0 + \frac{1}{2}Y_0^2\right) = 0. \end{align*} Integrating: \begin{align*} \xi\,Y_0 + \frac{1}{2}Y_0^2 = a, \end{align*} where $a$ is an integration constant. This is a quadratic equation in $Y_0$. By the quadratic formula ($\frac{1}{2}Y_0^2 + \xi\,Y_0 - a = 0$): \begin{align*} Y_0 = \frac{-\xi \pm \sqrt{\xi^2 + 2a}}{1} = -\xi \pm \sqrt{\xi^2 + 2a}. \end{align*} There are two branches; both the constant $a$ and the correct branch must be determined by matching. **Step 4: Branch selection by matching.** To match with the outer solution as $\xi \to +\infty$ (i.e. $x \to +\infty$ in the inner scale), we need the large-$\xi$ behaviour of each branch. For $\xi > 0$, factor out $\xi$: \begin{align*} Y_0 = -\xi \pm \xi\sqrt{1 + \frac{2a}{\xi^2}}. \end{align*} Expanding the square root for large $\xi$ using $\sqrt{1 + u} \sim 1 + \frac{1}{2}u + \cdots$ for $u \to 0$: \begin{align*} Y_0 \sim -\xi \pm \xi\!\left(1 + \frac{a}{\xi^2} + \cdots\right) = -\xi \pm \xi \pm \frac{a}{\xi} + \cdots \end{align*} For the $+$ branch: $Y_0 \sim -\xi + \xi + a/\xi = a/\xi$. For the $-$ branch: $Y_0 \sim -\xi - \xi - a/\xi = -2\xi - a/\xi$. The outer solution has $y_0 = 1 + 1/x$, so for $x$ small and positive, $y \sim 1/x = \epsilon^{-1/2}/\xi$. Since $y = \epsilon^{-1/2}Y$, this means $Y_0 \sim 1/\xi$ as $\xi \to +\infty$. Only the $+$ branch matches (the $-$ branch grows like $-2\xi$). Furthermore, comparing $a/\xi$ with $1/\xi$ gives $a = 1$. Therefore: \begin{align*} Y_0(\xi) = -\xi + \sqrt{\xi^2 + 2}. \end{align*} **Step 5: Higher-order inner solution.** Substituting $Y \sim Y_0 + \epsilon^{1/2}Y_1 + \cdots$ into $(\xi + Y)Y' + Y = \epsilon^{1/2}$: \begin{align*} (\xi + Y_0 + \epsilon^{1/2}Y_1 + \cdots)(Y_0' + \epsilon^{1/2}Y_1' + \cdots) + (Y_0 + \epsilon^{1/2}Y_1 + \cdots) = \epsilon^{1/2}. \end{align*} Expanding the product: $(\xi + Y_0)Y_0' + \epsilon^{1/2}[(\xi + Y_0)Y_1' + Y_1 Y_0'] + O(\epsilon)$. Adding the remaining terms: \begin{align*} \underbrace{(\xi + Y_0)Y_0' + Y_0}_{= 0 \text{ (leading order)}} + \epsilon^{1/2}\underbrace{\bigl[(\xi + Y_0)Y_1' + Y_1 Y_0' + Y_1 - 1\bigr]}_{O(\epsilon^{1/2})} + \cdots = 0. \end{align*} At $O(\epsilon^{1/2})$: \begin{align*} \xi\,Y_1' + Y_0\,Y_1' + Y_1\,Y_0' + Y_1 = 1. \end{align*} Observe that $\xi\,Y_1' + Y_1 = \frac{d}{d\xi}(\xi\,Y_1)$ (product rule) and $Y_0\,Y_1' + Y_1\,Y_0' = \frac{d}{d\xi}(Y_0\,Y_1)$ (product rule). So the equation becomes $\frac{d}{d\xi}[(\xi + Y_0)\,Y_1] = 1$. Integrating: \begin{align*} (\xi + Y_0)\,Y_1 = \xi + b, \end{align*} where we used $\int 1\,d\xi = \xi$ and $b$ is an integration constant. Substituting $Y_0 = -\xi + \sqrt{\xi^2 + 2}$: \begin{align*} Y_1 = \frac{\xi + b}{\xi + Y_0} = \frac{\xi + b}{\sqrt{\xi^2 + 2}}. \end{align*} **Step 6: Van Dyke matching at $\epsilon^{1/2} : 1$.** The inner expansion contains two unknowns ($a$ and $b$) at orders $\epsilon^{-1/2}$ and $1$ respectively (since $y = \epsilon^{-1/2}Y_0 + Y_1 + \cdots$). To capture both, we match the inner expansion to order $1$ and the outer to order $\epsilon^{1/2}$ (which turns out to require only the leading outer term). *Inner $\to$ outer.* Define $W_2(x;\epsilon) := y_{\text{in}}^{(\text{ord}(1))}(x;\epsilon)$ — the inner expansion (expressed in terms of $y$, not $Y$) evaluated at $\xi = x/\epsilon^{1/2}$: \begin{align*} W_2(x;\epsilon) = \epsilon^{-1/2}\!\left(-\frac{x}{\epsilon^{1/2}} + \sqrt{\frac{x^2}{\epsilon} + 2a}\right) + \frac{x/\epsilon^{1/2} + b}{\sqrt{x^2/\epsilon + 2a}}. \end{align*} Expand $W_2$ for small $\epsilon$ with $x > 0$ fixed. For the $\epsilon^{-1/2}Y_0$ term, factor out $x/\epsilon^{1/2}$: \begin{align*} \epsilon^{-1/2}\!\left(-\frac{x}{\epsilon^{1/2}} + \frac{x}{\epsilon^{1/2}}\sqrt{1 + \frac{2a\epsilon}{x^2}}\right) &= \epsilon^{-1/2} \cdot \frac{x}{\epsilon^{1/2}}\!\left(-1 + 1 + \frac{a\epsilon}{x^2} + \cdots\right) \\ &= \epsilon^{-1/2} \cdot \frac{a\epsilon^{1/2}}{x} + O(\epsilon^{1/2}) = \frac{a}{x} + O(\epsilon^{1/2}). \end{align*} For the $Y_1$ term: substitute $\xi = x/\epsilon^{1/2}$ into $Y_1 = (\xi + b)/\sqrt{\xi^2 + 2a}$: \begin{align*} Y_1 = \frac{x/\epsilon^{1/2} + b}{\sqrt{x^2/\epsilon + 2a}} = \frac{x/\epsilon^{1/2} + b}{(x/\epsilon^{1/2})\sqrt{1 + 2a\epsilon/x^2}}. \end{align*} For small $\epsilon$ with $x$ fixed, $\sqrt{1 + 2a\epsilon/x^2} = 1 + O(\epsilon)$. The numerator is $x/\epsilon^{1/2} + b$ and the denominator is $x/\epsilon^{1/2} + O(\epsilon^{1/2})$, so: \begin{align*} Y_1 = \frac{x/\epsilon^{1/2} + b}{x/\epsilon^{1/2}(1 + O(\epsilon))} = \left(1 + \frac{b\epsilon^{1/2}}{x}\right)(1 + O(\epsilon)) = 1 + \frac{b\epsilon^{1/2}}{x} + O(\epsilon). \end{align*} Combining: $W_2(x;\epsilon) = \epsilon^{-1/2}Y_0 + Y_1 + \cdots = \frac{a}{x} + O(\epsilon^{1/2}) + 1 + \frac{b\epsilon^{1/2}}{x} + O(\epsilon)$. So $W_2$, expanded to order $\epsilon^{1/2}$: \begin{align*} \frac{a}{x} + 1 + \frac{\epsilon^{1/2}\,b}{x} + O(\epsilon). \end{align*} *Outer $\to$ inner.* Define $W_1(\xi;\epsilon) := y_{\text{out}}^{(\epsilon^{1/2})}(\epsilon^{1/2}\xi;\epsilon)$ — the outer expansion evaluated at $x = \epsilon^{1/2}\xi$. We only need the outer expansion to order $\epsilon^{1/2}$. The leading outer term is $y_0 = 1 + 1/x$ (order $1$). Is there a term at order $\epsilon^{1/2}$? Substituting $y \sim y_0 + \epsilon^{1/2}\,\tilde{y}$ into $(x + \epsilon y)y' + y = 1$ and collecting at $O(\epsilon^{1/2})$: \begin{align*} x\,\tilde{y}' + \tilde{y} = 0, \qquad \tilde{y}(1) = 0. \end{align*} This is the homogeneous equation $(x\,\tilde{y})' = 0$, giving $\tilde{y} = c/x$; the boundary condition forces $c = 0$, so $\tilde{y} = 0$. Therefore the outer expansion to order $\epsilon^{1/2}$ is simply $y_0 = 1 + 1/x$. So: \begin{align*} W_1(\xi;\epsilon) = 1 + \frac{1}{\epsilon^{1/2}\xi}. \end{align*} Expand $W_1$ for small $\epsilon$ with $\xi$ fixed, to order $1$: \begin{align*} W_1(\xi;\epsilon) = \epsilon^{-1/2}\frac{1}{\xi} + 1. \end{align*} (There is nothing to expand — $W_1$ is already an exact expression in $\epsilon$ and $\xi$.) *Set equal.* Van Dyke's principle requires $W_1 = W_2$ (to the required orders). Both $W_1$ and $W_2$ approximate $y$ in the overlap region. To compare, write $W_2$ in terms of $x$ and $W_1$ in terms of $\xi$, noting $x = \epsilon^{1/2}\xi$: \begin{align*} \frac{a}{x} + 1 + \frac{\epsilon^{1/2}\,b}{x} = \frac{1}{x} + 1. \end{align*} At order $1/x$: $a = 1$ (confirming our earlier finding). At order $\epsilon^{1/2}/x$: $b = 0$. All constants are determined. The final inner solution is: \begin{align*} Y_0(\xi) = -\xi + \sqrt{\xi^2 + 2}, \qquad Y_1(\xi) = \frac{\xi}{\sqrt{\xi^2 + 2}}. \end{align*} [/example] ## Worked Example: Logarithmic Switchback The previous examples involved inner solutions with algebraic growth (constants, linear terms). When the inner solution instead grows *logarithmically*, the matching procedure encounters a fundamental complication: logarithmic terms sit between algebraic orders and force corrections at unexpected places in the expansion. This phenomenon is called **switchback**. [example: Linear BVP with a Degenerate Operator] Consider the boundary-value problem \begin{align*} (\epsilon + x)\frac{d^2 y}{dx^2} + \frac{dy}{dx} = 1, \qquad y(0) = 0, \quad y(1) = 2, \qquad \epsilon \searrow 0. \end{align*} **Step 1: Outer expansion — a degenerate operator.** Setting $\epsilon = 0$ gives $xy_0'' + y_0' = 1$, which is *still second order* — the order has not been reduced. So why does a boundary layer form? The operator $x\,d^2/dx^2 + d/dx$ is **degenerate at $x = 0$**: its leading coefficient vanishes there. The general solution of $xy_0'' + y_0' = 1$ is found by writing the equation as $(xy_0')' = 1$ (since $xy_0'' + y_0' = \frac{d}{dx}(xy_0')$). Integrating once: $xy_0' = x + c_1$, so $y_0' = 1 + c_1/x$. Integrating again: \begin{align*} y_0 = x + c_1 \ln x + c_2. \end{align*} For $c_1 \neq 0$, the solution diverges as $x \to 0$. If we try to impose $y_0(0) = 0$, we need $c_1 = 0$ (to avoid the logarithmic singularity), leaving $y_0 = x + c_2$. But then $y_0(0) = 0$ gives $c_2 = 0$ and $y_0(1) = 2$ gives $c_2 = 1$ — a contradiction. The outer expansion cannot satisfy both boundary conditions. We therefore expect a boundary layer at $x = 0$ (where the operator degenerates). The layer width is $O(\epsilon)$, since the $\epsilon y''$ term regularises the degeneracy on that scale. The outer expansion inherits the condition at $x = 1$: \begin{align*} y_0(1) = 2 \implies 1 + c_1 \ln 1 + c_2 = 2 \implies c_2 = 1. \end{align*} Writing $a_0 := c_1$ (to be determined by matching): \begin{align*} y_0(x) = x + 1 + a_0 \ln x. \end{align*} **Step 2: Inner expansion.** Define $\xi = x/\epsilon$ and $Y(\xi;\epsilon) = y(\epsilon\xi;\epsilon)$. By the chain rule: $y' = (1/\epsilon)Y'$, $y'' = (1/\epsilon^2)Y''$. Substituting into $(\epsilon + x)y'' + y' = 1$: \begin{align*} (\epsilon + \epsilon\xi) \cdot \frac{1}{\epsilon^2}Y'' + \frac{1}{\epsilon}Y' &= 1 \\[4pt] \frac{\epsilon(1 + \xi)}{\epsilon^2}Y'' + \frac{1}{\epsilon}Y' &= 1 \\[4pt] \frac{1}{\epsilon}(1 + \xi)Y'' + \frac{1}{\epsilon}Y' &= 1. \end{align*} Multiplying through by $\epsilon$: \begin{align*} (1 + \xi)Y'' + Y' = \epsilon. \end{align*} We pose $Y(\xi;\epsilon) \sim Y_0(\xi) + \cdots$ with $Y_0(0) = 0$. At $O(1)$ (the right-hand side is $O(\epsilon)$ and drops out): \begin{align*} (1 + \xi)Y_0'' + Y_0' = 0. \end{align*} Write this as $\frac{d}{d\xi}[(1+\xi)Y_0'] = 0$ (since $(1+\xi)Y_0'' + Y_0' = \frac{d}{d\xi}[(1+\xi)Y_0']$ by the product rule). Integrating: $(1+\xi)Y_0' = C$, so $Y_0' = C/(1+\xi)$. Integrating again: \begin{align*} Y_0 = C\ln(1+\xi) + D. \end{align*} Applying $Y_0(0) = 0$: $D = 0$. Writing $A_0 := C$: \begin{align*} Y_0(\xi) = A_0 \ln(1 + \xi). \end{align*} **Step 3: Naive matching fails.** We attempt to match using an intermediate variable $\eta = x/\epsilon^\alpha$ with $0 < \alpha < 1$, defining $Z(\eta;\epsilon) = y(\epsilon^\alpha\eta;\epsilon)$. *$Z$ from the outer expansion.* Define $Z_{\text{out}}(\eta;\epsilon) := y_{\text{out}}^{(1)}(\epsilon^\alpha\eta;\epsilon)$. Substitute $x = \epsilon^\alpha\eta$: \begin{align*} Z_{\text{out}}(\eta;\epsilon) &= \epsilon^\alpha\eta + 1 + a_0\ln(\epsilon^\alpha\eta). \end{align*} Now $\ln(\epsilon^\alpha\eta) = \alpha\ln\epsilon + \ln\eta$. So: \begin{align*} Z_{\text{out}}(\eta;\epsilon) = \alpha\,a_0\ln\epsilon + a_0\ln\eta + 1 + O(\epsilon^\alpha). \end{align*} *$Z$ from the inner expansion.* Define $Z_{\text{in}}(\eta;\epsilon) := Y_{\text{in}}^{(1)}(\epsilon^{\alpha-1}\eta;\epsilon)$. Substitute $\xi = \epsilon^{\alpha-1}\eta$: \begin{align*} Z_{\text{in}}(\eta;\epsilon) &= A_0\ln(1 + \epsilon^{\alpha-1}\eta). \end{align*} Since $\alpha < 1$, we have $\alpha - 1 < 0$, so $\epsilon^{\alpha-1} \to +\infty$ as $\epsilon \to 0$. Therefore $1 + \epsilon^{\alpha-1}\eta \approx \epsilon^{\alpha-1}\eta$ for large values: \begin{align*} \ln(1 + \epsilon^{\alpha-1}\eta) = \ln(\epsilon^{\alpha-1}\eta) + \ln\!\left(1 + \frac{1}{\epsilon^{\alpha-1}\eta}\right) = (\alpha - 1)\ln\epsilon + \ln\eta + O(\epsilon^{1-\alpha}). \end{align*} So: \begin{align*} Z_{\text{in}}(\eta;\epsilon) = A_0(\alpha - 1)\ln\epsilon + A_0\ln\eta + O(\epsilon^{1-\alpha}). \end{align*} *Comparison.* Setting $Z_{\text{out}} = Z_{\text{in}}$: Matching $\ln\eta$ coefficients: $a_0 = A_0$. Matching $\ln\epsilon$ coefficients: $\alpha\,a_0 = (\alpha - 1)A_0$. Substituting $a_0 = A_0$ into the second equation: $\alpha\,A_0 = (\alpha - 1)A_0$, i.e. $A_0 = 0$. But $A_0 = 0$ makes the inner solution trivially zero — it cannot satisfy the boundary condition $y(0) = 0$ in a nontrivial way. **Naive matching has failed.** Matching the constant terms: $1 = 0$. This is also a contradiction. **Step 4: Diagnosing the problem — the ordering-matching circularity.** The contradiction in Step 3 is not an accident of this particular equation — it is a manifestation of the circularity described in the remark on logarithmic matching above. Here is exactly where the reasoning breaks down. The naive matching in Step 3 proceeds as follows: expand $Z_{\text{out}}$ and $Z_{\text{in}}$ as $\epsilon \to 0$ with $\eta$ fixed, then compare term by term. When we expand $Z_{\text{out}}$, the outer term $a_0\ln(\epsilon^\alpha\eta) = a_0\alpha\ln\epsilon + a_0\ln\eta$ splits into two pieces. The first piece, $a_0\alpha\ln\epsilon$, depends on $\epsilon$. The matching procedure then requires us to decide what asymptotic order this term belongs to — is it $O(1)$? $O(\ln\epsilon)$? $o(1)$? — and place it accordingly in the hierarchy. But we do not know $a_0$ yet; it is the very constant we are trying to determine by matching. **The asymptotic order of $a_0\ln\epsilon$ depends on the size of $a_0$:** - If $a_0 = O(1)$, then $a_0\ln\epsilon \to -\infty$ as $\epsilon \to 0$. It is *larger* than any $O(1)$ term and must be balanced by something equally large on the inner side. The matching equations in Step 3 implicitly made this assumption: they treated $a_0\ln\epsilon$ and $A_0\ln\epsilon$ as the dominant terms (matching $\ln\epsilon$ coefficients separately from $O(1)$ terms). The result — $A_0 = 0$ and the contradiction $1 = 0$ — shows that this assumption is inconsistent. - If $a_0 = O(1/\ln\epsilon)$, then $a_0\ln\epsilon = O(1)$. The $\ln\epsilon$ factor is absorbed into the constant, and the term is genuinely $O(1)$ — it does not blow up. In this case, $a_0\ln\epsilon$ should *not* be separated from the other $O(1)$ terms. It sits at the same algebraic order and must be matched together with them. The correct answer turns out to be the second case. The inner solution grows as $A_0\ln\xi$ for large $\xi$, and for this to match the bounded outer solution, $A_0$ must compensate for the logarithmic growth by being logarithmically small: $A_0 \sim -1/\ln\epsilon$. Since $a_0 = A_0$ (from matching the $\ln\eta$ coefficients, which is valid regardless), we have $a_0 = -1/\ln\epsilon$. When we write $a_0 = A_0 = -1/\ln\epsilon$ and look at the outer expansion $y_0 = x + 1 + a_0\ln x = x + 1 - (\ln x)/\ln\epsilon$, the last term is $O(1/\ln\epsilon)$ — a logarithmic correction to an $O(1)$ expansion. This "unexpected" term at an intermediate order is the **switchback**: matching at one order has forced a correction at a nominally higher order. **Step 5: The fix — generalised expansions.** Rather than tracking every logarithmic sub-order, we group terms by algebraic power and allow the coefficients to depend on $\ln\epsilon$. We pose: \begin{align*} y(x;\epsilon) &\sim y_0(x;\epsilon) + \epsilon\,y_1(x;\epsilon), \\[4pt] Y(\xi;\epsilon) &\sim Y_0(\xi;\epsilon) + \epsilon\,Y_1(\xi;\epsilon), \end{align*} where each term may depend on $\ln\epsilon$ (but not on powers of $\epsilon$). *Leading outer:* the equation is unchanged: $xy_0'' + y_0' = 1$, $y_0(1;\epsilon) = 2$. The solution is: \begin{align*} y_0(x;\epsilon) = x + 1 + a_0(\epsilon)\ln x, \end{align*} where $a_0(\epsilon)$ is now allowed to depend on $\ln\epsilon$. *$O(\epsilon)$ outer:* substituting $y \sim y_0 + \epsilon y_1$ into $(\epsilon + x)y'' + y' = 1$ and collecting at $O(\epsilon)$: the $\epsilon y_0''$ term comes down from the $(\epsilon + x)y''$ piece. Writing $(\epsilon + x)y'' = xy'' + \epsilon y''$, the full equation is: \begin{align*} x(y_0'' + \epsilon y_1'' + \cdots) + \epsilon(y_0'' + \cdots) + (y_0' + \epsilon y_1' + \cdots) = 1. \end{align*} At $O(1)$: $xy_0'' + y_0' = 1$ (already solved). At $O(\epsilon)$: \begin{align*} xy_1'' + y_1' = -y_0'', \qquad y_1(1;\epsilon) = 0. \end{align*} Computing $y_0' = 1 + a_0/x$, so $y_0'' = -a_0/x^2$. Thus $xy_1'' + y_1' = a_0/x^2$. Writing this as $(xy_1')' = a_0/x^2$ and integrating: $xy_1' = -a_0/x + c$, so $y_1' = -a_0/x^2 + c/x$. Integrating again: \begin{align*} y_1 = \frac{a_0}{x} + c\ln x + d. \end{align*} Applying $y_1(1;\epsilon) = 0$: $a_0 + 0 + d = 0$, so $d = -a_0$. Writing $a_1(\epsilon) := c$: \begin{align*} y_1(x;\epsilon) = a_0(\epsilon)\!\left(\frac{1}{x} - 1\right) + a_1(\epsilon)\ln x. \end{align*} *Leading inner:* the equation is $(1+\xi)Y_0'' + Y_0' = 0$, $Y_0(0;\epsilon) = 0$, with solution (as before): \begin{align*} Y_0(\xi;\epsilon) = A_0(\epsilon)\ln(1 + \xi). \end{align*} *$O(\epsilon)$ inner:* substituting $Y \sim Y_0 + \epsilon Y_1$ into $(1+\xi)Y'' + Y' = \epsilon$: \begin{align*} (1+\xi)(Y_0'' + \epsilon Y_1'' + \cdots) + (Y_0' + \epsilon Y_1' + \cdots) = \epsilon. \end{align*} At $O(1)$: $(1+\xi)Y_0'' + Y_0' = 0$ (already solved). At $O(\epsilon)$: \begin{align*} (1+\xi)Y_1'' + Y_1' = 1, \qquad Y_1(0;\epsilon) = 0. \end{align*} Writing $(1+\xi)Y_1'' + Y_1' = [(1+\xi)Y_1']'$ and integrating: $(1+\xi)Y_1' = \xi + C_1$. So $Y_1' = (\xi + C_1)/(1+\xi) = 1 - (1 - C_1)/(1+\xi)$. Integrating: \begin{align*} Y_1 = \xi - (1 - C_1)\ln(1+\xi) + C_2. \end{align*} Applying $Y_1(0;\epsilon) = 0$: $0 - 0 + C_2 = 0$, so $C_2 = 0$. Writing $A_1(\epsilon) := C_1$, we have $-(1 - C_1) = A_1 - 1$. The solution becomes: \begin{align*} Y_1(\xi;\epsilon) = A_1(\epsilon)\ln(1+\xi) + \xi - \ln(1+\xi), \end{align*} where $A_1(\epsilon)$ is a free constant (allowed to depend on $\ln\epsilon$). *Verification:* $Y_1' = A_1/(1+\xi) + 1 - 1/(1+\xi) = (A_1 - 1)/(1+\xi) + 1$. Then $(1+\xi)Y_1' = A_1 - 1 + (1+\xi) = \xi + A_1$. Differentiating: $[(1+\xi)Y_1']' = 1$. ✓ Also $Y_1(0) = 0$. ✓ So the two-term inner expansion is: \begin{align*} Y(\xi;\epsilon) \sim A_0(\epsilon)\ln(1+\xi) + \epsilon\bigl[A_1(\epsilon)\ln(1+\xi) + \xi - \ln(1+\xi)\bigr]. \end{align*} **Step 6: Van Dyke matching at $\epsilon : \epsilon$ (up to logarithmic factors).** *Direction 1 (Outer $\to$ inner).* Define $W_1(\xi;\epsilon) := y_{\text{out}}^{(\epsilon)}(\epsilon\xi;\epsilon)$: \begin{align*} W_1(\xi;\epsilon) &= \epsilon\xi + 1 + a_0\ln(\epsilon\xi) + \epsilon\!\left[a_0\!\left(\frac{1}{\epsilon\xi} - 1\right) + a_1\ln(\epsilon\xi)\right]. \end{align*} Now $\ln(\epsilon\xi) = \ln\epsilon + \ln\xi$. Expanding term by term: *$y_0$ piece:* $\epsilon\xi + 1 + a_0(\ln\epsilon + \ln\xi) = \epsilon\xi + 1 + a_0\ln\epsilon + a_0\ln\xi$. *$\epsilon y_1$ piece:* $\epsilon \cdot a_0/({\epsilon\xi}) - \epsilon a_0 + \epsilon a_1(\ln\epsilon + \ln\xi) = a_0/\xi - \epsilon a_0 + \epsilon a_1\ln\epsilon + \epsilon a_1\ln\xi$. Combining, and keeping all terms to algebraic order $\epsilon$ (including logarithmic factors): \begin{align*} W_1(\xi;\epsilon) = \epsilon\xi + (a_0 + \epsilon a_1)\ln\xi + \frac{a_0}{\xi} + (1 + a_0\ln\epsilon - \epsilon a_0 + \epsilon a_1\ln\epsilon). \end{align*} *Direction 2 (Inner $\to$ outer).* Define $W_2(x;\epsilon) := Y_{\text{in}}^{(\epsilon)}(x/\epsilon;\epsilon)$: \begin{align*} W_2(x;\epsilon) = A_0\ln(1 + x/\epsilon) + \epsilon\bigl[A_1\ln(1 + x/\epsilon) + x/\epsilon - \ln(1 + x/\epsilon)\bigr]. \end{align*} For $x > 0$ fixed and $\epsilon \to 0$: $x/\epsilon \to +\infty$, so $\ln(1 + x/\epsilon) = \ln(x/\epsilon) + \ln(1 + \epsilon/x)$. Expanding $\ln(1 + \epsilon/x) = \epsilon/x + O(\epsilon^2)$: \begin{align*} \ln(1 + x/\epsilon) = \ln(x/\epsilon) + \frac{\epsilon}{x} + O(\epsilon^2) = -\ln\frac{\epsilon}{x} + \frac{\epsilon}{x} + O(\epsilon^2). \end{align*} And $\ln(x/\epsilon) = \ln x - \ln\epsilon$, so $\ln(1 + x/\epsilon) = -\ln\epsilon + \ln x + \epsilon/x + O(\epsilon^2)$. Substituting into each piece of $W_2$: *$A_0$ piece:* $A_0(-\ln\epsilon + \ln x + \epsilon/x + \cdots)$. *$\epsilon A_1$ piece:* $\epsilon A_1(-\ln\epsilon + \ln x + \cdots)$. *$\epsilon \cdot x/\epsilon$ piece:* $x$. *$-\epsilon\ln(1+x/\epsilon)$ piece:* $-\epsilon(-\ln\epsilon + \ln x + \cdots) = \epsilon\ln\epsilon - \epsilon\ln x + \cdots$. Combining: \begin{align*} W_2(x;\epsilon) &= -A_0\ln\epsilon + A_0\ln x + \frac{A_0\epsilon}{x} + \epsilon A_1(-\ln\epsilon + \ln x) + x + \epsilon\ln\epsilon - \epsilon\ln x + O(\epsilon^2). \end{align*} Now re-express in terms of $\xi = x/\epsilon$ (since the final comparison should be in a common variable). Using $x = \epsilon\xi$ and $\ln x = \ln\epsilon + \ln\xi$: \begin{align*} W_2 &= -A_0\ln\epsilon + A_0(\ln\epsilon + \ln\xi) + \frac{A_0}{\xi} + \epsilon A_1(-\ln\epsilon + \ln\epsilon + \ln\xi) + \epsilon\xi + \epsilon\ln\epsilon - \epsilon(\ln\epsilon + \ln\xi) + \cdots \\ &= A_0\ln\xi + \frac{A_0}{\xi} + \epsilon A_1\ln\xi + \epsilon\xi - \epsilon\ln\xi + \cdots \\ &= \epsilon\xi + (A_0 + \epsilon A_1 - \epsilon)\ln\xi + \frac{A_0}{\xi} + \cdots \end{align*} (Note how the $\ln\epsilon$ terms cancelled in the $A_0$ piece — this is the key simplification that makes generalised expansions work.) *Comparison.* Setting $W_1 = W_2$ and comparing term by term: *Coefficient of $\epsilon\xi$:* $1 = 1$. ✓ *Coefficient of $\ln\xi$:* $a_0 + \epsilon a_1 = A_0 + \epsilon A_1 - \epsilon$. *Coefficient of $1/\xi$:* $a_0 = A_0$. *Remaining constant:* $1 + a_0\ln\epsilon - \epsilon a_0 + \epsilon a_1\ln\epsilon = 0$ (since the inner side has no $\xi$-independent constant term at this order). From the $1/\xi$ coefficient: $a_0 = A_0$. From the $\ln\xi$ coefficient and $a_0 = A_0$: $\epsilon a_1 = \epsilon A_1 - \epsilon$, so $a_1 = A_1 - 1$. From the constant term: $1 + a_0\ln\epsilon - \epsilon a_0 + \epsilon a_1\ln\epsilon = 0$. This equation involves terms at multiple logarithmic sub-orders, so we solve it iteratively: *At leading order* (keeping only terms that are $O(1)$ or larger — the $\epsilon a_0$ and $\epsilon a_1\ln\epsilon$ terms are both algebraically $O(\epsilon)$ and drop out): \begin{align*} 1 + a_0\ln\epsilon = 0 \implies a_0 = A_0 = -\frac{1}{\ln\epsilon}. \end{align*} *At order $\epsilon\ln\epsilon$:* substitute $a_0 = -1/\ln\epsilon$ and keep the next-largest terms. The equation becomes $1 - 1 + \epsilon a_1\ln\epsilon + \epsilon/\ln\epsilon = 0$. The $\epsilon/\ln\epsilon$ term is algebraically smaller than $\epsilon\ln\epsilon$, so at this sub-order: \begin{align*} \epsilon a_1\ln\epsilon = 0 \implies a_1 = 0, \qquad A_1 = a_1 + 1 = 1. \end{align*} *At order $\epsilon$:* now keep the $\epsilon/\ln\epsilon$ term. The equation becomes $\epsilon a_1\ln\epsilon = -\epsilon/\ln\epsilon$, giving: \begin{align*} a_1 = -\frac{1}{\ln^2\epsilon}, \qquad A_1 = 1 - \frac{1}{\ln^2\epsilon}. \end{align*} This iterative extraction of logarithmic corrections is characteristic of switchback problems. **Step 7: Final result.** Substituting back: \begin{align*} y_{\text{out}}(x;\epsilon) &\sim 1 + x - \frac{\ln x}{\ln\epsilon} + \epsilon\!\left[-\frac{1}{\ln\epsilon}\!\left(\frac{1}{x} - 1\right) - \frac{\ln x}{\ln^2\epsilon}\right], \\[8pt] Y_{\text{in}}(\xi;\epsilon) &\sim -\frac{\ln(1+\xi)}{\ln\epsilon} + \epsilon\!\left[\xi - \frac{\ln(1+\xi)}{\ln^2\epsilon}\right]. \end{align*} The characteristic signature of logarithmic switchback is visible: the coefficients involve **inverse powers of $\ln\epsilon$** rather than powers of $\epsilon$. The outer expansion contains the term $-(\ln x)/\ln\epsilon$, which is technically $o(1)$ but larger than $\epsilon$ — it sits at the intermediate logarithmic order that forced the switchback. [/example] ## Worked Example: Logarithmic Switchback from Nonlinearity The switchback example above involved a degenerate operator — the coefficient of $y''$ vanished at $x = 0$, producing a logarithmic inner solution. But logarithmic complications can also arise from nonlinearity, even when the leading-order outer solution is perfectly regular. The following example demonstrates this second mechanism. [example: Nonlinear BVP with Logarithmic Correction] Consider the nonlinear boundary-value problem \begin{align*} x^2 \frac{d^2 y}{dx^2} - \epsilon \, y \frac{dy}{dx} = 0, \qquad y(0) = 1, \quad y(1) = 0, \qquad \epsilon \searrow 0. \end{align*} **Step 1: Leading-order outer expansion.** We pose $y(x;\epsilon) \sim y_0(x) + \epsilon\,y_1(x) + \cdots$ with $x$ fixed. At $O(1)$: \begin{align*} x^2 y_0'' = 0. \end{align*} Since $x > 0$ in the outer region, this gives $y_0'' = 0$, so $y_0$ is linear. Both boundary conditions can be applied: \begin{align*} y_0(0) = 1, \quad y_0(1) = 0 \implies y_0(x) = 1 - x. \end{align*} The leading-order outer solution is completely regular — no logarithmic singularity, no degenerate operator. The expansion appears to be standard. **Step 2: The $O(\epsilon)$ outer correction introduces logarithms.** At $O(\epsilon)$, the equation is \begin{align*} x^2 y_1'' = y_0\,y_0', \qquad y_1(0) = 0, \quad y_1(1) = 0. \end{align*} Computing: $y_0 = 1 - x$ and $y_0' = -1$, so $y_0\,y_0' = -(1-x) = x - 1$. The equation becomes \begin{align*} y_1'' = \frac{x - 1}{x^2} = \frac{1}{x} - \frac{1}{x^2}. \end{align*} Integrating once: $y_1' = \ln x + 1/x + C_1$. Integrating again: \begin{align*} y_1 = x\ln x - x + \ln x + C_1 x + C_2. \end{align*} The boundary condition $y_1(1) = 0$ gives $0 - 1 + 0 + C_1 + C_2 = 0$, so $C_2 = 1 - C_1$. Writing $a_1 := C_1$: \begin{align*} y_1(x) = a_1 x + 1 - a_1 + x\ln x + \ln x - x. \end{align*} Now try $y_1(0) = 0$: as $x \to 0$, $x\ln x \to 0$ but $\ln x \to -\infty$. The term $\ln x$ diverges — **the $O(\epsilon)$ correction cannot satisfy the boundary condition at $x = 0$**. This divergence signals a boundary layer near $x = 0$, even though the leading-order solution was perfectly smooth there. **Step 3: Boundary layer scaling.** The logarithmic divergence in $y_1$ tells us that the outer expansion breaks down near $x = 0$. With $y = O(1)$ (since $y(0) = 1$), a scaling analysis of $x^2 y'' - \epsilon yy' = 0$ suggests: if $x = O(\delta)$, then $x^2 y'' = O(\delta^2/\delta^2) = O(1)$ and $\epsilon yy' = O(\epsilon/\delta)$. Balancing requires $\epsilon/\delta \sim 1$, giving $\delta = \epsilon$. Define $\xi = x/\epsilon$ and $Y(\xi;\epsilon) = y(\epsilon\xi;\epsilon)$. By the chain rule: $y' = (1/\epsilon)Y'$, $y'' = (1/\epsilon^2)Y''$. Substituting: \begin{align*} (\epsilon\xi)^2 \cdot \frac{1}{\epsilon^2}Y'' - \epsilon \cdot Y \cdot \frac{1}{\epsilon}Y' = 0, \end{align*} which simplifies to \begin{align*} \xi^2 Y'' - YY' = 0, \qquad Y(0;\epsilon) = 1. \end{align*} This is the full nonlinear equation — no simplification has occurred. However, the key simplification comes from the *size* of the inner solution. Both the boundary value $y(0) = 1$ and the limit of the leading-order outer solution as $x \to 0$ are equal to $1$. This means $Y = 1 + O(\epsilon)$ in the inner region, suggesting the inner expansion \begin{align*} Y(\xi;\epsilon) \sim 1 + \epsilon\,Y_1(\xi;\epsilon) + \cdots \quad \text{as } \epsilon \to 0, \; \xi \text{ fixed}, \end{align*} where $Y_1$ is allowed to depend logarithmically on $\epsilon$. **Step 4: Leading-order inner correction.** Substituting $Y = 1 + \epsilon Y_1 + \cdots$ into $\xi^2 Y'' - YY' = 0$: \begin{align*} \xi^2(\epsilon Y_1'' + \cdots) - (1 + \epsilon Y_1 + \cdots)(\epsilon Y_1' + \cdots) = 0. \end{align*} At $O(\epsilon)$ (the $O(1)$ equation is trivially satisfied since the leading term is constant): \begin{align*} \xi^2 Y_1'' - Y_1' = 0, \qquad Y_1(0;\epsilon) = 0. \end{align*} This is a linear second-order ODE. To solve it, write $Y_1' = W$, so $\xi^2 W' - W = 0$, i.e. $W'/W = 1/\xi^2$. Integrating: \begin{align*} \ln|W| = -\frac{1}{\xi} + \text{const} \implies W = Y_1' = C\,e^{-1/\xi}. \end{align*} Integrating $Y_1' = C\,e^{-1/\xi}$ from $\xi = 0$ (where $Y_1 = 0$) is delicate because $e^{-1/\xi} \to 0$ as $\xi \to 0^+$ (faster than any power). We write the solution as \begin{align*} Y_1(\xi;\epsilon) = A(\epsilon)\int_0^\xi C\,e^{-1/s}\,ds, \end{align*} but it is more convenient to express the solution using integration from $\infty$: \begin{align*} Y_1(\xi;\epsilon) = A(\epsilon)\int_\infty^{1/\xi} u^{-2}\,e^{-u}\,du = -A(\epsilon)\int_{1/\xi}^\infty u^{-2}\,e^{-u}\,du, \end{align*} where the substitution $u = 1/s$ gives $ds = -du/u^2$ and $e^{-1/s} = e^{-u}$. (The sign and limits are chosen so that the boundary condition $Y_1(0) = 0$ is automatic: as $\xi \to 0^+$, $1/\xi \to +\infty$ and $\int_{1/\xi}^\infty u^{-2}e^{-u}\,du \to 0$.) The constant $A(\epsilon)$ is allowed to depend logarithmically on $\epsilon$ and will be determined by matching. **Step 5: Asymptotics of the inner solution for matching.** For matching, we need the behaviour of $Y_1$ as $\xi \to +\infty$ (or more precisely, the $\epsilon \to 0$ limit of the composition $Y_1(\delta^{\alpha-1}\eta)$ for $\eta > 0$, which probes the large-$\xi$ behaviour of $Y_1$). This requires the expansion of \begin{align*} I(\rho) := \int_\rho^\infty u^{-2}\,e^{-u}\,du \quad \text{as } \rho \to 0^+. \end{align*} We derive this expansion using [integration by parts](/theorems/210). Start by splitting the [integral](/page/Integral) at $u = 1$ (or any convenient fixed point): \begin{align*} I(\rho) = \int_\rho^1 u^{-2}\,e^{-u}\,du + \underbrace{\int_1^\infty u^{-2}\,e^{-u}\,du}_{\text{a fixed constant}}. \end{align*} The second integral is a numerical constant. For the first integral, expand $e^{-u} = 1 - u + \frac{1}{2}u^2 - \cdots$ for $u$ near $0$: \begin{align*} \int_\rho^1 u^{-2}\,e^{-u}\,du &= \int_\rho^1 u^{-2}(1 - u + \tfrac{1}{2}u^2 - \cdots)\,du \\ &= \int_\rho^1 (u^{-2} - u^{-1} + \tfrac{1}{2} - \cdots)\,du \\ &= \left[-\frac{1}{u} - \ln u + \frac{1}{2}u - \cdots\right]_\rho^1 \\ &= \left(-1 - 0 + \frac{1}{2} - \cdots\right) - \left(-\frac{1}{\rho} - \ln\rho + \frac{1}{2}\rho - \cdots\right) \\ &= \frac{1}{\rho} + \ln\rho + \left(-1 + \frac{1}{2} + \int_1^\infty \cdots\right) + O(\rho). \end{align*} Combining with the second integral and collecting all numerical constants into a single quantity, we obtain \begin{align*} I(\rho) = \frac{1}{\rho} + \ln\rho + (\gamma - 1) + O(\rho) \quad \text{as } \rho \to 0^+, \end{align*} where $\gamma = 0.5772\ldots$ is the Euler–Mascheroni constant. (The value $\gamma - 1$ can be verified by computing $I(\rho)$ numerically for small $\rho$ and comparing.) Setting $\rho = 1/\xi$ (so $\rho \to 0$ as $\xi \to \infty$): \begin{align*} \int_{1/\xi}^\infty u^{-2}\,e^{-u}\,du = \xi + \ln\frac{1}{\xi} + (\gamma - 1) + O(1/\xi) = \xi - \ln\xi + (\gamma - 1) + O(1/\xi). \end{align*} Therefore the inner solution has the large-$\xi$ behaviour: \begin{align*} Y_1(\xi;\epsilon) \sim -A(\epsilon)\left[\xi - \ln\xi + (\gamma - 1)\right] \quad \text{as the matching limit is taken.} \end{align*} **Step 6: Generalised outer expansion.** Since $y_1$ contains $\ln x$, we employ the generalised approach: allow $a_1(\epsilon)$ to depend logarithmically on $\epsilon$. The outer expansion to algebraic order $\epsilon$ is \begin{align*} y(x;\epsilon) \sim (1 - x) + \epsilon\left[a_1(\epsilon)\,x + 1 - a_1(\epsilon) + x\ln x + \ln x - x\right] \quad \text{as } \epsilon \to 0, \; x \text{ fixed.} \end{align*} **Step 7: Van Dyke matching at $\epsilon : \epsilon$ (up to logarithmic factors).** *Direction 1 (Outer $\to$ inner).* Evaluate the outer expansion at $x = \epsilon\xi$: \begin{align*} y_{\text{out}}^{(\epsilon)} &= (1 - \epsilon\xi) + \epsilon\left[a_1\,\epsilon\xi + 1 - a_1 + \epsilon\xi\ln(\epsilon\xi) + \ln(\epsilon\xi) - \epsilon\xi\right]. \end{align*} Expand as $\epsilon \to 0$ with $\xi$ fixed. Note $\ln(\epsilon\xi) = \ln\epsilon + \ln\xi$. The terms $\epsilon\xi\ln(\epsilon\xi) = \epsilon\xi\ln\epsilon + \epsilon\xi\ln\xi$ and $a_1\,\epsilon\xi$ are all $O(\epsilon)$ or smaller (they involve $\epsilon\xi$ times something). Keeping terms to algebraic order $\epsilon$: - From $(1 - \epsilon\xi)$: $1 - \epsilon\xi$. - From $\epsilon(1 - a_1)$: $\epsilon(1 - a_1)$. - From $\epsilon\ln(\epsilon\xi) = \epsilon\ln\epsilon + \epsilon\ln\xi$: these are $O(\epsilon\ln\epsilon)$ and $O(\epsilon)$. Combining all terms to algebraic order $\epsilon$ (grouping by powers of $\xi$): \begin{align*} W_1(\xi;\epsilon) = 1 - \epsilon\xi + \epsilon\left[1 - a_1(\epsilon) + \ln\epsilon + \ln\xi\right] + O(\epsilon^2\ln\epsilon). \end{align*} *Direction 2 (Inner $\to$ outer).* The inner expansion is \begin{align*} Y_{\text{in}}^{(\epsilon)} = 1 - \epsilon\,A(\epsilon)\int_{1/\xi}^\infty u^{-2}\,e^{-u}\,du. \end{align*} Evaluate at $\xi = x/\epsilon$, so $1/\xi = \epsilon/x$: \begin{align*} Y_{\text{in}}^{(\epsilon)} = 1 - \epsilon\,A(\epsilon)\int_{\epsilon/x}^\infty u^{-2}\,e^{-u}\,du. \end{align*} Now expand as $\epsilon \to 0$ with $x$ fixed, using the asymptotic expansion of $I(\rho)$ at $\rho = \epsilon/x \to 0^+$: \begin{align*} \int_{\epsilon/x}^\infty u^{-2}\,e^{-u}\,du = \frac{x}{\epsilon} + \ln\frac{\epsilon}{x} + (\gamma - 1) + O(\epsilon/x). \end{align*} Substituting $\ln(\epsilon/x) = \ln\epsilon - \ln x$: \begin{align*} Y_{\text{in}}^{(\epsilon)} &= 1 - \epsilon\,A\left[\frac{x}{\epsilon} + \ln\epsilon - \ln x + (\gamma - 1) + O(\epsilon)\right] \\ &= 1 - Ax - \epsilon A\ln\epsilon + \epsilon A\ln x - \epsilon A(\gamma - 1) + O(\epsilon^2). \end{align*} Re-express in terms of $\xi = x/\epsilon$ (using $x = \epsilon\xi$, $\ln x = \ln\epsilon + \ln\xi$): \begin{align*} W_2(\xi;\epsilon) &= 1 - \epsilon A\xi - \epsilon A\ln\epsilon + \epsilon A(\ln\epsilon + \ln\xi) - \epsilon A(\gamma - 1) + O(\epsilon^2\ln\epsilon) \\ &= 1 - \epsilon A\xi + \epsilon A\ln\xi - \epsilon A(\gamma - 1) + O(\epsilon^2\ln\epsilon). \end{align*} Note the key cancellation: $-\epsilon A\ln\epsilon + \epsilon A\ln\epsilon = 0$. The $\ln\epsilon$ terms cancel completely — this is characteristic of well-matched generalised expansions. *Comparison.* Setting $W_1 = W_2$ and comparing term by term: *Constant term:* $1 = 1$. ✓ *Coefficient of $\epsilon\xi$:* $-1 = -A$, giving \begin{align*} A(\epsilon) = 1. \end{align*} Note that in the end, $A$ does not depend on $\ln\epsilon$ — it is a genuine constant. *Coefficient of $\epsilon\ln\xi$:* $1 = A = 1$. ✓ (This is automatically consistent.) *$\xi$-independent terms at $O(\epsilon)$:* $1 - a_1(\epsilon) + \ln\epsilon = -(\gamma - 1) = 1 - \gamma$, giving \begin{align*} a_1(\epsilon) = \ln\epsilon + \gamma. \end{align*} Here $a_1$ *does* depend on $\ln\epsilon$ — this is the switchback. The logarithmic dependence was inherited from the nonlinear term $y_0 y_0' = -(1-x)$, which produced $\ln x$ in $y_1$, which upon rescaling generated $\ln\epsilon$. **Step 8: Final result.** Substituting $a_1(\epsilon) = \ln\epsilon + \gamma$ into the outer expansion: \begin{align*} y_{\text{out}}(x;\epsilon) \sim (1 - x) + \epsilon(\ln\epsilon + \gamma)\,x + \epsilon(1 - \ln\epsilon - \gamma) + \epsilon(x\ln x + \ln x - x). \end{align*} Collecting terms by algebraic sub-order: \begin{align*} y_{\text{out}}(x;\epsilon) \sim (1 - x) + \epsilon\ln\epsilon \cdot (x - 1) + \epsilon\left[(x+1)\ln x + (\gamma - 1)x + (1 - \gamma)\right]. \end{align*} The inner expansion is \begin{align*} Y_{\text{in}}(\xi;\epsilon) \sim 1 - \epsilon\int_{1/\xi}^\infty u^{-2}\,e^{-u}\,du. \end{align*} Several features are worth noting. First, the inner solution has $A = 1$ — a genuine constant, not logarithmically dependent on $\epsilon$. The switchback appears only in the outer expansion, through $a_1 = \ln\epsilon + \gamma$. Second, the term $\epsilon\ln\epsilon\,(x-1)$ in the outer expansion is larger than $O(\epsilon)$ but smaller than $O(1)$ — it is the intermediate-order correction forced by the logarithmic matching. Third, the inner solution involves a special function (the exponential integral variant $\int_{1/\xi}^\infty u^{-2}e^{-u}\,du$), whose small-argument asymptotics required the Euler–Mascheroni constant $\gamma$. This is typical of logarithmic problems: the matching constants often involve transcendental numbers that arise from asymptotic expansions of special functions. [/example] ## References 1. M. Van Dyke, *Perturbation Methods in Fluid Mechanics*, annotated edition, Parabolic Press, 1975. 2. E. J. Hinch, *Perturbation Methods*, Cambridge University Press, 1991. 3. C. M. Bender and S. A. Orszag, *Advanced Mathematical Methods for Scientists and Engineers*, Springer, 1999. 4. O. Schnitzer, *Asymptotic Analysis: Matched Asymptotics and Logarithms*, Imperial College London lecture notes, 2025.