Proposition. Let $n$ be a positive odd integer. There is a one-to-one correspondence between factorizations of $n$ in the form $n=ab$ , where $a\geq b>0$, and the representations of $n$ in the form $t^2-s^2$, where $s$ and $t$ are nonnegative integers. The correspondence is given by the equation
$$t=\frac{a+b}{2},\ s=\frac{a-b}{2};\ a=t+s,\ b=t-s$$

Proof. Given $n=ab$, we can write
$$n=ab=\left(\frac{a+b}{2}\right)^2-\left(\frac{a-b}{2}\right)^2$$
Conversely, given $n=t^2-s^2$, $n$ can be factored as $n=(t+s)(t-s)$.

If $n=ab$ with $a$ and $b$ close together, then $s=\frac{a-b}{2}$ is small, and so $t=\frac{a+b}{2}$ is slightly larger than $\sqrt{n}$. In that case, we can find $a$ and $b$ by trying all values of $t$ starting with $[\sqrt{n}]+1$, until we find one for which $t^2-n=s^2$ is a perfect square. This method is called the Fermat factorization.

Example. Factor $200819$.

Solution. $[\sqrt{200819}]+1=449$.
\begin{align*} 449^2-200819&=782,\ \mbox{not a perfact square}\\ 450^2-200819&=1681=41^2 \end{align*}
Hence,
$$200819=450^2-41^2=(450+41)(450-41)=491\cdot 409$$

If $a$ and $b$ are not close together for $n=ab$, although the Fermat factorization method will eventually find $a$ and $b$, one will have to try a large number of $t=[\sqrt{n}]+1, [\sqrt{n}]+2,\cdots$. There is a generalization of Fermat factorization. Choose small $k$, successively set $t=[\sqrt{kn}]+1, [\sqrt{kn}]+2,\cdots$, etc. until we find a $t$ for which $t^2-kn=s^2$ is a perfect square. Then $(t+s)(t-s)=kn$ and a nontrivial common factor of $n$ can be found by calculating $(t+s,n)$.

Example. Factor 141467.

Solution. First we factor 141467 by simple Fermat factorization. Since $\sqrt{141467}\approx 376.1209911717239$, we successively try $t=377,378,\cdots$ until we find $t^2-141467$ is a perfect square. We find
$$t^2-141467=414^2-141467=29929=173^2$$
Thus,
$$141467=414^2-173^2=(414+173)(414-173)=587\cdot 241$$
Now, this time we factor 141467 by generalized Fermat factorization with $k=3$. We try $t=[\sqrt{3n}]+1=652, 653,\cdots$ and find
$$t^2-3n=655^2-3\cdot 141467=4624=68^2=s^2$$
$(t+s,n)=(723,141467)=241$ and so $141467=241\cdot 587$.

The following two propositions tell us that it is not a good idea to choose an even number for $k$ in generalized Fermat factorization.

Proposition. If $k=2$, or if $k$ is any integer divisible by 2 but not by 4, then we cannot factor a large odd $n$ using generalized Fermat factorization with this choice of $k$.

Proof. $n$ is an odd integer, so $n=2m+1$ for some $m\in\mathbb{Z}$. For $k=2$, $kn=4m+2\equiv 4\mod 4$. If $k$ is an integer divisible by 2 but not by 4, $k=2(2l+1)=4l+2$ for some $l\in\mathbb{Z}$ and $kn\equiv 2\mod 4$. $t^2-s^2=kn\equiv 2\mod 4$, but the difference of two squares cannot be 2 modulo 4. (Each of the integers $t$ and $s$ is one of the forms $4u$, $4u+1$, $4u+2$, or $4u+3$ for some $u\in\mathbb{Z}$, so there are 16 different cases of $t$ and $s$ and in each case you can easily check if $t^2-s^2$ can be 2 modulo 4. For instance if $t=4u_1+1$ and $s=4u_2+2$, then
$t^2-s^2\equiv 1^2-2^2\equiv 1\mod 4$.)

Proposition. If $k=4$, and if generalized Fermat factorization works for a certain $t$, then simple Fermat factorization (with $k=1$) would have worked equally well.

Proof. $t^2-s^2=kn=4n\equiv 4\mod 8$ which can hold only if both $s$ and $t$ are even. In that case, $\left(\frac{s}{2}\right)^2=\left(\frac{t}{2}\right)^2-n$, so simple Fermat factorization would have worked equally well.

The Monte Carlo method is a computational simulation scheme, originally introduced by Stanisław Ulam, that solves a wide variety of problems arising in chemistry, economics, finance, mathematics, physics, etc. In this note, we discuss an application of the Monte Carlo method in factoring of integers. It is also called rho method and was introduced by J. M. Pollack.

The first step is to choose an easily evaluated map from $\mathbb{Z}_n$ to itself. A popular choice is $f(x)=x^2+1$. Next, one chooses an initial value $x=x_0$, and then computes the successive iterates of $f$: $x_1=f(x_0)$, $x_2=f(x_1)$, $x_3=f(x_2)$, etc. i.e. $x_{j+1}=f(x_j)$, $j=0,1,2,\cdots$. Compare the $x_j$’s, hoping to find two which are different modulo $n$ but the same modulo some divisor of $n$. Once we find such $x_i$, $k_k$, we have $(x_k-x_j,n)$ equal to a proper divisor of $n$

Example. Let us factor $91$ by choosing $f(x)=x^2+1$, $x_0=1$.
\begin{align*} x_1&=f(x_0)=2\\ x_2&=f(2)=5\\ x_3&=f(5)=26\\ &\vdots \end{align*}
Since $(x_3-x_2,n)=(21,91)=7$, 7 is a factor.

The method works by successively computing $x_k=f(x_{k-1})$ and comparing $x_k$ with the earlier $x_j$ until we find a pair satisfying $(x_k-x_j,n)=r>1$. But as $k$ gets larger, it becomes more time consuming to compute $(x_k-x_j,n)$ for each $j<k$. Note that once there is a pair $(k_0,j_0)$ such that $x_{k_0}\equiv x_{j_0} \mod r|n$, we have the same relation $x_k\equiv x_j\mod r$ for any pair $(j,k)$ such that $k-j=k_0-j_0$: Set $k=k_0+m$ and $j=j_0+m$, and apply $f$ to both sides of the congruence $x_{k_0}\equiv x_{j_0}\mod r$ repeatedly $m$ times.

The previous algorithm can be modified so that we need to calculate the gcd only once for each $k$. This significantly reduces the required computational burden. Here is the modified algorithm.

We successively compute the $x_k$. For each $k$, we proceed as follows. Suppose $k$ is an $(h+1)$-bit integer, i.e. $2^h\leq k<2^{h+1}$. Let $j$ be the largest $h$-bit integer: $j=2^h-1$. We compare $x_k$ with this particular $x_j$, i.e. compute $(x_k-x_j,n)$. If this gcd gives a nontrivial factor of $n$, stop. Otherwise continue on to $k+1$.

Example. $n=91$, $f(x)=x^2+1$, $x_0=1$.
\begin{align*} x_1&=f(1)=2\\ x_2&=f(2)=5;\ (x_2-x_1,n)=(5-2,91)=1\\ x_3&=f(5)=26;\ (x_3-x_1,n)=(24,91)=1\\ x_4&=f(26)=26^2+1\equiv 40\mod 91;\ (x_4-x_3,n)=(14,91)=7 \end{align*}

Example. Factor 4087 using $f(x)=x^2+x+1$ and $x_0=2$.
\begin{align*} x_1&=f(2)=7;\ (x_1-x_0,n)=(7-2,4087)=1\\ x_2&=f(7)=57;\ (x_2-x_1,n)=(57-7,4087)=1\\ x_3&=f(57)=3307;\ (x_3-x_1,n)=(3307-7,4087)=1\\ x_4&=f(3307)\equiv\mod 4087;\ (x_4-x_3,n)=(2745-3307,4087)=1\\ x_5&=f(2745)\equiv 1343\mod 4087; (x_5-x_3,n)=(1343-3307,4087)=1\\ x_6&=f(1343)\equiv 2626\mod 4087;\ (x_6-x_3, n)=(2626-3307,4087)=1\\ x_7&=f(2626)\equiv 3734\mod 4087;\ (x_7-x_3,n)=(3734-3307,4087)=61 \end{align*}
Hence, $61$ is a factor of $4087$ and $4087=61\cdot 67$.

Here, we study some computing problems that are related to Divisibility and the Division Algorithm.

Maxima

GCD: Maxima has a function to compute the GCD (Greatest Common Divisor) gcd(a,b).

Example. Find the GCD of 158 and 188.

(%i1) gcd(158,188);

(%o1) 2

LCM: Maxima has a function to compute the LCM (Least Common Multiple) lcm(express_1,...express_n). This function needs to be loaded by running the command load("functs").

Example. Find the LCD of 6, 8, 15, 32.

(%i1) load("functs")$

(%i2) lcm(6,8,15,32);

(%o2) 480

Bézout’s Lemma (see here) says that given two integers $a,b$, there is an integer solution $(x,y)$ to the equation $ax+by=(a,b)$. The Maxima command gcdex implements the Euclidean algorithm. gcdex(a,b) returns [x, y, u] where $u=(a,b)$ and $ax+by=u$. So, one can find a solution of the equation $ax+by=(a,b)$ using gcdex.

Example. Solve $2=158x+188y$.

We saw earlier that $(158,188)=2$.

(%i2) gcdex(158,188);

(%o2) [25, - 21, 2]

Hence, $x=25$ and $y=-21$.

One can also use the command igcdex. gcdex and igcdex both implement the Euclidean algorithm. The only difference between them is that the argument of igcdex must be integers while the arguments of gcdex can be also polynomials. For instance,

(%i7) gcdex(x^2 + 1,x^3 + 4);

(%o7)/R/ [$\displaystyle\frac{x^2+4x-1}{17}$, $\displaystyle\frac{x+4}{17}$, 1]

Python

GCD: Algorithm

1. The name of the Python code is gcd(a,b) and it implements the Euclidean algorithm to find the GCD of two integers $a$ and $b$.
2. If b == 0 then gcd(a,b) returns 0 because every nonzero integer divides 0.
3. Otherwise, it returns gcd(b,a % b). % is the modulo or remainder operator. a % b returns the remainder when $a$ is divided by $b$ if $a\geq b$ or $a$ if $a<b$ because $a=b\cdot 0+a$. For example, 12345 % 2 returns 1.
4. The loop continues until b==0.

Use your favorite editor to write the following python code (it is called a Python module) and save it as, for example, gcd.py. The Python interface IDLE has own editor. I also recommend Emacs and Geany for Python editors if you are running Python in command shell with iPython.

# Computing the GCD of two integers a and b recursively using the Euclidean algorithm

def gcd(a,b):
    if (b == 0):
        return a
    return gcd(b, a%b)

On IDLE, the Python module can be run from its pull-down (or drop-down) menu under Run. Its shortcut key is F5. To run it on iPython, first change your current folder location to the one that contains the Python module gcd.py using cd (change directory) command in command shell and then start iPython by typing ipython (you can first start iPython and then change your current folder location as well using cd within iPython) and pressing the enter key. Once iPython is loaded, do

$ ipython
Python 3.7.9 (default, Feb 4 2021, 01:17:56)
Type 'copyright', 'credits' or 'license' for more information
IPython 7.19.0 -- An enhanced Interactive Python. Type '?' for help.
In [1]: run gcd.py
In [2]: gcd(98,56)
Out[2]: 14

The Euclidean Algorithm: The Euclidean algorithm itself is a perfect computer algorithm. A Python code of the Euclidean algorithm is

# This Python module runs the Euclidean Algorithm to find the GCD of two integers a and b
def euclid(a,b):
    r = a % b
    while r > 0:
        print (a, b, r)
        a, b, r = b, r, b % r
    return b

Running this Python module (named euclid.py) on iPython shows

$ ipython
Python 3.7.9 (default, Feb 4 2021, 01:17:56)
Type 'copyright', 'credits' or 'license' for more information
IPython 7.19.0 -- An enhanced Interactive Python. Type '?' for help.
In [1]: run euclid.py
In [2]: euclid(98,56)
98 56 42
56 42 14
Out[2]: 14

LCM: The LCM $[a,b]$ can be obtained by $\frac{ab}{(a,b)}$ (see here).

Algorithm: The name of the Python code is lcm(a,b). It computes the LCM of two integers $a$ and $b$.

1. Find the GCD of $a$ and $b$ by running the Python program gcd(a,b).
2. Compute the LCM of $a$ and $b$ using the formula a*b/gcd(a,b).

Its Python code is

# Compute the LCM of two integers a and b
Recursive function to return the GCD of a and b

def gcd(a,b):
    if (b == 0):
        return a
    return gcd(b, a%b)

# Function to return the LCM of a and b

def lcm(a,b):
    return a*b/gcd(a,b)

Running this Python code (named lcm.py) on iPython shows

In [3]: run lcm.py
In [4]: lcm(15,20)
Out[4]: 60.0

Linear Diophantine Equations: As seen here, a linear diophantine equation
\begin{equation}
\label{eq:lindiophant}
ax+by=c
\end{equation}
is solvable if and only if $(a,b)|c$. If $b|a$ then $(a,b)=b$ and
$$x=0,\ y=\frac{c}{b}$$
is a solution. If $b\not|a$, then write $a=bq+r$ and substitute into \eqref{eq:lindiophant}:
$$(bq+r)x+by=c$$
This can be rearranged as
$$b(qx+y)+rx=c$$
Setting $u=qx+y$ and $v=x$, we obtain the equation
\begin{equation}
\label{eq:lindiophant2}
bu+rv=c
\end{equation}
with smaller coefficients. If we can solve \eqref{eq:lindiophant2}, we recover a solution of the original equation as
$$x=v,\ y=u-qv$$
The process terminates (according to the Euclidean algorithm) with an equation where the second coefficient divides the first. This scheme is embodied in the following Python code

# Solve a linear diophantine equation using the Euclidean algorithm

# The function divmod returns a pair of numbers, the quotient and remainder, which we store in q, r.

def diosolve(a,b,c):
    q, r = divmod(a,b)
    if r == 0:
        return( [0,c/b] )
    else:
        sol = diosolve( b, r, c )
        u = sol[0]
        v = sol[1]
        return( [ v, u-q*v ] )

Running this Python code (named diosolve.py) on iPython to solve the linear diophantine equation $12345x+54321y=3$:

In [9]: run diosolve.py
In [10]: diosolve(12345, 54321, 3)
Out[10]: [3617.0, -822.0]

Note that $(12345,54321)=3$, so of course the diophantine equation is solvable.

Refernces.

[1] Jim Carlson, A Short Course in Python for Number Theory, Draft of May 21, 2004

[2] Maximal 5.29.0 Manual, 29. Number Theory

We begin by discussing a method of computing $a^d\mod n$, called the repeated squaring method, which is efficient enough to require only a handheld calculator using the following example.

Example. Let us compute the least residue of $848^{187}\mod 1189$ by using only a handheld calculator. First we convert 187 to the base 2:
\begin{align*} 187&=(10111011)_2\\ &=1\cdot 2^7+0\cdot 2^6+1\cdot 2^5+1\cdot 2^4+1\cdot 2^3+0\cdot 2^2+1\cdot 2+1\\ &=128+32+16+8+2+1 \end{align*}
Let $a=848$ and $n=1189$.
$$\begin{array}{|c|c|}
\hline
k & a^k\mod n\\
\hline
1 & 848\\
\hline
2 & 848^2=719104\equiv 948\mod 1189\\
\hline
4 & 948^2=898704\equiv 1009\mod 1189\\
\hline
8 & 1009^2= 1018081\equiv 297\mod 1189\\
\hline
16 & 297^2=88209\equiv 223\mod 1189\\
\hline
32 & 223^2=49729\equiv 980\mod 1189\\
\hline
64 & 980^2=960400\equiv 877\mod 1189\\
\hline
128 & 877^2=769129\equiv 1035\mod 1189\\
\hline
\end{array}$$
Hence,
\begin{align*} 848^{187}&=848^{128+32+16+8+2+1}\\ &=848^{128}\cdot 848^{32}\cdot 848^{16}\cdot 848^8\cdot 848^2\cdot 848\\ &\equiv 1035\cdot 980\cdot 223\cdot 297\cdot 948\cdot 848\mod 1189\\ &\equiv 190\mod 1189 \end{align*}

Here is our first primality test. Let us recall Fermat Little Theorem: If $n$ is a prime and $n\not|a$, then $a^{n-1}\equiv 1\mod n$. The contrapositive of the above statement can serve as a primality test: if $a^{n-1}\not\equiv 1\mod n$, then $n$ is not a prime or $n|a$.

Example. Let us show that $33$ is not a prime using the above primality test. (Of course it is not prime because $33=3\cdot 11$.) $33\not|2$ and
\begin{align*} 2^{33-1}&=2^{32}\\ &=(2^5)^62^2\\ &=(32)^62^2\\ &\equiv(-1)^62^2\mod 33\\ &\equiv 4\mod 33 \end{align*}
Thus, $33$ is not a prime.

Here is another primality test.

Theorem. If the integer $n>1$ has no prime divisor less than or equal to $\sqrt{n}$, then $n$ is prime.

Proof. Suppose that $n$ is a composite number. Then $n=d_1d_2$ for some $d_1>1$ and $d_2>1$. Suppose that all prime divisors are greater than $\sqrt{n}$. Then $d_1>\sqrt{n}$ and $d_2>\sqrt{n}$. But $n=d_1d_2>(\sqrt{n})^2=n$, which is impossible. This proves the theorem. If $d_1\leq\sqrt{n}$, then $d_1$ is either prime or else has a prime divisor less than or equal to $\sqrt{n}$.

Question: Is the statement “If $2^{n-1}\equiv 1\mod n$ then $n$ is a prime.” true? The ancient Chinese believed so, but there is a counterexample. $2^{340}\equiv 1\mod 341$. But, $341=11\cdot 31$ (this is left as an exercise), so it is not a prime.

Definition. We call $n$ a pseudoprime if $2^{n-1}\equiv 1\mod n$ but $n$ is composite. More generally, a composite number $n$ such that $a^{n-1}\equiv 1\mod n$ is called a pseudoprime to base $a$.

The smallest pseudoprime is $341$, and it was not discovered until 1819. Bases other than $2$ can be used to identify composite numbers. For example,
$$3^{340}\equiv 56\mod 341$$
Thus, $341$ is not a prime number. These pseudoprimes are rarer then primes. Ever rarer are pseudoprimes to multiples bases. It is known, for example, that there are only 1770 integers below $25\cdot 10^9$ that are simultaneously pseudoprimes to the bases $2$, $3$, $5$ and $7$. So, the primality of numbers less than $25\cdot 10^9$ can be determined by testing Fermat’s congruence $a^{n-1}\equiv 1\mod n$ with these four bases, then comparing any number passing all 4 tests with a list of the 1770 exceptions. If the number is not any of those 1770 exceptions, it must be a prime.

Next we discuss Mersenne numbers and Fermat numbers, and primality tests those numbers.

Definition. If $k$ is a positive integer, we call $M_k=2^k-1$ a Mersenne number.

$$\begin{array}{|c|c|c|}
\hline
k & M_k=2^k-1 & \mbox{prime?}\\
\hline
1 & 1 & \mbox{N}\\
\hline
2 & 3 & \mbox{Y}\\
\hline
3 & 7 & \mbox{Y}\\
\hline
4 & 15=3\cdot 5 & \mbox{N}\\
\hline
5 & 31 & \mbox{Y}\\
\hline
6 & 63=9\cdot 7 & \mbox{N}\\
\hline
7 & 127 & \mbox{Y}\\
\hline
8 & 255=3\cdot 5\cdot 7 & \mbox{N}\\
\hline
9 & 511=7\cdot 73 & \mbox{N}\\
\hline
10 & 1023=3\cdot 11\cdot 31 & \mbox{N}\\
\hline
\end{array}$$
From this table, one may expect that $M_k$ is prime when $k$ is prime. However, $M_{11}=2047=23\cdot 89$ is not a prime.

Theorem.
If $d|k$, then $M_d|M_k$, so if $M_k$ is prime, then $k$ is prime.

Proof. Since $d|k$, $k=dd_1$ for some $d_1\in\mathbb{Z}$.
$$M_k=2^k-1=2^{dd_1}-1=(2^d-1)(2^{(d(d_1-1)}+2^{(d(d_1-2)}+\cdots+2^d+1)$$

$$\begin{array}{|c|c|c|}
\hline
k & 2^k+1 & \mbox{prime?}\\
\hline
1 & 3 & \mbox{Y}\\
\hline
2 & 5 & \mbox{Y}\\
\hline
3 & 9=3\cdot 3 & \mbox{N}\\
\hline
4 & 17 & \mbox{Y}\\
\hline
5 & 33=3\cdot 11 & \mbox{N}\\
\hline
6 & 65=5\cdot 13 & \mbox{N}\\
\hline
7 & 129=3\cdot 43 & \mbox{N}\\
\hline
8 & 257 & \mbox{Y}\\
\hline
9 & 513=3\cdot 3\cdot 3\cdot 19 & \mbox{N}\\
\hline
10 & 1025=5\cdot 5\cdot 41 & \mbox{N}\\
\hline
\end{array}$$

Theorem. If $k$, $a$, and $b$ are positive integers such that $k=ab$, where $a$ is odd, then
$$2^b+1|2^k+1$$
In particular, if $2^k+1$ is prime, then $k$ is $0$ or a power of $2$.

Proof. Under the assumption, we want to show that
$$2^k+1\equiv 0\mod 2^b+1$$
or
$$2^k\equiv -1\mod 2^b+1$$
Since, $2^b\equiv -1\mod 2^b+1$,
$$2^k=2^{ab}=(2^b)^a\equiv (-1)^a\mod 2^b+1\equiv -1\mod 2^b+1$$
If $k>0$ is not a power of $2$, then we can take $a>1$ so that
$$1<2^b+1<2^k+1$$
This means that $2^k+1$ has a positive divisor other than 1 and itself.

Definition. $F_r=2^{2^r}+1$ is called a Fermat number.

Example. $F_0=3$, $F_1=5$, $F_2=17$, $F_3=257$, $F_4=65537$ are primes, so Fermat conjectured that Fermat numbers are prime. But $F_5=2^{32}+1$ is not a prime. It is divisible by $641$. $641$ can be written
\begin{align*} 641&=640+1=2^7\cdot 5+1\\ 641&=625+16=5^4+2^4 \end{align*}
Thus, we have $2^7\cdot 5\equiv -1\mod 641$ and $2^4\equiv -5^4\mod 641$.
Now,
\begin{align*} F_5&=2^{32}+1\\ &=(2^7)^4\cdot 2^4+1\\ &\equiv(2^7)^4(-5^4)+1\mod 641\\ &\equiv -(2^7\cdot 5)^4+1\mod 641\\ &\equiv 0\mod 641 \end{align*}
No other Fermat numbers that are prime have been found yet.

One may expect that for any composite number $n$, there exists a base $a$ for which Fermat’s theorem can be used to show that $n$ is composite. However, there are composite numbers called Carmichael numbers, which are pseudo primes to every base. That is, $n$ is composite, but $a^{n-1}\equiv 1\mod n$, whenever $(a,n)=1$. The smallest Carmichael number is $561=3\cdot 11\cdot 17$. It was proved in 1994 that there are infinitely many Carmichael numbers.

There is a primality test just for Mersenne numbers. It is called the Lucas-Lehmer test. Define a sequence $S_1,S_2,\cdots$ by
\begin{align*} S_1&=4,\\ S_n&=S_{n-1}^2-2\ \mbox{for}\ n>1. \end{align*}
For example, $S_2=4^2-2=14$, $S_3=14^2-2=194$. Then the Lucas-Lehmer test is

Theorem. If $p$ is an odd prime then $M_p=2^p-1$ is prime if and only if $S_{p-1}\equiv 0\mod M_p$.

For example, take $p=7$ and show that $M_p=2^7-1=127$ is a prime.
\begin{align*} S_3&=194\equiv 67\mod 127\\ S_4&\equiv 67^2-2\mod 127\equiv 42\mod 127\\ S_5&\equiv 42^2-2\mod 127\equiv 111\mod 127\\ S_6&\equiv 111^2-2\mod 127\equiv 0\mod 127. \end{align*}

There is also a primality test just for Fermat numbers. It is called Pépin’s test.

Theorem. If $n>0$, the Fermat number $F_n=2^{2^n}+1$ is prime if and only if $3^{\frac{F_n-1}{2}}\equiv -1\mod F_n$.

Example. Use Pépin’s test to show that $F_3=257$ is prime.

Solution. $3^{\frac{257-1}{2}}=3^{128}=3^{2^7}$.
\begin{align*} 3^2&=9\\ 3^4&=9^2=81\\ 3^8&=81^2\equiv 136\mod 257\\ 3^{16}&\equiv 136^2\mod 257\equiv 249\mod 257\\ 3^{32}&\equiv 249^2\mod 257\equiv 64\mod 257\\ 3^{64}&\equiv 64^2\mod 257\equiv 241\mod 257\\ 3^{128}&\equiv 241^2\mod 257\equiv 256\mod 257\equiv -1\mod 257 \end{align*}