A more general form of the problem is: $$\Delta u = g$$ with $g > 0$, $u \geq 0$ and $u(x',0) = 0$ for all $x' \in \mathbb{R}^{n-1}$. Here are some things we can say:
If $0 < g < K$ then by $W^{2,p}$ estimates and Sobolev embedding we have an interior $C^{1,\alpha}$ estimate (for any $\alpha$), giving in particular the (scaling-invariant) estimate $$|\nabla u|_{B_{r/2}(x)} \leq \frac{C}{r}\|u\|_{L^{\infty}(B_r(x))}$$ for all $B_r(x) \subset \Omega$, where $C$ depends only on $K$.
If $0 \leq u \leq M$ for some $M$ apriori, then using barriers of the form $-C(|x-(x_0',-1)|^{2-n}-1)$ for any $x_0'$ and $C$ large depending on $M$ we easily obtain $$0 \leq u \leq C(M)x_n.$$ Of course, we must have this apriori bound since quadratic growth like $x_n^2$ is allowed otherwise.
Thus, 2) gives a boundary gradient bound, and when coupled with the scaling invariant interior estimate of 1) we get a full gradient bound.
- As a remark, if $g \geq \delta > 0$ then $u$ must have quadratic growth. Indeed, if not, the function $\frac{\delta}{2n}(|x-e_n|^2 - \frac{1}{2})$ would lie above $u$ on the boundary of a large ball and on $\{x_n = 0\}$ and also be an upper barrier, contradicting that $u \geq 0$.