Heat Equation with lower order terms and separation of variables. Let $u_t=a^2u_{xx}+hu_x$ with the following conditions: $u(0,t)=0=u(L,t)$ and $u(0,t)=f(x)$ such that $x \in (0,L)$ and $t>0$.
I'm using the method of separation of variables.  Let $u(x,t)=X(x)T(t)$.  Then:
$$XT'=a^2X''T+hX'T$$
Shifting stuff around we get:
$$a^2X''+hX'-\lambda X=0\text{ and }T'-T\lambda=0.$$
What's throwing me off here is trying to apply separation of variables with the additional $u_x$ term and then moving on to finding eigenvalues.
 A: If you want to solve directly, and have orthogonality of eigenfunctions in $x$, you'll need to write the equation for $X$ in Sturm-Liouville form. An integrating factor will help:
$$
                  e^{hx/a^2}a^2X''+e^{hx/a^2}hX'=\lambda e^{hx/a^2}X \\
                  \frac{d}{dx}\left(e^{hx/a^2}a^2\frac{dX}{dx}\right)=\lambda e^{hx/a^2}X
$$
The reason for doing this is to discover an orthogonality relation between eigenfunctions in $X$, which you'll need in order to expand in a series. From this you see that, if $X_1,X_2$ are solutions of the equation that vanish at both endpoints of the interval $[0,L]$,
$$
         \int_{0}^{L}X_1(x)X_2(x)e^{hx/a^2}dx=0.
$$
This is because of the weight function multiplying $\lambda$. Then you can solve the original equation directly, and use this orthogonality relation to isolate the Fourier coefficients in the expansion with respect to the eigenfunctions of the equation in $X$. The equation
$$
              X''+\frac{h}{a^2}X'-\frac{\lambda}{a^2} X = 0
$$
is Euler's equation, which has known solutions. So you solve for a solution $X$ such that $X(0)=0$, $X'(0)=1$, and then you determine the set of $\lambda$ for which that solution satisfies $X(L)=0$, and you'll have the eigenvalues, along with the eigenfunctions. Then, once you have the eigenfunctions
$$
           X_1, X_2, X_3, \cdots
$$
with corresponding eigenvalues $\lambda_k$, the general solution is
$$
              u(x,t) = \sum_{k=1}^{\infty}C_kX_k(x)e^{\lambda_k t},
$$
and you'll be able to solve $u(x,0)=f(x)$ using orthogonality with respect to the weight $e^{hx/a^2}$:
$$
             f(x) = \sum_{k=1}^{\infty}C_k X_k(x) \\
        \int_{0}^{L}f(x)X_k(x)e^{hx/a^2}dx = C_k\int_{0}^{L}X_k(x)^2e^{hx/a^2}dx \\
        C_k = \frac{\int_{0}^{L}f(x)X_k(x)e^{hx/a^2}dx}{\int_{0}^{L}X_k(x)^2e^{hx/a^2}dx}
$$
Note: Because of the way you chose $\lambda$ instead of $-\lambda$, the eigenvalues will be mostly negative. You know this from the $T(t)=e^{\lambda t}$ solutions.
A: How to deal  if there is  $u_{x}$ 
We have : $$a^2X''+hX'-\lambda X=0 $$ let $ X=e^{mx} $ which implies : 
$$a^2m^2+hm-\lambda=0 \Rightarrow m=\frac{-h\pm\sqrt{h^2+4a^2 \lambda}}{2a^2}$$
You can see that the only case that gives non trivial solution is when $h^2+4a^2\lambda<0$
Now let $-\gamma^2=h^2+4a^2\lambda$ which implies :
$$X=e^{\frac{-hx}{2a^2}}[c_{1}cos(\frac{\gamma x}{2a^2})+c_{2}sin(\frac{\gamma x}{2a^2})] $$
since $X(0)=X(L)=0$ we have $c_{1}=0$ and 
$$ sin(\frac{\gamma L}{2a^2})=0 \Rightarrow \gamma =\frac{2a^2n\pi}{L} 
\Rightarrow -\gamma^2=-\frac{4a^2n^2\pi^2}{L^2}\Rightarrow \lambda_{n}=-(\frac{h^2}{4a^2}+\frac{n^2\pi^2}{L^2})$$
Now we have :
$$X_{n}(x)T_n(t)=e^{\frac{-hx}{2a^2}}sin(\frac{n\pi x}{L})e^{-(\frac{h^2}{4a^2}+\frac{n^2\pi^2}{L^2})t} $$
Finally we have :
$$ u(x,t)=\sum_{n=1}^{\infty}c_n   e^{\frac{-hx}{2a^2}}sin(\frac{n\pi x}{L})e^{-(\frac{h^2}{4a^2}+\frac{n^2\pi^2}{L^2})t}$$
and from $u(x,0) =f(x)$ you can get $c_{n}$.
