Solving a 2-variable Second Order Linear Partial Differential Equation Part 1: Initial Problem
I am trying to solve the following partial differential equation. 
$$py + q= a\frac{\partial f}{\partial x} +  by\frac{\partial f}{\partial y} +  c\frac{\partial^2f}{\partial y^2}$$
with boundary condition: $f(d, y) = 0$ and $f(0,y) = 0$. 
where $a,b,c,d, p$ and $q$ are all constants. 
I have tried making the assumption that the function $f$ is separable by making the substitution $f(x,y) = g(x) h(y)$, but to no avail due to the constant $d$.
Appreciate any help in this, thank you! 

Part 2 : By applying Merkh's suggestion
I took Laplace Transform with respect to $x$ for the equation given below
and obtained
$$\frac{1}{s}(py + q) = asF + by\frac{\partial F}{\partial y} + c\frac{\partial^2 F}{\partial y^2}$$
Now I am trying to solve this using Variation of Parameters but do not know what form does $y_1$ and $y_2$ in the article take. Is there a way to figure this out?
Found some Variation of Parameters Examples here. Now, I am just going by trial and error to figure out the form that $y_1$ and $y_2$ takes. 
I am keen to find out about alternative methods as well, thanks!
 A: Have you considered taking the Fourier transform?
Treating $x$ as the time variable, you can take the Fourier transform with respect to $y$ to turn the PDE into an ODE with respect to $x.$  Solve the ODE and transform back to get your solution. 

Part 2)
$y_1$ and $y_2$ are the solutions to the homogeneous equations.  However, we should really be calling them $F_1(y;s)$ and $F_2(y;s)$ to avoid confusion.  I use the semicolon notation to indicate that we are treating $s$ as an independent parameter, at least for now.  You have $c F'' + by F' + as F = 0.$  First step is to simplify as much as possible,
\begin{align}
c F'' + by F' + as F = F'' + \frac{by}{c} F' + \frac{as}{b} F,
\end{align}
motivated by the form $f'' + y f' = (f' y)',$ this may be rearranged,
\begin{align}
F'' + \frac{by}{c} F' + \frac{as}{b} F &= \frac{d}{dy} \bigg( F'(y) \exp{\big( \frac{by^2}{2c}\big)} \bigg) + \frac{a}{c}\exp{\big(\frac{by^2}{2c}\big)}F(y) = 0.
\end{align}
So you are left with a Sturm-Liouville type equation, where there exists entire books based on solving such problems.  Once solved, you have found your $F_1$ and $F_2$ to be used in the variation of parameters approach,
\begin{align}
u_1'(y) F_1(y;s) + u_2'(y) F_2(y;s) &= 0\\
u_1' F_1' + u_2' F_2' &= \frac{py + q}{s}.
\end{align}
However, my original suggestion was to Fourier transform with respect to $y,$ not $x.$  That will take care of the $\frac{\partial^2 }{\partial y^2}$ term and give you are first order ODE in $x.$ I haven't done this out, there could be other difficulties, but that is my first instinct when dealing with an infinite spatial domain (treating $y$ as space).
Taking the laplace transform in the $x$ variable where you have a finite domain $x \in[0,d]$ would lead to some pretty bad complications in applying boundary conditions properly.
