I want to compute a solution of the congruence $$ a^2+b^2\equiv -1\pmod{p}. $$ I know that a solution always exists (this is a simple argument using the pidgeonhole principle). I also know that if $p\equiv 1\pmod{4}$, $-1$ is a quadratic residue mod $p$ and the problem can be solved by calculating a square root $a$ of $-1$ and letting $b=0$.
However, for general $p$, I know of no better method than simply trial and error (which is actually not so bad either since of course the chance to "stumble upon" quadratic residues is about a half). Does anyone know or can anyone think of a more efficient method?
I asked this question already about an hour ago and it was closed as a duplicate, but all the answers in the linked questions only proved the solvability of the congruence or determined the number of solution pairs, which is not what I was looking for. So let me be clear about this: What I am interested in is an algorithm to compute a solution pair $(a,b)$ which is better than just randomly trying values for $a$ and checking if $-1-a^2$ is a square mod $p$.
(I don't need to find all solutions - any particular one will do.)