Is there an algorithm for writing a number as a sum of three squares?

Question

By Gauss's Theorem, every positive integer $n$ is a sum of three triangular numbers; these are numbers of the form $\frac{m(m+1)}2$. Clearly $$ n = \frac{m_1^2+m_1}2 + \frac{m_2^2+m_2}2 + \frac{m_3^2+m_3}2, $$ so multiplying through by $4$ and completing the squares gives $$ 8n+3 = (2m_1+1)^2 + (2m_2+1)^2 + (2m_3+1)^2. $$ Thus writing $n$ as a sum of three triangular numbers is equivalent to writing $8n+3$ as a sum of three (necessarily odd) squares.

My question is;

Is there an algorithm for writing a positive integer as a sum of three squares?

Answer 1

A representation of $n$ as a sum of three triangular numbers is equivalent to representing $8n+3$ as a sum of three odd squares. The question of computing representations as a sum of three squares has been much discussed here, see Efficient computation of integer representation as a sum of three squares

Answer 2

This problem is discussed in my paper with Rabin, Randomized algorithms in number theory, Commun. Pure Appl. Math. 39, 1985, S239 - S256. We give an algorithm that, assuming a couple of reasonable conjectures, will produce a representation as a sum of three squares in random polynomial time.

Answer 3

One point that I do not see in the answer to which Igor links is size. Your target number is some $k \equiv 3 \pmod 8.$ So we take some odd $z$ and find out whether $k - z^2$ is the sum of two squares by factoring. My advice is to take $z$ as large as possible to begin, the decrease $z$ by 2 at each failure. There are two reasons for this.

First, the numbers $j \equiv 2 \pmod 8$ that actually are the sum of two squares are more frequent the smaller the approximate size of $j.$ Combining all congruence classes $\pmod 8,$ the number of integers up to some real positive $x$ is about $$ \frac{0.7642 \; x}{\sqrt{\log x}}, $$ so they get less frequent near $x$ as $x$ gets bigger.

Second, deciding whether $k-z^2$ is the sum of two squares is just factoring, and factoring is quicker for smaller numbers: powers of $2$ are irrelevant, any positive integer $j$ is the sum of two squares if and only if, when factoring $j,$ the exponent of any prime divisor $q \equiv 3 \pmod 4$ is even. Indeed, if there are any such, what you actually do is divide out all the appropriate $q^{2a}$ to arrive at a smaller number $j_0,$ write that as $x_0^2 + y_0^2 = j_0$ by solving that for each remaining prime power $p^w$ with $p \equiv 1 \pmod 4,$ which involves finding a square root of $-1 \pmod p$ and then screwing around. Combining pieces comes from $$ (a^2 + b^2)(c^2 + d^2) = (ad-bc)^2 + (ac + bd)^2. $$ Oh, when yopu are done with $x_0^2 + y_0^2 = j_0,$ you put back each $q \equiv 3 \pmod 4$ with $(q^a x_0)^2 + (q^a y_0)^2 = q^{2a} j_0.$

Well, there is more to it, as you can see. But start with large $z.$ Size Matters.