We assume that the stock price process $\{S_t,\,t>0\}$ satisfies, under the real-world probability measure $P$, an SDE of the form
\begin{align*}
dS_t=S_t\big((\mu-q)dt+\sigma dW_t\big),
\end{align*}
where $\{W_t, \, t >0\}$ is a standard Brownian motion.
Here, we need to consider the total return asset $e^{qt}S_t$, that is, the asset with the dividend payments invested in the same underlying stock.
We consider a locally risk-free self-financing portfolio of the form
\begin{align*}
\pi_t =\Delta_t^1 \big(e^{qt}S_t\big) + \Delta_t^2 V_t,
\end{align*}
where $V_t$ is the option price. Then,
\begin{align*}
d\pi_t &= \Delta_t^1 d\big(e^{qt}S_t\big) + \Delta_t^2 dV_t\\
&= \Delta_t^1 e^{qt}\big(q S_t dt + dS_t \big) + \Delta_t^2\left(\frac{\partial V}{\partial t}dt + \frac{\partial V}{\partial S}dS_t + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2 dt\right)\\
&=\left[\mu\Delta_t^1 e^{qt} S_t + \Delta_t^2\left(\frac{\partial V}{\partial t} + (\mu-q) S_t \frac{\partial V}{\partial S} + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2 \right)\right]dt \\
&\qquad\qquad\qquad\qquad\qquad\quad + \left(\sigma\Delta_t^1 e^{qt}S_t + \sigma \Delta_t^2 S_t \frac{\partial V}{\partial S}\right)dW_t.
\end{align*}
Since $\pi_t$ is locally risk-free, we assume that $\pi_t$ earns the risk-free interest rate $r$, that is,
\begin{align*}
d\pi_t = r \pi_t dt,
\end{align*}
Then,
\begin{align*}
&\left[\mu \Delta_t^1 e^{qt} S_t + \Delta_t^2\left(\frac{\partial V}{\partial t} + (\mu-q) S_t \frac{\partial V}{\partial S} + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2 \right)\right]dt \\
&\qquad\qquad\qquad\qquad\qquad + \left(\sigma\Delta_t^1 e^{qt} S_t + \sigma \Delta_t^2 S_t \frac{\partial V}{\partial S}\right)dW_t= r \pi_t dt.
\end{align*}
Consequently,
\begin{align*}
\sigma\Delta_t^1 e^{qt}S_t + \sigma \Delta_t^2 S_t \frac{\partial V}{\partial S}=0, \tag{1}
\end{align*}
and
\begin{align*}
\mu e^{qt} \Delta_t^1 S_t + \Delta_t^2\left(\frac{\partial V}{\partial t} + (\mu-q) S_t \frac{\partial V}{\partial S} + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2 \right) = r(\Delta_t^1 e^{qt}S_t + \Delta_t^2 V_t).
\end{align*}
From $(1)$,
\begin{align*}
\Delta_t^1 = -e^{-qt} \Delta_t^2 \frac{\partial V}{\partial S}.
\end{align*}
Then,
\begin{align*}
-\mu \Delta_t^2 S_t \frac{\partial V}{\partial S}+ \Delta_t^2\left(\frac{\partial V}{\partial t} + (\mu-q) S_t \frac{\partial V}{\partial S} + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2 \right) = r\Big(-\Delta_t^2 S_t\frac{\partial V}{\partial S} + \Delta_t^2 V_t\Big),
\end{align*}
or
\begin{align*}
\Delta_t^2\left(\frac{\partial V}{\partial t} -q S_t \frac{\partial V}{\partial S} + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2\right) &= r\Delta_t^2\Big(-\frac{\partial V}{\partial S} S_t + V_t\Big). \tag{2}
\end{align*}
Canceling the term $\Delta_t^2$ from both sides of $(2)$, we obtain the Black–Scholes equation of the form
\begin{align*}
\frac{\partial V}{\partial t} + (r-q) S_t \frac{\partial V}{\partial S} + \frac{1}{2}\frac{\partial^2 V}{\partial S^2} \sigma^2S_t^2 -rV = 0.
\end{align*}