This is a follow up of Three Cases of Primes From Squares which is itself a follow up of prior blogs articles. I started this journey by recalling the theorem that all primes congruent to 1 mod 4 can be written as a sum of two squares. I then recreated the algorithm for efficiently computing the two squares using the Euclidean algorithm on the Gaussian integers. This then led me to examine finding primes and their decompositions into integer squares for other d in the expression p = d·a2 + b2 where p is a prime. The cases d = 1 and d = 4 are equivalent, so I started with the cases of d = 2 and d = 3 where the quadratic imaginary extension produced by adding the square root of -d was still Euclidean. At that point I thought that I was done. I did not, at first, see a practical way to find the a and b for arbitrary values of d. But, as what is typical for me, I did not know how to stop picking at the problem. One day, when I figured out that there was actually an efficient algorithm for finding square roots mod p, I saw a way forward. But the main breakthrough was when I discovered that the Euclidean algorithm could be extended to non Euclidean fields in a useful way. Investigating this approach, I ran into another roadblock problem that I describe in Approximating Integers with Multiples of Rational Numbers. After eventually finding an efficient solution to that problem, I had all the pieces for breaking down primes with 100s of digits into d times a square plus a square for d up to about 10,000 in less than a minute on my laptop.
However, once I had code that could do arbitrary d, I wondered about finding primes that could do many d simultaneously. For example, what is the first prime that can be written as two times a square plus a square, three times a square plus a square, all the way up to 10 times a square plus square? The answer, which you can get if you execute the code in the link to the Groovy Console below, is 1009. Note that this prime is not as large as you might think it should be. In researching this general problem, I found some other interesting behaviors, which I will describe in blog articles following this one.
The Groovy code for the arbitrary and many simultaneous d can be found at the Gist Arbitrarily split primes into a multiple of a square plus a square. A preloaded and executable version can be found at Groovy Web Console. When executed it may take time for the code to be compiled. Results for simpler cases usually return in under 10 seconds. And as described in earlier blog articles, the Result tab is where you can find the primes and their decompositions. If you are truly interested in reading and understanding the code, you should probably start with prior iterations in the earlier blog articles.
When doing larger d, one issue immediately comes up. The solutions you produce are for multiples of the prime, not just for the prime itself. As an example, if you examine the case of d = 5, the prime 7 cannot be written a 5 times a square plus a square. Instead it is 14 = 2·7 that can be written as one squared times 5 plus three squared. I deal with that by filtering out primes that have such solutions and only keeping the ones where p = d·a2 + b2 with no multiple of p. But the code has a flag where this filtering can be turned off.