Solution (source code)

= Solution
{c}

Choose $\delta=\sqrt\epsilon$ in the <WKB approximation> and substitute $y=\exp(S_0/\delta+S_1+\cdots)$. The leading equation is
$$
(S_0')^2=x,
$$
so take $S_0=\pm2x^{3/2}/3$. At the next order,
$$
2S_0'S_1'+S_0''+(x-1)S_0'=0,
$$
giving
$$
S_1=-\tfrac14\log x-\tfrac14x^2+\tfrac12x
$$
up to an absorbed constant. Label the recessive and dominant branches respectively by
$$
y_1=x^{-1/4}e^{-x^2/4+x/2}e^{-2x^{3/2}/(3\sqrt\epsilon)},\qquad
y_2=x^{-1/4}e^{-x^2/4+x/2}e^{2x^{3/2}/(3\sqrt\epsilon)}.
$$
The formal outer solution is $c_1y_1+c_2y_2$. The point zero is a simple <wave turning point>: $S_0'$ vanishes and the WKB amplitude diverges.

The boundary condition at infinity selects \b[$c_2=0$], with $c_1$ to be fixed by matching. This statement concerns the actual recessive branch. One cannot literally send $x\to\infty$ in the truncated small-$\epsilon$ formulas: their ordering fails when $\epsilon x$ is no longer small, and their common Gaussian factor would then spuriously suggest that both branches decay.

To justify the selection, remove the first derivative exactly with $y=e^{-x^2/4+x/2}v$. Then
$$
v''=\left[\frac x\epsilon+\frac{(x-1)^2}4+\frac12\right]v.
$$
The coefficient is positive for $x\geq0$, so it has distinct dominant and recessive branches. For very large $x$, the dominant branch restores $y$ growth of exponential rate $1/\epsilon$, while the recessive branch decays. Matching to this recessive branch removes $c_2$.