The problem is set as follows: $$\frac{\partial W(\tau,z)}{\partial \tau}=-\frac{1}{\varepsilon} \frac{\partial(z_1 W(\tau,z))}{\partial z_1}-\frac{\partial(z_2 W(\tau,z))}{\partial z_2}+h^2\frac{\partial^2 W(\tau,z)}{\partial z_1^2}$$ where $\varepsilon>0,\, h>0$ with boundary conditions on $\Omega=\{z\in \mathbb{R}^2\vert z_1 \leqslant a,\,z_2 \in \mathbb{R} \}$, where $a>0$: $$W(\tau,z)\vert_{\partial \Omega}=0. $$ and initial conditions $$W(0,z)=\delta(z).$$ Is it possible to find a classical solution for that problem? The Dirac delta function for the initial condition can be substituted with some smooth enough function in order to simplify analysis. Can the solution of the stationary problem be obtained?
$\textbf{Edited on November 28}$
There is a mistake in the initial system. The corrected one looks as follows: $$\frac{\partial W(\tau,z)}{\partial \tau}=\frac{1}{\varepsilon} \frac{\partial(z_1 W(\tau,z))}{\partial z_1}+\frac{\partial(z_2 W(\tau,z))}{\partial z_2}+h^2\frac{\partial^2 W(\tau,z)}{\partial z_1^2}.$$
Not a solution, just a few thoughts about a product ansatz, which the solution is most likely not of the form. If $$ W(t,z_1,z_2)=w_1(t)w_2(z_1)w_3(z_2) $$ then one obtains the three equations \begin{align} \left\{ \frac{{\rm d}}{{\rm d}t} - c_1\right\} w_1(t) &= 0 \\ \left\{ \frac{{\rm d}^2}{{\rm d}z_1^2} w_2(t) - \frac{z_1}{h^2 \epsilon} \frac{\rm d}{{\rm d}z_1} + \frac{c_2}{h^2\epsilon} \right\} w_2(z_1) &= 0 \\ \left \{ \frac{\rm d}{{\rm d}z_2} + \frac{\epsilon c_1 + \epsilon + c_2 + 1}{z_2 \epsilon} \right\} w_3(z_2) &= 0 \end{align} with the solutions \begin{align} w_1(t) &= C_1 \, {\rm e}^{c_1 t} \tag{1} \\ w_2(z_1) &= C_{21} \, z_1 \, M\left(\frac{1-c_2}{2} , \frac{3}{2} , \frac{z_1^2}{2h^2\epsilon} \right) + C_{22} \, z_1 \, U\left(\frac{1-c_2}{2} , \frac{3}{2} , \frac{z_1^2}{2h^2\epsilon} \right) \tag{2} \\ w_3(z_2) &= C_3 \, z_2^{-\left(c_1 + 1 + \frac{c_2+1}{\epsilon}\right)} \tag{3} \end{align} where $c_1,c_2,C_1,C_{21},C_{22},C_3$ are arbitrary constants and $M,U$ are the two Kummer solutions.
You say that the solutions must vanish at $z_1=a, z_1=-\infty$ and $z_2=\pm \infty$:
If the first argument in the two Kummer solutions of $w_2(z_1)$ was a negative integer, they were polynomials of some degree and thus would not go to zero in the limit $z_1\rightarrow -\infty$. Thus $M$ and $U$ are not polynomials in which case we must have $C_{21}=0$ since $$M(a,b,z)\sim \frac{\Gamma(b)}{\Gamma(a)} \, {\rm e}^{z} \, z^{a-b}$$ diverges exponentially. For the second term we have $$ z_1 \, U\left(\frac{1-c_2}{2} , \frac{3}{2} , \frac{z_1^2}{2h^2\epsilon} \right) \sim |z_1|^{c_2} \qquad {\rm as} \quad z_1 \rightarrow -\infty $$ and so $c_2<0$. From (1) we must as well require $c_1<0$ and from (3) it follows $c_1+1+\frac{c_2+1}{\epsilon} > 0$ because of the boundary condition at $z_2=\pm \infty$.
Turning back to $U$ however it seems we need to find zeros in order to fulfill the boundary condition at $z_1=a$. But the total number of zeros for $U(a,b,z)$ is given by $T(a,b)=T(a-b+1,2-b)$ if $b\geq 1$ which here becomes $T\left(-\frac{c_2}{2},\frac{1}{2}\right)=0$ since the first argument $-c_2 >0$ (cf. https://dlmf.nist.gov/13.9).
Thus is seems there is no solution $\neq 0$ with this product-ansatz, since the boundary condition at $z_1=a$ can not be fulfilled.