This is a factorization algorithm similar to Pollard p-1 factorization,
but it is faster, because there is no use of a gcd in every loop.
The algorithm has the same limits as the Pollard p-1 algorithm.
You can combine it with the Helmes p+1 algorithm in the adjoined field.
Implementation in Mupad
powermod_I:=proc (zahl, exp, modul)
begin
bit:=exp mod 2;
res_quadrat:=zahl;
if bit = 1 then
res:=zahl;
else
res:=1;
end_if;
exp:=(exp-bit)/2;
while TRUE do
bit:=exp mod 2;
res_quadrat:=res_quadrat*res_quadrat;
res_re:=Re (res_quadrat) mod modul;
res_im:=Im (res_quadrat) mod modul;
res_quadrat:=res_re+res_im*I;
if bit = 1 then
res:=res*res_quadrat;
res_re:=Re (res) mod modul;
res_im:=Im (res) mod modul;
res:=res_re+res_im*I;
end_if;
exp:=(exp-bit)/2;
if exp=0 then
return (res);
end_if;
end_while;
end;
while TRUE do
bb:=b;
ee:=e;
for i from 1 to 1000 do
b:=powermod_I (b, e, p);
e:=e+1;
end_for;
if Im (b)=0 then
// Repetition
for i from 1 to 1000 do
bb:=powermod_I (bb, ee, p);
f:=gcd (Im (bb), p);
if f > 1 then
print (p,"=", f, "*", p/f);
return (0);
end_if;
ee:=ee+1;
end_for;
else
f:=gcd (Im (b), p);
if f > 1 then
print (p,"=", f, "*", p/f);
return (0);
end_if;
end_if;
end_while;
