p=131;a=2;b=29;m=Ceiling[(Sqrt[p-1])]
d=PowerMod[a,m,p]

Do[If[Mod[b*PowerMod[a,-r,p],p]==1,Print[" r = ",r]],{r,0,m-1}]

Do[If[Mod[b*PowerMod[a,-r,p],p]==PowerMod[d,q,p],Print["x=",q*m+r]],{r,0,m-1},{q,1,m}]