Special arithmetic progressions involving perfect squares

Starting from the equations in my previous answer, we get, by multiplying them in pairs, $$(x-y)x(x+y)(x+2y) + (x-y)x + (x+y)(x+2y) + 1 = (z_1 z_6)^2\,,$$ $$(x-y)x(x+y)(x+2y) + (x-y)(x+y) + x(x+2y) + 1 = (z_2 z_5)^2\,,$$ $$(x-y)x(x+y)(x+2y) + (x-y)(x+2y) + x(x+y) + 1 = (z_3 z_4)^2\,.$$ Write $u = z_1 z_6$, $v = z_2 z_5$, $w = z_3 z_4$ and take differences to obtain $$3 y^2 = u^2 - v^2 \qquad\text{and}\qquad y^2 = v^2 - w^2\,.$$ The variety $C$ in ${\mathbb P}^3$ described by these two equations is a smooth curve of genus 1 whose Jacobian elliptic curve is 24a1 in the Cremona database; this elliptic curve has rank zero and a torsion group of order 8. This implies that $C$ has exactly 8 rational points; up to signs they are given by $(u:v:w:y) = (1:1:1:0)$ and $(2:1:0:1)$. So $y = 0$ or $w = 0$. In the first case, we do not have an honest AP ($y$ is the difference). In the second case, we get the contradiction $abcd + ad + bc + 1 = 0$ ($a,b,c,d$ are supposed to be positive). So unless I have made a mistake somewhere, this proves that there are no such APs of length 4.

Addition: We can apply this to rational points on the surface. The case $y = 0$ gives a bunch of conics of the form $$x^2 + 1 = z_1^2, \quad z_2 = \pm z_1, \quad \dots, \quad z_6 = \pm z_1\,;$$ the case $w = 0$ leads to $ad = -1$ or $bc = -1$. The second of these gives $ad + 1 < 0$, and the first gives $ac + 1 = (a^2 + 1)/3$, which cannot be a square. This shows that all the rational points are on the conics mentioned above; in particular, (weak) Bombieri-Lang holds for this surface.

Already for three-term progressions it's somewhat surprising that there are infinitely many solutions, because the usual probabilistic guess for the expected number of solutions leads to a convergent sum: a random number of size $M$ is a square with probability about $M^{-1/2}$, so we're summing something like $1/(abc)$ over all three-term progressions $(a,b,c)$, etc. To be sure such a guess cannot account for non-random patterns arising from polynomial identities, but it does suggest that past a certain point such identities will be the only source of solutions.

Now a mindless exhaustive search over progressions $(x,x+y,x+2y)$ with $0 < x,y < 10^4$ finds only the first six examples $$ (1,7),\phantom+ (4,26),\phantom+ (15,97),\phantom+ (56,362),\phantom+ (209,1351),\phantom+ (780,5042) $$ of an infinite family associated with the solutions $(2,1)$, $(7,4)$, $(26,15)$, $(97,56)$, $(362,209)$, $(1351,780)$, etc. of the Pell equation $x^2-3y^2=1$. If it can be proved that these are the only solutions then it will immediately follow that there are no four-term arithmetic progressions with the same property. But that seems like a very hard problem.

Here's the gp code; with a bound of $10^4$ it takes only a few minutes. One can surely do better with a more intelligent search procedure (e.g. start by finding all solutions of $ab+1=r^2$ by factoring $r^2-1$).

H = 10^4
progsq(x,y,n) = sum(i=0,n-2,sum(j=i+1,n-1,issquare((x+i*y)*(x+j*y)+1)))
for(x=1,H,for(y=1,H,if(progsq(x,y,3)==3,print([x,y]))))

According to Magma, the projective closure of the variety associated to the problem (given by the equations $$x(x-y) + 1 = z_1^2, \quad (x+y)(x-y) + 1 = z_2^2, \quad (x+2y)(x-y) + 1 = z_3^2,$$ $$(x+y)x + 1 = z_4^2, \quad (x+2y)x + 1 = z_5^2, \quad (x+2y)(x+y) + 1 = z_6^2 \quad)$$ is an irreducible surface in ${\mathbb P}^8$ with 34 isolated singularities. Since it is a complete intersection of six quadrics, it should be of general type (and it has trivial rational points with $x = y = 0$ and slightly less trivial ones with $y = 0$, so reduction methods will not work), which makes it very likely that this is a hard question.

Added later: You may want to look at Question 73346 for an explanation by Noam Elkies of the reasoning behind this.