upto(n) = { my(i, res = List()); maxp = prime(logint(n, 2) + 1); t = 0; primepimap = Map(); forprime(p = 2, maxp, t++; mapput(primepimap, p, t); ); inc = 10^6; printat = inc; for(i = 3, n, if(i > printat, print1(i\inc", "); printat += inc; ); if(isforUpto(i), listput(res, i); print(); print(); print(res); print(); print(); ) ); res } isforUpto(n) = { my(f = factor(n, maxp)); if(f[#f~, 1] > maxp, return(0); ); my(v = vector(mapget(primepimap, f[#f~, 1])), i); if(#v == 1, return(0)); for(i = 1, #f~, v[#v + 1-mapget(primepimap, f[i, 1])] = f[i, 2] ); if(#v == 2, return(denominator((n-v[2])/v[1]) == 1); , if(#v == 3, my(D = v[2]^2 + 4*n*v[1] - 4*v[1]*v[3]); if(issquare(D), B = (-v[2] + sqrtint(D)) / (2*v[1]); return(denominator(B) == 1) , return(0) ) , u = sqrtnint(n, #v - 1); ) ); for(i = 2, u, if(fromdigits(v, i) == n, return(1) ) ); 0 }