<s>[INST] Diophantine equation ${x ^ 2}+{y ^ 2}+{z ^ 2}={m ^ 2}+{n ^ 2} $ has infinite solutions. What is its general solution formula? Diophantine equation ${x ^ 2}+{y ^ 2}+{z ^ 2}={m ^ 2}+{n ^ 2} $ has infinite solutions. For examples,
$\begin{array}{l}
{1^2}+ {3^2} + {4^2}  = {1^2} + {5^2}\\
{1^2} + {2^2} + {6^2} = {4^2} + {5^2}\\
{322^2} + {562^2} + {597^2} = {601^2} + {644^2}\\
{938^2} + {1063^2} + {4722^2} = {536^2} + {4901^2}\\
 \vdots 
\end{array}$
What is its general solution formula? How to get the general solution formula?
 [/INST] When $z\neq 0,$ that means we want all rational solutions to:
$$m^2+n^2-x^2-y^2=1\tag 1$$
We pick one solution, $\mathbf s_0=(m_0,n_0,x_0,y_0)=(-1,0,0,0)$ and any directional vector $\mathbf p=(p_1,p_2,p_3,p_4)$ with the $p_i$ integers with $\gcd(p_1,p_2,p_3,p_4)=1.$
Then $\mathbf s_0+t\mathbf p$ yields a solution for $t=0$ trivially, and we get another value $t$ where you get another solution.
Specifically, you want $$(-1+p_1t)^2+(p_2t)^2-(p_3t)^2-(p_4t)^2=1$$ which, when simplified and factoring out $t,$ gives you:
$$-2p_1+(p_1^2+p_2^2-p_3^2-p_4^2)t=0\\
t=\frac{2p_1}{p_1^2+p_2^2-p_3^2-p_4^2}$$ where now we assume $p_1^2+p_2^2\neq p_3^2+p_4^2.$
Then our integer solution formula is $$\begin{align}m&=p_1^2-p_2^2+p_3^2+p_4^2\\
n&=2p_1p_2\\
x&=2p_1p_3\\
y&=2p_1p_4\\
z&=p_1^2+p_2^2-p_3^2-p_4^2.
\end{align}$$
This won't give primitive solutions, in that we won't always have $\gcd(m,n,x,y,z)=1,$ but we'll get one from every equivalence class of solutions with $z\neq 0.$
This will also give solutions with $z=0,$ but we don't know yet if it will give all (equivalence classes) of solutions for $z=0.$ But we can show if $z=0,$ the above solution is a scalar multiple of $(m,n,x,y,z)=(p_1,p_2,p_3,p_4,0),$ which gives all solutions.
The fact that three of the numbers are even makes it clear that we can't get all the primitive solutions with $z\neq 0.$ For example, we can get values $m$ odd, so we should be able to get $n$ odd.

For example, $(p_i)=(2,1,1,1)$ gives a solution $(m,n,x,y,z)=(5,4,4,4,3).$
This means we need a solution with $(m,n)$ switched, and possibly scaled .
The corresponding rational solution to $(1)$ is $(m,n,x,y)=(\frac43,\frac53,\frac43,\frac43)$ which is in the direction $p=(7,5,4,4)$ from $s_0.$
That gives $m=56, n=70, x=56,y=56, z=42.$ And that is $14$ times our goal $(m,n,x,y,z)=(4,5,4,4,3).$

This technique works for all homogenous quadratic Diophantine equations: once you find one non-zero solution, you can find an infinite class. But the formula won't general all solutions - it will only generate all solutions of the corresponding rational equation, which means it will generate all solutions modulo scalar multiplication, and with one fixed variable chosen to be non-zero.
In general, we would first assume $z\neq 0,$ apply the above technique, then assume $z=0$ and solve the resulting smaller Diophantine equation using the same technique. But often, the resulting formula for $z\neq 0$ would give (up to equivalence) all solutions for $z=0.$
Sometimes, we get an obvious variable to assume non-zero. $x^2+xy+y^2=z^2$ has only the trivial solution when $z=0,$ so we can assume $z\neq 0$ for all others.
</s>