I would like to find the first two terms of the asymptotic approximation for $\epsilon^2 y'' = axy$ on $0\leq x <
\infty$ where $y(\infty) = 0$ and $a>0$.
Here is my work so far, with the standard WKB method, we plug in
$$y(x) = \exp\left[\frac{1}{\delta} \sum_{n=0}^\infty \delta^n S_n(x)\right]$$ and divide by it; our ODE simplifies to
$$\epsilon^2\left[\frac{1}{\delta} S_0'' + S_1 + \cdots + \frac{1}{\delta^2} S_0'^2 + 2 \frac{1}{\delta} S_0'S_1' + \cdots \right] -ax = 0 \quad \quad (\star)$$
Our largest term $\frac{\epsilon^2}{\delta^2} S_0'^2 $must match with $-ax$, so we make $\delta = \epsilon$, and get
$$S_0'^2 \sim ax$$
thus
$$S_0 \sim \pm \sqrt{a} \frac{2}{3} x^{3/2}.$$
And by matching $O(\epsilon)$ terms in $(\star)$, we get
$$2S_0'S_1' + S_0'' = 0$$
and we get $S_1 \sim – \frac{1}{4} \log(ax)$.
So our approximation
$$y (x)\sim C_1 (ax)^{-1/4} \exp \left[\frac{1}{\epsilon} \sqrt{a} \frac{2}{3} x^{3/2}\right] + C_2 (ax)^{-1/4} \exp \left[-\frac{1}{\epsilon} \sqrt{a} \frac{2}{3} x^{3/2}\right]$$
And from our condition $y(\infty) = 0$, we get $C_1 = 0$.
But the above approximation does not make sense at $x=0$, how can we take care of this problem?
I know boundary layer theory but often it has problems when the coefficient of $y'$ is constant zero. I also worked out the WKB method for a boundary layer problem here, I am not sure how would I apply to this case.
Edit: If I add a boundary layer at $x=0$, we let $\xi = x/\delta$, then the ODE in terms of delta will become
$$\frac{\epsilon^2}{\delta^2} (y_0'' + \epsilon y_1'' + \cdots ) = a\xi\delta(y_0 + \epsilon y_1, + \cdots)$$
which gives
$$\frac{\epsilon^2}{\delta^2} y_0'' + \frac{\epsilon^3}{\delta^2} y_1'' + \cdots – \delta a\xi y_0 + \cdots = 0$$
If I let $\delta = \epsilon$, I will get $y_0 = a\xi + b$, how would I match this with $y(x)$ from WKB method
Best Answer
Your asymptotic approximation is correct : $$y (x)\sim C_2 (ax)^{-1/4} \exp \left[-\frac{1}{\epsilon} \sqrt{a} \frac{2}{3} x^{3/2}\right] \tag 1$$ according to the boundary condition $y(\infty)=0$.
Note : A second order ODE requires at least two boundary conditions to expect having a unique solution. Since the wording of the question specifies only the condition $y(\infty)=0$, the solution remains undetermined : The coefficient $C_2$ cannot be computed.
Another key point is that an asymptotic formula is valid only for $x\to\infty$, not for $x\to 0$.So, don't expect Eq.$(1)$ makes sens at $x=0$. It is normal that the asymptotic approximation be far to be correct for small values of $x$.
If you want an approximate around $x=0$, expand $y(x)$ in power series. $$y(x)\simeq c_0+c_1x+c_2x^2+...$$ Of course, this formula is not valid for large $x$, so the boundary condition at infinity cannot be fulfilled. Moreover, two boundary conditions should be required close to $x=0$ in order to determine the coefficients $c_0$ , $c_1$, $c_2$ , …
In fact, as Claude Leibovici quite rightly wrote in comments, the exact solving involves the Airy functions. $$y(x)=c_1\text{Ai}\left(\sqrt[3]{\frac{a}{\epsilon^2}}\:x\right)+c_2\text{Bi}\left(\sqrt[3]{\frac{a}{\epsilon^2}}\:x\right)$$ Ai and Bi are the Airy functions http://mathworld.wolfram.com/AiryFunctions.html
The condition $y(\infty)=0$ is satisfied any $c_1$ and $c_2$, which shows that the condition $y(\infty)=0$ is not sufficient to fully determine a solution.
On physicists viewpoint, defining the properties of a boundary layer is certainly a way to introduce some boundary conditions in the range around $x=0$. Supposing that $y(x)\simeq y_0+y'_0 x$ in a layer of a given thickness $\delta$ allows to compute $c_1$ and $c_2$ as a function of $y_0$ , $y'_0$ and $\delta$, thanks to the series expansion of the Airy functions. Then the asymptotic expansion of the Airy functions will leads to the coefficient $C_2$ in the asymptotic approximation $(1)$.
http://functions.wolfram.com/Bessel-TypeFunctions/AiryAi/06/
http://functions.wolfram.com/Bessel-TypeFunctions/AiryBi/06/