(PARI) lista(nn) = forstep(n=183, nn, [3, 5, 7, 5, 3, 1, 4, 7, 5, 3, 5, 7, 5, 3, 5, 7, 5, 3, 5, 7, 4, 1], my(s = eval(concat(Str(n), Str(n+1)))); if(issquare(s), print1(s, ", "))) \\
Iain Fox, Dec 27 2017
(PARI) eea(x, y) = my(a=max(x, y), b=min(x, y), s=0, so=1, st, r=b, ro=a, rt, q, t); while(r, q=ro\r; rt=r; r=ro-q*r; ro=rt; st=s; s=so-q*s; so=st); t=(ro-so*a)\b; if(x>y, [so, t], [t, so]) \\ Extended Euclidean Algorithm
lista(nn) = my(res=Set(), b, f2, c, s); for(d=3, nn, b=10^d+1; fordiv(b, f, if(f!=1 && f!=b, f2=b/f; if(gcd(f, f2)==1, c=eea(f, f2); if(c[1]<0, s=f*(f2+2*c[1])*f2*(f-2*c[2])+1, s=f*(2*c[1])*f2*(-2*c[2])+1); if(#digits(s)==d*2, res=setunion(res, Set(s))))))); Vec(res) \\ (Will find all values of length nn*2 or shorter)
Iain Fox, Oct 16 2021