<s>[INST] Fermat's infinite descent for finding the squares that sum to a prime Fermat's theorem on sum of two squares states that an odd prime $p = x^2 + y^2 \iff p \equiv 1 \pmod 4$
Applying the descent procedure I can get to $a^2 + b^2 = pc$ where $c \in \mathbb{Z} \gt 1$
I want $c = 1$, so how do I proceed from here? How do I apply the procedure iteratively?
Example:
$$
p = 97
$$
$$97 \equiv 1 \pmod 4 \implies \left(\frac{-1}{97}\right) = 1 \implies x^2 \equiv -1 \pmod {97}$$ has a solution
$$x^2 + 1 \equiv 0 \pmod {97}$$
$$x^2 + 1 = 97m$$
We find an $x,m$ that solves the equation.
$$x = 75, m = 58$$
Now, we pick an $a,b$ such that $\frac{-m}{2} \leq a,b \leq \frac{m}{2}$
$$a \equiv x \pmod m = 17$$
$$b \equiv y \pmod m = 1$$
Observations:


*

*$ a^2 + b^2 \equiv x^2 + 1 \equiv 0 (\mod m)$

*$ (a^2 + b^2) = mc$

*$ (x^2 + 1) = mp$
Plugging in $a,b,m$ for 2, we get $c = 5$

By this identity, we know that
$(a^2 + b^2)(c^2 + d^2) = (ac + bd)^2 + (ad - bc)^2$
**$(a^2 + b^2)(x^2 + 1^2) = (ax + b)(a - bx) = m^2pc$
Dividing ** by $m^2$, $pc = (\frac{ax+b}{m})^2 + (\frac{a-bx}{m})^2$
Plugging in $a,b,m,p,c$ we get that $22^2 + (-1)^2 = 97*5$
So we have two squares that add up to 5 times our $p$. How do we turn the 5 into a 1? What is the next step in the descent?
 [/INST] We need to know that if $q$ is prime and $q\equiv 3\pmod 4$ then $q|(x^2+y^2)\to (q|x\land q|y)$. Now for prime $p\equiv 1\pmod 4$ we have $c p=a^2+b^2$ where $p/2>a>0<b<p/2$. This implies $0<c<p$. We ask for the smallest possible positive $c$. We show that if $c>1$ then $c$ is not the minimum. First,if $q$ is prime and $q|c$, where $q\equiv 3\pmod 4$, then $q|c\to q|(a^2+b^2) \to (q|a\land q|b) \to q^2|c\to (c/q^2)p=(a/q)^2+(b/q)^2$; and $0<c/q^2<c$ so $c$ is not the minimum. Second, suppose that every prime $r$ such that $r<p$ and $r\not \equiv 3\pmod 4$ is the sum of two squares. (This is the Inductive Hypothesis.) Then if the prime $r$ divides $c$ and $r\not \equiv 3\pmod 4$ then $r=(s^2+t^2) $ for some positive integers $s,t,$....Recall that $r$ also divides $a^2+b^2$. Let $s\equiv i t\pmod r$ and $a\equiv j b\pmod r$  where $i^2\equiv j^2\equiv -1\pmod r$ and $j\equiv \pm i\pmod r. $ In the case $j\equiv i\pmod r$ we have, with $c'=c/r,$ $$ p c' r^2=p c  r=c'(a^2+b^2)(s^2+t^2)=c'((a s +b t)^2+(a t-b s)^2).$$ Now $(a s +b t)\equiv j i b t+b t\equiv i^2 b t+ b t\equiv 0\pmod r$ and similarly $a t-b s\equiv 0\pmod r$. Therefore for some $u,v$ we have $$0\ne (r u)^2+(r v)^2=(a s+b t)^2+(a t-b s)^2=p c' r^2$$ from which $$0\ne p c'=u^2+v^2$$ Since $0<c'<c$, we see that $c$ is not minimum.....In the case $j\equiv -i\pmod r$ we use the equality $(a^2+b^2)(s^2+t^2)=(a s-b t)^2+(a t+b s)^2.$
</s>