18
$\begingroup$

Recently I've been studying the problem of integer representation as sum of three squares. Most of the articles that I've found study the function $r_m(n)$ which counts the number of representations of $n$ as the sum of $m$ squares. However, this is not what I am interested in. What I'm looking for is an efficient way (for some given $n$) to find $x$, $y$ and $z$ such that $n = x^2 + y^2 + z^2$. I need to find at least one such representation. Can you recommend me some articles that study this problem?

P.S. I believe that Emil Grosswald's book "Representation of integers as Sums of Squares" contains the answer. However, I could not find this book on my university's web-site.

Improve this question
$\endgroup$
8
  • 2
    $\begingroup$ For a fast probabilistic approach, based on Cramer's conjecture, simply take small values of $z$ in the appropriate arithmetic progression, and check if $n-z^2$ is a prime number which is $1$ mod $4$. Then, to get $p=x^2+y^2$, use Lagrange's algorithm. $\endgroup$
    – Dror Speiser
    Commented Aug 9, 2012 at 6:59
  • $\begingroup$ Sounds reasonable, but I have a couple of questions. First, notice that your proposition is true only for $m \equiv 1 \pmod{4}$. In that case $x$ and $z$ are odd, $y$ is even, and therefore $m - z^2 = x^2 + y^2 \equiv 1 \pmod{4}$. Is it always possible to find such a prime $p$ so that $p = x^2 + y^2$? Do you know any proof for this? Second, for $m \equiv 3 \pmod{8}$ we have $m - z^2 = x^2 + y^2 \equiv 2 \pmod{8}$. So in this case, I guess, we should look for $p$ such that $2p = x^2 + y^2$. Again, why does this $p$ always exists? Also, could you please give a link to an article about algorithm? $\endgroup$
    – Anton
    Commented Aug 9, 2012 at 7:33
  • 1
    $\begingroup$ I believe there are no references in the literature to this problem, nor does the current knowledge on primes allow to prove theorems such as those above. I also have a copy of Grosswald's book, and I can assure you that it does not contain an answer to this problem. I hate to sound so pessimistic, but hey, on the bright side, if you do find anything better than exponential complexity, you'll probably be the first to do so! $\endgroup$
    – Dror Speiser
    Commented Aug 9, 2012 at 8:51
  • $\begingroup$ @Dror factoring is subexponential, so finding a $z$ that does work will give subexponential complexity or am I missing something? This approach is far from efficient though. $\endgroup$
    – joro
    Commented Aug 9, 2012 at 9:06
  • $\begingroup$ You are missing two things: the first is that what you said isn't exactly true, you don't know how long you'll search for a $z$ that actually works, i.e. that $n-z^2$ is a sum of two squares - sure, once it works you have $x,y$ from factorisation, but how many $z$ must you try before one works? second, subexponential factoring are probabilistic algorithms, without proofs, just like the method that has been proposed. $\endgroup$
    – Dror Speiser
    Commented Aug 9, 2012 at 9:15

2 Answers 2

Reset to default
25
$\begingroup$

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 polynomial time.

Improve this answer
$\endgroup$
2
  • $\begingroup$ Looks like knowing one solution allows full? parametrization over Q. $\endgroup$
    – joro
    Commented Oct 22, 2012 at 12:44
  • $\begingroup$ @JeffreyShallit Is there a paywall free copy available? $\endgroup$
    – Turbo
    Commented Jun 21, 2019 at 6:42
11
$\begingroup$

A modification of Dror's comment.

This probabilistic algorithm worked for me.

The main idea is to pick some $z$, compute $m=n - z^2$, factor $m$ with trial division and express it as a sum of two squares if possible. The probability of finding prime $m=4a+1$ or $2p$ is high enough for practical purposes.

The algorithm:

  1. z:=0
  2. z:=z+1
  3. m:=n-z^2
  4. if can't trial factor m goto 2
  5. if m=x^2+y^2 (the factorization is known) then x^2+y^2+z^2=n. Done
  6. goto 2

Added Experimental results: 100 random integers of size 1000 bits were solved at average time 3.5 seconds per solution. 10 random integers of size 10000 bits were solved at average time 5.8 minutes.

Here is a pari/gp program and example:

/*
? n=nextprime(10^220+30000)*nextprime(2^1000+40000000);n%8
%136 = 3
? t=threesquares(n);

? ##
  ***   last result computed in 2,210 ms.
? round(log(n))
%138 = 1200
*/

pl=10^6;\\ bound for trial division, may need change
default(primelimit,pl);
{
twosquares(n)=
local(K,i,v,p,c1,c2);
K=bnfinit(x^2+1);
v=bnfisintnorm(K,n);
for(i=1,#v,p=v[i];c1=polcoeff(p,0);c2=polcoeff(p,1);if(denominator(c1)==1&&denominator(c2)==1,return([c1,c2])) );
return([]);
}

{
threesquares(n)=
local(m,z,i,x1,y1,j,fa,g);
if(n/4^valuation(n,4)==7,return([]););
for(z=1,n,
\\forstep(z=sqrtint(n),1,-1,
m=n-z^2;
if(m%4==3,next);
print1(z," ",);
fa=factor(m,pl);
g=1;
for(i=1,#fa~,if(!ispseudoprime(fa[i,1])/*||!isprime(fa[i,1])*/||(fa[i,2]%2==1&&fa[i,1]%4==3),g=0;break; ));
if(!g,next);
print("\nfound ",z," "," m=",m,factor(m));
j=twosquares(m);
print("j=",j);
x1=abs(j[1]);
y1=abs(j[2]);
return([x1,y1,z]);
);
}
Improve this answer
$\endgroup$
8
  • $\begingroup$ OK, this works. The question is, does this work in a polynomial time? There exist two algorithms (by Serret and by Schoof) that compute $x$, $y$ for a given $p=x^2+y^2$ in a polynomial time. But how many values of z will we check before finding a proper one? I need this for cryptographic purposes, so efficiency is crucial. $\endgroup$
    – Anton
    Commented Aug 9, 2012 at 8:52
  • $\begingroup$ By the way, here's the online tool for solving the two squares problem: numbertheory.org/php/main_pell.html. Also, I found a book that contains an algorithm by Serret that solves the two squares problem in $O(\log^4{p})$ (based on ERH): plouffe.fr/simon/math/…. $\endgroup$
    – Anton
    Commented Aug 9, 2012 at 8:57
  • $\begingroup$ The algorithm is probabilistic and I don't know how many values of z are needed. You can check the program for your data (probably lowering 10^8 to 10^6, check for yourself). As the program showed 1332 bit integer was solved in about a minute. $\endgroup$
    – joro
    Commented Aug 9, 2012 at 9:00
  • $\begingroup$ That sounds like a pretty good result. Ok, I'll do some more research on how fast this works and how many $z$ should be checked and write about it later. By the way, I think that it is more efficient to start not with $z=0$, but with $z= \lfloor n / \sqrt{3} \rfloor$. Then we'll get the smallest value of $n−z^2$ and it would be more probable for this number to be prime. $\endgroup$
    – Anton
    Commented Aug 9, 2012 at 9:14
  • 3
    $\begingroup$ Here's some heuristics for predicting how many values of $z$ are needed: the solutions to $x^2 + y^2 + z^2 = n$ are roughly equidistributed on the sphere, and the number of solutions is proportional to the class number $h(-4n)$, which grows like $\sqrt{n}$. So consider $C \cdot \sqrt{n}$ points equidistributed on the sphere of area $4 \pi n$. What is the statistical expectation of the smallest (positive) $z$-coordinate among these points? $\endgroup$
    – Marty
    Commented Aug 9, 2012 at 22:07

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .