Some functions related to primes 
Math 302 Cryptography
February 26, 2025

T. Ratliff, Wheaton College
ratliff_thomas@wheatoncollege.edu

A mrWitness[ ] function

mrWitness[a_,n_]:=Module[{temp,k,q,atable,i},

(* Factor out the powers of 2 from n-1 and find q *)
temp=n-1;
k=0;
While[Mod[temp,2]==0,temp=temp/2;k++];
q=(n-1)/2^k;

(* If a^q = 1 mod n, then not a MR witness *)
If[PowerMod[a,q,n]==1,
Return[False]
];

(* Now check second condition *)
atable = Table[ PowerMod[a, (2^i) q, n], {i,0,k-1}];
(*Print[atable];  You can uncomment this if you want to see the values for Condition (ii) *) 
If[ MemberQ[atable, n-1]==True, 
Return[False], 
(* else *) Return[True]
];

]

A primeQ[ ] function that calls mrWitness[ ]

primeQ[n_]:=Module[{isPrime, i},

isPrime=True;
i=0;

(* If 10 values are NOT MR witnesses, then strong evidence that n is a pseudo-prime *)
While[ (isPrime == True) && (i<10),
a = RandomInteger[{2,n-1}];
If[ mrWitness[a,n]==True, isPrime=False, i++]
];
Return[isPrime]
]

A Pollard’s ρ  pFactor[ ] function

This function will apply the algorithm discussed in class to find one factor of n. 
I hope it’s fairly self-explanatory.

pFactor[n_]:=Module[{x=1,y=1,g=1,i=1, f},
f[x_]:=Mod[x^2+1,n];

While[g==1, 
x=f[x];
y=f[f[y]];
g=GCD[Abs[x-y],n];
(*Print[i," ", x," ", y, " ", g];  You can uncomment this to see each stage in the sequences*)
i++
];
Print["Collision at i = ",i-1,", x = ",x,", y = ",y,", g = ", g]
]

Let’s look for a k-bit prime

(* A function to tell how many values we should expect to pick before getting a k-bit prime *)
expectedNum[k_]:=Ceiling[N[ (2^k-2^(k-1))/(2^k/Log[2^k] - 2^(k-1)/Log[2^(k-1)])]]

How many values do we expect to need to pick, on average, before finding a 20-bit prime?

expectedNum[20]

Now let’s look for that prime

(* k bit prime *)
k=20

(* Count how many times we need to look for a prime *)
j=0;
gotOne=False;

While[ gotOne==False, 
n=RandomInteger[{2^(k-1), 2^k-1}]; (* Pick a random k-bit integer *)
j++;
If[ primeQ[n]==True, 
gotOne=True]
];

Print["Pseudoprime = ", n, ", found in ", j," steps"]

Let’s look for Diffie-Hellman parameters (p, α) where p is a k-bit prime

(* bits in p *)
k=20

(* Count how many times we need to look for a prime *)
j=0;
gotOne=False;

While[ gotOne==False, 
(* Look for a (k-1)-bit prime q *)
q=RandomInteger[{2^(k-2), 2^(k-1)-1}];
j++;
If[ primeQ[q]==True, 
If[primeQ[2*q+1]==True, (* Check if 2q + 1 is prime*)
gotOne=True]
];
];

Print["Pseudoprime = ", 2*q+1, ", found in ", j," steps"]

Use Pollard’s ρ for factoring

n=2465;
pFactor[n]