Four squares such that the difference of any two is a square?

Question

I. This post asks to find $4$ integers $a,b,c,d$ such that the difference between any two is a square. As mentioned by my answer, it is equivalent to finding $3$ squares such that the difference of any two is also a square.

With the positive answer to that question, the OP of that post muses if we can also find $FIVE$ integers $a,b,c,d,e$ such that the difference of any two is a square. Equivalently, we are to solve System 1,

$$p^2+w^2 = x^2\\ q^2+w^2 = y^2\\ r^2+w^2 = z^2\\ s^2+x^2 = y^2\\ t^2+x^2 = z^2\\ u^2+y^2 = z^2$$

If this has a positive solution, then it involves three special hypotenuse $\color{blue}{x,y,z}$ expressible as Pythagorean triples in $1$, $2$, and $3$ ways,

$$p^2+w^2 = \color{blue}{x^2}$$ $$q^2+w^2 = s^2+x^2 = \color{blue}{y^2}$$ $$r^2+w^2 = t^2+x^2 = u^2+y^2 = \color{blue}{z^2}$$

which seems doable.

---

II. Alternatively, by solving the system for $p,q,r,s,t,u$, then we can reduce the number of variables to finding just four squares $w^2< x^2<y^2<z^2$ such that the difference between any two is a square, or System 2,

$$-w^2+x^2 =\square_1\\-w^2+y^2=\square_2\\-w^2+z^2=\square_3\\-x^2+y^2=\square_4\\-x^2+z^2=\square_5\\-y^2+z^2=\square_6$$

For the special case $w = 0$, the smallest of infinitely many solutions is,

$$w,x,y,z = 0, 153, 185, 697$$

Q: More generally, can we find four squares $w^2,x^2,y^2,z^2$ of System 2 such that $w\neq0$?

Update: Using zwim's data here, we find that for $w,x,y,z = 448, 952, 1073, 1105$, then,

$$-w^2+x^2 =840^2\\-w^2+y^2=975^2\\-w^2+z^2\neq\square_3\,\,\\-x^2+y^2=495^2\\-x^2+z^2=561^2\\-y^2+z^2=264^2$$

Almost, but not quite. But I believe a higher search range will yield something.

Answer 1

I could manage to find many almost solutions (with the aim to derive full solutions) by a systematic search:

My first trials, code and ideas are posted here, where I received very helpful input that directed me to a promising approach, proposed by Peter: Search for 6-tuples $(s,t,u,t+u,t+u−s,t−s)$ consisting of perfect squares. I implemented a Python script which systematically searched for those 6-tuples and it found a lot of instances. Unfortunatelly at a certain point it became slow and did not scale anymore.

With the help of Arty, who developed an enhaced version, I could manage to find a first allmost matches. This performance breakthrough origins from this question, I asked on SO.

On a machine with 128GB RAM, an Intel Xeon CPU E5-2630 v4, 2.20GHz processor (and two graphic cards of type Tesla V100 with 16GB RAM) I generated a large textfile of these 6-tuples $(s,t,u,t+u,t+u−s,t−s)$.

To make the final step, namely to identify these four squares, I implemented the Mathematica Script pythagorean.nb, which outputs an "almost solution": [w=40579, x=-58565, y=-65221].

From the data set I selected the corresponding row:

42228, 51060, 185472, 1783203984, 2607123600, 34399862784, 37006986384, 35223782400, 823919616

Hence $s=42228^2,t=51060^2,u=185472^2$. To calculate the last variable $z$, I need to use the equation $u+y^2=z^2$ leading to $z=196605,294$.

Now my idea is to multiply all integers $w,x,y,z$ with a suitable integer, for example with $10^{10}$ to bring us closer towards a solution:

• $z^2-y^2=1854720000103154^2$
• $z^2-x^2=1876800000101940.2^2$
• $z^2-w^2=1923720000099454^2$
• $y^2-x^2=287040000000000^2$
• $y^2-w^2=510600000000000^2$
• $x^2-w^2=422280000000000^2$

An excerpt of the list containing such "almost" solutions, I could generate so far is (raw output of Mathematica):

• $w=40579,x=-58565,y=-65221,z=196605.2940$
• $w=53900,x=-71540,y=-91756,z=157415.4376$
• $w=534240,y=-540600,z=-554115,x=537302.7941$
• $w=67375,x=-89425,y=-114695,z=196769.2970$
• $w=81158,x=-117130,y=-130442,z=393210.5880$

The Mathematica Code is shown below:

arr = Import[
   "C:/Users/esultano/git/pythagorean/pythagorean_stu_Arty_.txt", 
   "CSV", "HeaderLines" -> 0];
f[i_] := Part[arr[[i]], 1 ;; 3];
len = Length[arr];
For[i = 1, i < len, i++, {
  triplet = f[i];
  s = triplet[[1]];
  t = triplet[[2]];
  u = triplet[[3]];
  ins = FindInstance[
    x*x - w*w == s*s && y*y - w*w == t*t && w != 0, {w, x, y}, 
    Integers];
  If[Length[ins] > 0, 
   Print[Append[First[ins], 
     z -> N[Sqrt[u^2 + First[ins][[3, 2]]^2], 10]]], Continue];
  ins = FindInstance[
    y*y - w*w == t*t && z*z - y*y == u*u && w != 0, {w, y, z}, 
    Integers];
  If[Length[ins] > 0, 
   Print[Append[First[ins], 
     x -> N[Sqrt[s^2 + First[ins][[1, 2]]^2], 10]]], Continue];
  }
 ]