Some DLP functions 
Math 302 Cryptography
February 12, 2025

T. Ratliff, Wheaton College
ratliff_thomas@wheatoncollege.edu

A Pollard’s ρ function

f[x_,g_,h_,p_]:=  Mod[Piecewise[
{{g*x, 0 ≤ x < p/3},
{x^2, p/3 ≤ x < 2p/3},
{h*x, 2p/3 ≤ x < p}}] , p];

nextAlpha[ai_,xi_,p_]:= Mod[Piecewise[
{{ai+1, 0 ≤ xi < p/3},
{2ai, p/3 ≤ xi < 2p/3},
{ai, 2p/3 ≤ xi < p}}] , p-1];

nextBeta[bi_,xi_,p_]:= Mod[Piecewise[
{{bi, 0 ≤ xi < p/3},
{2bi, p/3 ≤ xi < 2p/3},
{bi+1, 2p/3 ≤ xi < p}}] , p-1];

pollardsRho[g_,h_,p_]:= Module[{x=1,α=0,β=0,y=1,γ=0,δ=0,i,u,v,v2,w,d,ordg,options,intOptions,gpowers,xloc},

For[i=1,(i==1)||(x≠y),i++,

α=nextAlpha[α,x,p];
β=nextBeta[β,x,p];
x=f[x,g,h,p];

γ=nextAlpha[γ,y,p];
δ=nextBeta[δ,y,p];
y=f[y,g,h,p];

γ=nextAlpha[γ,y,p];
δ=nextBeta[δ,y,p];
y=f[y,g,h,p];

If[Mod[i,100000]==0, Print[i]];
];

Print["\nCOLLISION AT\ni = ",i-1,":    (Subscript[x, i]) = ",x,", (Subscript[α, i]) = ",   α,", (Subscript[β, i]) = ",β,":     (Subscript[y, i]) = ",y,", (Subscript[γ, i]) = ",γ,", (Subscript[δ, i]) = ",δ] ;

u=Mod[α-γ,p-1];
v=Mod[δ-β,p-1];

v2 = ExtendedGCD[v,p-1][[2,1]];

w=Mod[u*v2,p-1];
d=Mod[v*v2,p-1];

Print["\nu = ", u,", v = ", v,",   g^",w," ≡ h^",d," mod ",p];

If[d==1,
Print["\nSolution x =",w],
(*else d != 1 and we need to take a root mod p*)
ordg = MultiplicativeOrder[g,p];
options = Table[1/d (w- k *ordg),{k, Ceiling[w/ordg -d], Floor[w/ordg]}];
Print["Options for x: ", options];
intOptions = Select[options,IntegerQ];
gpowers = PowerMod[g,intOptions,p];
Print[TableForm[Transpose[{intOptions,gpowers}]]];
xloc = Position[gpowers,h][[1,1]];
Print["\nSolution x = ", intOptions[[xloc]]]
]
];

A Shanks Babystep-Giantstep function

shanks[g_,h_,p_]:=Module[{n,babySteps, giantSteps,i,m,math,collision,bsPos, gsPos,x},
n=Floor[Sqrt[MultiplicativeOrder[g,p]]]+1;
babySteps = Table[ PowerMod[g,i,p],{i,0,n}];
giantSteps = Table[ Mod[h*PowerMod[g,-n*j,p],p],{j,0,n}];
match = Intersection[babySteps,giantSteps];
collision = match[[1]];
bsPos = Position[babySteps,collision];
gsPos = Position[giantSteps,collision];
i=bsPos[[1,1]] - 1;
j =gsPos[[1,1]] - 1;

x=i+j*n;
Print["Lists are n + 1 = ", n+1," long\nSolution x = ", x]
]

Using the functions

pollardsRho[6965,5254,21599]

shanks[6965,5254,21599]