Exactly, we can sum $S(p)$. Sum of roots of $z^2\equiv -x^2-y^2(mod p)$ is $\sum\limits_{x=0}^{p-1} \sum\limits_{y=0}^{p-1}((\frac{-x^2-y^2}{p})+1)$ If $x=0$ $\sum\limits_{y=0}^{p-1}(\frac{-y^2}{p})=(p-1)(\frac{-1}{p})$ Else, $\sum\limits_{y=0}^{p-1}(\frac{-x^2-y^2}{p})=(\frac{-1}{p}) \sum\limits_{y=0}^{p-1}(\frac{x^2+y^2}{p})=(\frac{-1}{p}) \sum\limits_{r=0}^{p-1} (1+(\frac{r}{p}))(\frac{r+x^2}{p})=(\frac{-1}{p})(\sum\limits_{r=0}^{p-1}(\frac{r+x^2}{p})+ \sum\limits_{r=1}^{p-1}(\frac{1+x^2r^{-1}}{p}))=-(\frac{-1}{p})$ Hence, $S(p)=\sum\limits_{x=0}^{p-1} \sum\limits_{y=0}^{p-1}((\frac{-x^2-y^2}{p})+1)=p^2+(p-1)(\frac{-1}{p})-(p-1)(\frac{-1}{p}))=p^2$

If you want to settle for, albeit weak, the original bound; here is a less technical solution. For $z=0$, consider the triple $(x,y,z)=(0,0,0)$. For any $z\in\{1,2,\dots,p-1\}\triangleq [p-1]$, let $S_z=\{-z^2-y^2\pmod{p}:0\leqslant y\leqslant \frac{p-1}{2}\}$. We claim $S_z$ contains at least one quadratic residue: Recall that there are exactly $\frac{p+1}{2}$ quadratic residues modulo $p$, call this set by $Q$; and suppose $S_z\cap Q=\varnothing$. Note that since $S_z\cup Q\subseteq \mathbb{Z}_p$, it follows that $p=|\mathbb{Z}_p|\geqslant |S_z\cup Q|=|S_z|+|Q|=p+1$, a contradiction (where the penultimate equality uses the disjointness assumption). Hence the claim. Now, for any $z\in[p-1]$, as I have established above, there is a $y_0$ such that $-z^2-y_0^2\equiv x_0^2\pmod{p}$ for some $x_0$. Now, if $x_0=y_0=0$, then $z=0$, which is not possible. Hence, either $x_0$ or $y_0$ is non-zero. Suppose this is $x_0$, then $(x_0,y_0,z)$ and $(p-x_0,y_0,z)$ are both pairs enjoying the condition. Thus, for each distinct $z\in[p-1]$, we obtain at least two different triples. This yields a total of at least $2(p-1)+1=2p-1$ triples (note that each triple is distinct: their third coordinate are different, and $p-x_0\neq x_0$). This yields the desired bound.