is(n, z=0, m=39, c=[1,16,22], q=[1,1,10]) = {
	if (setsearch(c, n%m) && isprime(n),
		forstep (x=0, +oo, +1,
			my (none=1);
			my (y0 = -q[2]/(2*q[3]));
			for (y=ceil(y0), oo,
				my (v=q[1]*x^2+q[2]*x*y+q[3]*y^2);
				if (v>n,
					break,
					none = 0;
					if (v==n,
						return (z);
					);
				);
			);
			if (y0,
				forstep (y=floor(y0), -oo, -1,
					my (v=q[1]*x^2+q[2]*x*y+q[3]*y^2);
					if (v>n,
						break,
						none = 0;
						if (v==n,
							return (z);
						);
					);
				);
			);
			if (none,
				return (!z);
			);
		);
	);
	return (0);
}

for (n=1, oo, if (is(n), print (k++ " " n); if (k==1 000, break)));

quit