What's with the restriction on $p$ ?

Let $a + 8p = x^2, a = y^2$ . We get $x+y = 4+2^b$ and $x-y = \frac{x^2-y^2}{x+y} = \frac{8p}{4+2^b}.$

From this we get that $x$ and $y$ are, respectively, $2+2^{b-1} \pm \frac{p}{1+2^{b-2}}.$

Since $x$ and $y$ are rational numbers whose squares are integers, they must be integers themselves.

If $b =1$ then they are $3 \pm \frac{2p}3$ , so $p = 3$ leads to a solution $(3,1,1)$ .

If $b = 2$ then they are $4 \pm \frac{p}2$ , so $p = 2$ leads to a solution $(2,9,2)$ .

Otherwise, the denominator is an odd integer $> 1$ , which must equal $p$ for $x$ and $y$ to be integral. But the only odd primes that are one more than powers of $2$ are the Fermat primes, i.e. $b-2$ is a power of $2$ . (This is a standard result.)

There are five odd Fermat primes less than $100000$ : $3,5,17,257,65537$ . So there are a total of $\fbox{7}$ solutions in all.

For now I'll just list out the solution triples:

$(2,9,2), (3,1,1), (3,25,3), (5,81,4), (17,1089,6), (257,263169,10), (65537,17180131329,18)$ .

The sequence of primes $3, 5, 17, 257, 65537$ should look familiar.