Find a two term asymptotic expansion of the following problem I want to find a two term asymptotic expansion, for small $\epsilon$, of the solution of the following problem:
$$
y'' - \epsilon y' - y = 1\tag{1}, \, y(0) = 0, \,y(1) = 1
$$
My approach: I assume that the solution has an asymptotic expansion that looks like this:
$$
y\sim y_0 + \epsilon^\alpha y_1 + \epsilon^\beta y_2 + \ldots \tag{2}
$$
with $0<\alpha <\beta<\ldots$.
If I substitute $(2)$ into $(1)$ I get the following equation:
$$
(y_0 + \epsilon y_1 + \ldots)'' - \epsilon(y_0 + \epsilon y_1 + \ldots)' - (y_0 + \epsilon y_1 + \ldots) = 1
$$
I will now try to find $y_0, y_1$ by inspecting different order terms:
$\mathcal{O}(1):$ In order to have balance in the equation, $y_0$ must satisfy
$$
y_0'' - y_0 = 1
$$
Solving this equation for $y_0$ gives 
$$y_0 = c_1\exp((\frac{1}{2} + \sqrt{5}/2)t) + c_2\exp((\frac{1}{2}
) - \sqrt{5}/2)t)$$
This solution needs to satisfy $y(0) = 0$. Substituting $t = 0$ and setting $y_0(0) = 0$ leads to the solution:
$$
y_0 = -\dfrac{1}{e^{1/2 - \sqrt{5}/2}}\exp((1/2 + \sqrt{5}/2)t) + \dfrac{1}{e^{1/2 - \sqrt{5}/2}}\exp((1/2 - \sqrt{5}/2)t)
$$
Furthermore, to have balance we need $\alpha = 1$.
$\mathcal{O}(\epsilon):$ In order to have balance in the equation $y_1$ must satisfy
$$
y_1'' - y_1 = -y_0'
$$
Or
$$
y_1'' - y_1 = \dfrac{d}{dt}(-\dfrac{1}{e^{1/2 - \sqrt{5}/2}}\exp((1/2 + \sqrt{5}/2)t) + \dfrac{1}{e^{1/2 - \sqrt{5}/2}}\exp((1/2 - \sqrt{5}/2)t))
$$
While I could plug this into wolfram to find a solution for $y_1$, I'm not sure whether this is the correct way to find a two term asymptotic expansion for the initial problem.
Question: Am I on the right track? It feels like there must be an easier way to do this. Could someone point me in the right direction?
 A: My calculation is different from yours. Let
$$
y=y_0 + \epsilon y_1  +O(\epsilon^2) 
$$
in the equation to get
$$ (y''_0 + \epsilon y''_1 +O(\epsilon^2))-\epsilon(y'_0 + \epsilon y'_1  +O(\epsilon^2))-(y_0 + \epsilon y_1  +O(\epsilon^2))=1. \tag{1} $$
Then the boundary conditions $y(0)=0,y(1)=1$ become
$$ y_0(0)=y_1(0)=0, y_0(1)=1,y_1(0)=0.$$
$O(1)$:
$$ y''_0-y_0=1, y_0(0)=0,y_0(1)=1 $$
which has the solution
$$ y_0= \frac{e^{-x} \left(e^x-e^{2 x}-e^{x+2}+2 e^{2 x+1}-2 e+e^2\right)}{e^2-1}. $$
$O(\epsilon)$:
$$ y''_1-y_1=y_0', y_1(0)=0,y_1(1)=0 $$
which has the solution
\begin{eqnarray*} y_1&=& \frac1{2 \left(e^2-1\right)^2}\bigg[e^{-x} (-e^{2 x+2} (x-4)+2 e^x-4 e^{x+2}+2 e^{x+4}+e^{2 x} (x-2)+2 e^{2 x+3} (x-2)\\
&&-2 e^{2 x+1} x+e^2 x-2 e x-e^4 (x+2)+2 e^3 (x+2))\bigg]. 
\end{eqnarray*}
Thus
\begin{eqnarray*} 
y&=&y_0+\epsilon y_1+O(\epsilon^2)\\
&=&\frac{e^{-x} \left(e^x-e^{2 x}-e^{x+2}+2 e^{2 x+1}-2 e+e^2\right)}{e^2-1}\\
&& +\frac{\epsilon}{2 \left(e^2-1\right)^2}\bigg[e^{-x} (-e^{2 x+2} (x-4)+2 e^x-4 e^{x+2}+2 e^{x+4}+e^{2 x} (x-2)+2 e^{2 x+3} (x-2)\\
&&-2 e^{2 x+1} x+e^2 x-2 e x-e^4 (x+2)+2 e^3 (x+2))\bigg]. 
\end{eqnarray*}
