OFFSET
1,1
PROG
(PARI) v=vector(M=10^4); for(a=1, M, v[a] = []; fordiv(2*a, x, if(x<(y=2*a/x) && issquare(x^2+y^2, &z) && 1==gcd([x, y, z]), v[a] = concat(z, v[a])))); concat(v) /* Michael Somos, Dec 21 2016 */
CROSSREFS
KEYWORD
nonn
AUTHOR
Lekraj Beedassy, Aug 17 2006
EXTENSIONS
a(20)=145 corrected by Philippe Guglielmetti, Dec 14 2016
a(43)=401 inserted by Michael Somos, Dec 21 2016
STATUS
approved