I will consider the case $k = 1$ and $p = 2$ (some arguments may generalize to $k \in \mathbb{N}$).
Let us use the norm $\|u\|^2 = \|u\|^2 + \|\nabla u\|^2$ in $H^1(\Omega)$ (both are $L^2$-norms).
The associated scalar product is denoted by $(\cdot,\cdot)$.
We denote by $K = \{v \in H^1(\Omega) : v \ge 0\}$ the positive cone in $H^1(\Omega)$.
Then, the projection of $u \in H^1(\Omega)$ onto $K$ is given by the minimizer $v$ of $\frac12 \, \| u - v \|^2$ over $K$.

Now, one can show, that there is a Lagrange multiplier $\lambda \in (H^1(\Omega))'$.
This multiplier lies in the polar cone $K^\circ$, that is
$$
\langle \lambda, w \rangle \le 0 \quad\mbox{for all } w \in K
$$
and is orthogonal to $v$, i.e., $\langle \lambda , v \rangle = 0$.
Moreover, we have the equation
$$
(v - y, w) + \langle \lambda, w \rangle = 0 \quad\mbox{for all } w \in H^1(\Omega),
$$
which couples all involved quantities. This equation is the weak formulation of a PDE.

Note that these optimality conditions are nothing more than an equivalent reformulation of the projection inequality
$$
(v - y, w - v) \ge 0 \quad\mbox{for all } w \in K.
$$

The inclusion $\lambda \in K^\circ$ can be interpreted as non-positivity of $\lambda$.
Moreover, it should be possible to show (at least in the case $n = 1$, you are interested in), that $\lambda$ is represented by some negative measure (along the lines that negative distributions are essentially measures).

Then, the complementarity $\langle \lambda , v \rangle = 0$ means, that $\lambda$ is only strictly negative, where $v = 0$ (i.e. where the projection is active).

As you can see,
the calculation of the projection requires the solution of a nonlinear (even non-differentiable) PDE.
This PDE is just the optimality condition (i.e., the projection inequality) of minimizing the distance in the Dirichlet energy.
Hence, there is no easy formula.

Finally, let me comment on the case $u \in H^2(\Omega)$.
That is, we are projecting a more regular element (but still w.r.t. the $H^1$-norm).
Then one can show (see, e.g., the book by Kinderlehrer and Stampacchia - I can provide an exact reference, if needed),
that (under some regularity of $\Omega$)
again $v \in H^2(\Omega)$ and $\lambda \in L^2(\Omega)$.
Due to this regularity, on can interpret the above relations pointwise and obtain
$\lambda(x) \ge 0$ and $v (x) \, \lambda(x) = 0$ for almost every $x \in \Omega$.
Finally, the PDE can be written as
$$
\max\Big( {-\Delta(u-v)} + (u-v), \; v \Big) = 0 \quad\mbox{almost everywhere in }\Omega.
$$
Again, one can see, that there is no easy (explicit) formula for the projection.

I found this 'high regularity case' interesting, because its similar to the case of projecting a $H^1$-function w.r.t. the $L^2$-norm, where the projection is again $H^1$.

projection. Of course, you can take the element in the positive cone with the smallest distance (which is unique for $p \in (1,\infty)$). But I don't think that this is your question. In the case $k = 1$, $p = 2$, this projection is given by the obstacle problem $\endgroup$