Solution to Differential Equation $\left( 1-2\lambda\frac{\partial}{\partial z}\right)w(x,y,z)-g(x,y,z+h)+2 \lambda h(x,y,z)=0$ I'm trying to solve the following Differential Equation:
$\left( 1-2\lambda\frac{\partial}{\partial z}\right)w(x,y,z)-g(x,y,z+h)+2 \lambda h(x,y,z)=0$
The unknown function is $w(x,y,z)$.
The function $h(x,y,z)$ is in turn related to $p(x,y,z)$ through the following differential equation:
$\left( 1-2\lambda\frac{\partial}{\partial z}\right)h(x,y,z)-p(x,y,z+h)=0$
this problem represents a slightly harder version of a differential equation already solved in this question:
Solution to differential equation $\left( 1-\lambda\frac{\partial}{\partial z}\right)w(x,y,z)-g(x,y,z+h)=0$
The authors report as solution the following expression:
$w(x,y,z)= \frac{1}{2\lambda}\int_{0}^{\infty}e^{-\frac{s}{2\lambda}}g(x,y,z+h+s)ds+\frac{1}{2\lambda}\int_{0}^{\infty}s \,e^{-\frac{s}{2\lambda}}p(x,y,z+h+s)ds$ 
While the first term of the RHS has  been shown in the link to be the solution of the simplified version of the problem, I still don't understand how the inclusion of the extra term $h(x,y,z)$ in the differential equation leads to second term in the RHS of the solution.
 A: As in my answer to the linked question, I'll consider only the $z$-dependence of each function (as such all functions in this derivation are written as capitalized e.g. $W(z)=w(x,y,z)$). Let the differential operator in this problem and the linked question be written as $\mathcal{D}=1-2\lambda D_z$. Then the results of the linked question may be summarized as: 

If $\mathcal{D}W=G$, then  $\displaystyle W=\mathcal{D}^{-1}G=\frac{1}{2\lambda}\int_{0}^\infty ds\,e^{-s/2\lambda}G(s+z)$.


With this in mind, we address the problem at hand. Writing the given equations using $\mathcal{D}$ and inverting, we have \begin{align}
\mathcal{D}H=P&\implies P=\mathcal{D}^{-1}H,\\
\mathcal{D}W=G-2\lambda P &\implies W=\mathcal{D}^{-1}G-2\lambda \mathcal{D}^{-1}P
=\mathcal{D}^{-1}G-2\lambda \mathcal{D}^{-2}H
\end{align}
Then $\displaystyle \mathcal{D}^{-1}G=\dfrac{1}{2\lambda}\int_{0}^\infty ds\,e^{-s/2\lambda}G(s+z)$ and 
\begin{align}
\mathcal{D}^{-2}H
&=\frac{1}{2\lambda}\int_{0}^\infty ds_1\,e^{-s_1/2\lambda}\mathcal{D}P(s_1+z)\\
&=\frac{1}{(2\lambda)^2}\int_{0}^\infty ds_1 \int_{0}^\infty ds_2\,e^{-(s_1+s_2)/2\lambda}P(s_1+s_2+z)\\
&=\frac{1}{(2\lambda)^2}\int_0^1 dt\int_{0}^\infty ds \,s e^{-s/2\lambda}P(s+z)\\
\end{align}
where we have change of variables to $s:=s_1+s_2,$ $t:=s_1/s$. (Note that $s=s_1+s_2\geq s_1\implies t=s_1/s\leq 1$). Since $\int_0^1 dt=1$, we may combine the integrals constituting $W$ to obtain
$$\boxed{\displaystyle W(z)=\dfrac{1}{2\lambda}\int_{0}^\infty ds\,e^{-s/2\lambda}\bigl(G(s+z)-sP(s+z)\bigr)}.$$
In terms of the functions used in the OP, this becomes
$$\boxed{\displaystyle w(x,y,z)=\dfrac{1}{2\lambda}\int_{0}^\infty ds\,e^{-s/2\lambda}\bigl(g(x,y,s+z+h)-s\,P(s+z+h)\bigr).}$$
