Tonelli-Shanks algorithm: Difference between revisions

Content added Content deleted
m (Convert to a task.)
m (→‎{{header|Phix}}: bigatom -> mpfr)
Line 1,726: Line 1,726:
=={{header|Phix}}==
=={{header|Phix}}==
{{trans|C#}}
{{trans|C#}}
{{libheader|mpfr}}
<lang Phix>include bigatom.e
<lang Phix>include mpfr.e
function ba_mod_exp(object base, exponent, modulus)
-- base/exponent/modulus can be integer/string/bigatom
⚫
function ts(string ns, ps)
-- returns mod(power(base,exponent),modulus) [aka b^e%m], but in bigatoms and faster.
bigatom res = BA_ONE
mpz n = mpz_init(ns),
base = ba_mod(base,modulus)
p = mpz_init(ps),
⚫
t = mpz_init(),
while ba_compare(exponent,0)!=0 do
if ba_mod(exponent,2)=BA_ONE then
r = mpz_init(),
res = ba_mod(ba_multiply(res,base),modulus)
pm1 = mpz_init(),
end if
pm2 = mpz_init()
mpz_sub_ui(pm1,p,1) -- pm1 = p-1
base = ba_mod(ba_multiply(base,base),modulus)
mpz_fdiv_q_2exp(pm2,pm1,1) -- pm2 = pm1/2
exponent = ba_idivide(exponent,2)
mpz_powm(t,n,pm2,p) -- t = mod(n^pm2,p)
⚫
end while
⚫
if mpz_cmp_si(t,1)!=0 then
return res
⚫
return "No solution exists"
end function

⚫
function ts(object n, p)
bigatom pm1 = ba_sub(p,1),
⚫
pm1o2 = ba_div(pm1,2)
⚫
if ba_mod_exp(n,pm1o2, p)!=BA_ONE then
return {0,0,false}
end if
end if
bigatom q = pm1
mpz q = mpz_init_set(pm1)
integer ss = 0
integer ss = 0
while ba_mod(q,2)=BA_ZERO do
while mpz_even(q) do
ss += 1
ss += 1
q = ba_div(q,2)
mpz_fdiv_q_2exp(q,q,1) -- q/=2
end while
end while
if ss=1 then
if ss=1 then
bigatom r1 = ba_mod_exp(n,ba_div(ba_add(p,1),4),p)
mpz_add_ui(t,p,1)
return {r1,ba_sub(p,r1),true}
mpz_fdiv_q_2exp(t,t,2)
mpz_powm(r,n,t,p) -- r = mod(n^((p+1)/4),p)
⚫
end if
integer z = 2
else
mpz z = mpz_init(2)
while ba_mod_exp(z,pm1o2,p)!=pm1 do
z += 1
while true do
mpz_powm(t,z,pm2,p) -- t = mod(z^pm2,p)
⚫
end while
if mpz_cmp(t,pm1)=0 then exit end if
bigatom c = ba_mod_exp(z,q,p),
r = ba_mod_exp(n,ba_div(ba_add(q,1),2),p),
mpz_add_ui(z,z,1) -- z+= 1
⚫
t = ba_mod_exp(n,q,p)
⚫
integer m = ss
⚫
while true do
if t=BA_ONE then
return {r,ba_sub(p,r),true}
⚫
end if
⚫
integer i = 0
bigatom zz = t
⚫
while zz!=BA_ONE and i<m-1 do
zz = ba_mod(ba_mul(zz,zz),p)
⚫
i += 1
end while
end while
bigatom b = c
mpz {b,c,zz} = mpz_inits(3)
integer e = m-i-1
mpz_powm(c,z,q,p) -- c = mod(z^q,p)
while e>0 do
mpz_add_ui(t,q,1)
b = ba_mod(ba_mul(b,b),p)
mpz_fdiv_q_2exp(t,t,1)
e -= 1
mpz_powm(r,n,t,p) -- r = mod(n^((q+1)/2),p)
mpz_powm(t,n,q,p) -- t = mod(n^q,p)
⚫
integer m = ss
while mpz_cmp_si(t,1) do -- t!=1
⚫
integer i = 0
⚫
mpz_set(zz,t)
⚫
while mpz_cmp_si(zz,1)!=0 and i<m-1 do
mpz_powm_ui(zz,zz,2,p) -- zz = mod(zz^2,p)
⚫
i += 1
⚫
end while
mpz_set(b,c)
integer e = m-i-1
⚫
while e>0 do
mpz_powm_ui(b,b,2,p) -- b = mod(b^2,p)
e -= 1
⚫
end while
⚫
mpz_mul(r,r,b)
mpz_mod(r,r,p) -- r = mod(r*b,p)
mpz_powm_ui(c,b,2,p) -- c = mod(b^2,p)
⚫
mpz_mul(t,t,c)
mpz_mod(t,t,p) -- t = mod(t*c,p)
⚫
m = i
end while
end while
⚫
end if
r = ba_mod(ba_mul(r,b),p)
c = ba_mod(ba_mul(b,b),p)
mpz_sub(p,p,r)
return mpz_get_str(r)&" and "&mpz_get_str(p)
⚫
t = ba_mod(ba_mul(t,c),p)
⚫
m = i
end while
end function
end function
constant tests = {{10,13},
constant tests = {{"10","13"},
{56,101},
{"56","101"},
{1030,10009},
{"1030","10009"},
{1032,10009},
{"1032","10009"},
{44402,100049},
{"44402","100049"},
{665820697,1000000009},
{"665820697","1000000009"},
{881398088036,1000000000039},
{"881398088036","1000000000039"},
{"41660815127637347468140745042827704103445750172002",
{"41660815127637347468140745042827704103445750172002",
ba_sprint(ba_add(ba_power(10,50),577))}}
sprintf("1%s577",repeat('0',47))}} -- 10^50+577

for i=1 to length(tests) do
for i=1 to length(tests) do
object {p1,p2} = tests[i]
string {p1,p2} = tests[i]
sequence sol = ts(p1,p2)
printf(1,"For n = %s and p = %s, %s\n",{p1,p2,ts(p1,p2)})
if not string(p1) then p1 = sprintf("%d",p1) end if
if not string(p2) then p2 = sprintf("%d",p2) end if
if sol[3] then
sol = ba_sprint(sol[1])&" and "&ba_sprint(sol[2])
else
⚫
sol = "No solution exits"
end if
printf(1,"For n = %s and p = %s, %s\n",{p1,p2,sol})
end for</lang>
end for</lang>
{{out}}
{{out}}
Line 1,816: Line 1,810:
For n = 56 and p = 101, 37 and 64
For n = 56 and p = 101, 37 and 64
For n = 1030 and p = 10009, 1632 and 8377
For n = 1030 and p = 10009, 1632 and 8377
For n = 1032 and p = 10009, No solution exits
For n = 1032 and p = 10009, No solution exists
For n = 44402 and p = 100049, 30468 and 69581
For n = 44402 and p = 100049, 30468 and 69581
For n = 665820697 and p = 1000000009, 378633312 and 621366697
For n = 665820697 and p = 1000000009, 378633312 and 621366697
Retrieved from "https://rosettacode.org/wiki/Tonelli-Shanks_algorithm"