Integer Solutions of $ x^2+y^2=161701$?

Question

I want to get all integer solutions of $ x^2+y^2=161701$ using the below theorem I have got the following pair $(401,30) , (399,50)$ because $161701$ is of the form $ 4k+1$ , But I want any analytic method to solve that kind of equation ?

Theorem:A positive integer $n$ is properly representable as a sum of two relatively prime squares if and only if the prime factors of $n$ are all of the form $4k+1$, except for the prime $2$ which may occur to at most the first power.

Answer 1

First of all, we factorize $161701$ in integers, it is: $$ 161701 = 101 \cdot 1601 \ . $$ Now we use the following result about the ring of the Gaussian integers, $$ R=\Bbb Z[i]\ . $$ A prime $p>0$ in the integers, seen now as an element in the unique factorization ring of gaussian integers, has the following chances:

• it is ramified, only in the case $p=2$, and here in the decomposition of $2$ as product of two conjugated (prime) elements $2=(1+i)(1-i)$ the two factors are associated. (They differ by a unit.)
• it splits as a product of two conjugated (prime) elements $p=(a+bi)(a-bi)$.
• it is inert, i.e. also a prime in the ring of gaussian integers.

The quadratic reciprocity laws tell us when we can split and when not. (Last two cases.) In our case, the primes are congruent to one modulo $4$, so we can split uniquely: $$ \begin{aligned} 101 &= 10^2+1^2 =(10+i)(10-i)\ ,\\ 1601 &= 40^2+1^2 =(40+i)(40-i)\ . \end{aligned} $$ So we have the prime factor decomposition in $R$: $$ 161701 = (10+i)(10-i)(40+i)(40-i)\ . $$ Each solution of $(x+iy)(x-iy)=161701$ corresponds (up to units) to grouping above the factors in two parts, in each part we have exactly one of $10\pm i$ and exactly one of $40\pm i$. This leads to the solutions obtained by expanding

• $x+iy=(10+i)(40+i)=399+50i$, and
• $x+iy=(10+i)(40-i)=401+30i$.

All integer solutions are obtained by also taking the conjugates, and also using the four units $1,i,-1,-i$. So we expect $16=4\cdot 4$ solutions with $x,y\in\Bbb Z$. (The computer confirms below.)

---

Some computer experiments:

sage: R.<j> = QuadraticField(-1)
sage: factor(161701)
101 * 1601
sage: factor(R(161701))
(j - 40) * (j - 10) * (j + 10) * (j + 40)
sage: (10+j)*(40+j)
50*j + 399
sage: (10+j)*(40-j)
30*j + 401

And all integer solutions:

sage: for x, y in cartesian_product( [ [-N..N], [-N..N] ] ):
....:     if x^2 + y^2 == 161701:
....:         print "x = %4s and y = %4s" % (x, y)
....:         
x = -401 and y =  -30
x = -401 and y =   30
x = -399 and y =  -50
x = -399 and y =   50
x =  -50 and y = -399
x =  -50 and y =  399
x =  -30 and y = -401
x =  -30 and y =  401
x =   30 and y = -401
x =   30 and y =  401
x =   50 and y = -399
x =   50 and y =  399
x =  399 and y =  -50
x =  399 and y =   50
x =  401 and y =  -30
x =  401 and y =   30
sage: