Finding the Two Squares

Date: 06/11/2003 at 14:42:19
From: Javier (Xavier in english)
Subject: Number Theory

One of Fermat's theorems says that every prime number that yields a 
remainder of 1 when divided by 4 can be expressed as the sum of two 
integer squares (e.g.: 97 = 4^2 + 9^2). This theorem was proven by 
Fermat.

How many methods are known for determining the two squares?

Recently I found a VERY easy method to calculate the two squares. It 
satisfies Fermat's equation for an infinite number of primes, as 
required, but not for all of them. Nevertheless, the concept could be 
interesting if it is still not known. I have tested my method with 
success for large primes. However, prior to publishing it, I would 
like to know if my observations are already known.  

I have prepared a example for a prime with more than 2000 digits:

5404039655023548476563582944026489940017046635919978393526155665887063
8005065726520071926444477856150183440664917310640681970481130525077555
1813089330186960026368746283439276196917097591174393565249415088734555
6578725238218090970466875266948231898935099582536395470695188077936119
2433771645028069183146642042538478132578942801308918185050611489995252
0908905971326635468074600142065310309528854077103841197427569314707284
4775685164355274259963758486561396521956859594280260399915587050758325
5400351927647415565632901827410376314922873985001988412462522813349024
2389210623459328019923365204682592825785893285607060446835655227048908
6574965426189138324445729498805245550926577038174305189357302351075179
0154520680233389589560258806467424680815075324158524134786072899580175
0068380594424435746385241825682145538749858981102306784640620872482092
5889676807213113544915873611251321932264546192328591369284190686405397
8796578316880512700669776845846557466093132336843994002159546498706062
1164704453834674852494088838569^2  + 
8743919838380357872240950333552884755248185799233734769136634226920403
9903706471911826674952059562985647958852613819239813261105188905461207
7875831229024747356357315107097914665106505179999291952312401355303994
6158383356240274174385610624223337333201384806920334421718904568214044
9493671859240701392246232222769417673828536565161132001453437639738150
0783674027108491505944224092684429561924608789102265054897025420125355
6399106401952751551971718373841859725066052370390227133112916072033624
3677709336762316401044741506168607114772696438521747666354137955673699
4592769178855046384251505557875385110054603399787713333296347735641649
9150561586213817623674946918384157119893250460621625898070273088370686
9866600591758453694615375633411248741292683150832411654904263544412460
4990511646043733048102767655959887431915411697178884044946116135815659
9841059268641957438756264022790459041565711695148082183594500446230387
1095440666251657225931565141507375703196852187822745740643739006237185
3290726979456675320377213569186^2  = 
1056597787330886165607283376335772947912004036112564579874375352759594
4392457501150256617516357718584986040635382499191505219101434089745828
2903824132176589968697437206589254353815861733407326847272527536408390
5833609186070912661512269116100758775431451837305406774829615930445817
6073261506856587953158671116932648237320048396566645752371112881790992
3244046629843158210452708553287703154960396488410169195658516335585159
1142327980744624330187879064855304798755771686797051547007947771882372
9526678934271436902243810677802914882915258441073230342319674534696073
8979826513948051491557075401036263453606330905092506968405840288315879
3747384557423242752292859137029927002634958591300831686286270148295738
1416422651851594740923740717278352278790933377212189709600451144671795
2125658282085764447893961587269325976762505165678443632053621488809924
6970935202477840794020001804824934094322942647313245015464146235410381
6909567885032645614056181741431685587944855905936219050748434454698615
2566196524585630781634530055632958969597562844061616974634095986796453
8987362157854244708043842221402221963277716301052548247805899502717820
0520330722553130716575868842408216811009769409987366723221135473664781
2607886111951361904437069608336082860180992488471757551982113395649211
4471403028698547136465705828789219477311237793325673906241774008541206
3517137807422012852658695334528240005246734404228893278003283623355475
1951150166619761910806328590241879930675318727225635622143110385488010
5596403020095508495581348319107302916061047962850483220408241536724453
5726451364346875102787079992815261713755256259556877623989107601134710
8694769329084944328439613846369487730397277862891300268053145113543905
8230938508307644282628222440206903725318897236157720164030865445790151
1491512269954396171770192346058842441982671364532811884245406494975736
4441477052680865532896708837038057967791144876558987357331177237668996
8432217235625731351030824652263016860210144187956704470629348157475577
780655133842312690878825193342235991049335995169792504550670357 =prime

Date: 06/12/2003 at 08:48:46
From: Doctor Jacques
Subject: Re: Number Theory

Hi Javier,

The method I know for solving that problem is based on a decomposition 
in Z[i], the Gaussian integers.

We want to solve:

  x^2 + y^2 = p

where p is a prime congruent to 1 (mod 4). We do this by finding the 
prime factors (x + yi) and (x - yi) of p in Z[i].

We start by finding a quadratic non-residue a (mod p). This is done by 
trial and error, using the Jacobi symbol (a/p) to check for the 
residue character. As half of the numbers are non-residues, the 
procedure should quickly find a residue.

Now, we compute:

  b = a^((p-1)/4) (mod p)

As a is a non-residue, we have:

  b^2 = a^((p-1)/2) = -1 (mod p)

i.e. p divides b^2 + 1 = (b + i)(b - i).

As p divides neither b+i nor b-i, each prime factor of p divides one 
of these two numbers. All we have to do is to compute the GCD of p and 
b+i in Z[i] (which is Euclidean), and this will give us the required 
factorization.

For example, with p = 13:

we take a = 2, which is a non-residue since 13 = 5 (mod 8).

We compute:

  b = 2^3 = 8 (mod 13)

(we can check that 8^2 = -1 (mod 13))

We compute the GCD of 13 and 8 + i:

  13 = (8+i)*2 + (-3 - 2i)
  8+i = (-3 - 2i)*(-2 + i)

which gives the GCD as (-3 - 2i), and we have the solution:

  13 = 3^2 + 2^2

Note that this procedure works for all primes (congruent to 1 mod 4, 
of course).

Is this like what you have in mind? Write back if you'd like to talk 
about this some more, or if you have any other questions.

- Doctor Jacques, The Math Forum
  http://mathforum.org/dr.math/