Tonelli-Shanks algorithm: Difference between revisions
Content added Content deleted
Thundergnat (talk | contribs) m (→{{header|Perl 6}}: Minor code simplifications, add missing failing test) |
(added FreeBASIC) |
||
Line 151: | Line 151: | ||
❌ error: not a square (mod p) (1032 10009) |
❌ error: not a square (mod p) (1032 10009) |
||
</pre> |
</pre> |
||
=={{header|FreeBASIC}}== |
|||
<lang FreeBASIC>' version 11-04-2017 |
|||
' compile with: fbc -s console |
|||
' maximum for p is 17 digits to be on the save side |
|||
' TRUE/FALSE are built-in constants since FreeBASIC 1.04 |
|||
' But we have to define them for older versions. |
|||
#Ifndef TRUE |
|||
#Define FALSE 0 |
|||
#Define TRUE Not FALSE |
|||
#EndIf |
|||
Function mul_mod(a As ULongInt, b As ULongInt, modulus As ULongInt) As ULongInt |
|||
' returns a * b mod modulus |
|||
Dim As ULongInt x, y = a Mod modulus |
|||
While b > 0 |
|||
If (b And 1) = 1 Then |
|||
x = (x + y) Mod modulus |
|||
End If |
|||
y = (y Shl 1) Mod modulus |
|||
b = b Shr 1 |
|||
Wend |
|||
Return x |
|||
End Function |
|||
Function pow_mod(b As ULongInt, power As ULongInt, modulus As ULongInt) As ULongInt |
|||
' returns b ^ power mod modulus |
|||
Dim As ULongInt x = 1 |
|||
While power > 0 |
|||
If (power And 1) = 1 Then |
|||
' x = (x * b) Mod modulus |
|||
x = mul_mod(x, b, modulus) |
|||
End If |
|||
' b = (b * b) Mod modulus |
|||
b = mul_mod(b, b, modulus) |
|||
power = power Shr 1 |
|||
Wend |
|||
Return x |
|||
End Function |
|||
Function Isprime(n As ULongInt, k As Long) As Long |
|||
' miller-rabin prime test |
|||
If n > 9223372036854775808ull Then ' limit 2^63, pow_mod/mul_mod can't handle bigger numbers |
|||
Print "number is to big, program will end" |
|||
Sleep |
|||
End |
|||
End If |
|||
' 2 is a prime, if n is smaller then 2 or n is even then n = composite |
|||
If n = 2 Then Return TRUE |
|||
If (n < 2) OrElse ((n And 1) = 0) Then Return FALSE |
|||
Dim As ULongInt a, x, n_one = n - 1, d = n_one |
|||
Dim As UInteger s |
|||
While (d And 1) = 0 |
|||
d = d Shr 1 |
|||
s = s + 1 |
|||
Wend |
|||
While k > 0 |
|||
k = k - 1 |
|||
a = Int(Rnd * (n -2)) +2 ' 2 <= a < n |
|||
x = pow_mod(a, d, n) |
|||
If (x = 1) Or (x = n_one) Then Continue While |
|||
For r As Integer = 1 To s -1 |
|||
x = pow_mod(x, 2, n) |
|||
If x = 1 Then Return FALSE |
|||
If x = n_one Then Continue While |
|||
Next |
|||
If x <> n_one Then Return FALSE |
|||
Wend |
|||
Return TRUE |
|||
End Function |
|||
Function legendre_symbol (a As LongInt, p As LongInt) As LongInt |
|||
Dim As LongInt x = pow_mod(a, ((p -1) \ 2), p) |
|||
If p -1 = x Then |
|||
Return x - p |
|||
Else |
|||
Return x |
|||
End If |
|||
End Function |
|||
' ------=< MAIN >=------ |
|||
Dim As LongInt b, c, i, k, m, n, p, q, r, s, t, z |
|||
For k = 1 To 7 |
|||
Read n, p |
|||
Print "Find solution for n ="; n; " and p =";p |
|||
If legendre_symbol(n, p) <> 1 Then |
|||
Print n;" is not a quadratic residue" |
|||
Print |
|||
Continue For |
|||
End If |
|||
If p = 2 OrElse Isprime(p, 15) = FALSE Then |
|||
Print p;" is not a odd prime" |
|||
Print |
|||
Continue For |
|||
End If |
|||
s = 0 : q = p -1 |
|||
Do |
|||
s += 1 |
|||
q \= 2 |
|||
Loop Until (q And 1) = 1 |
|||
If s = 1 And (p Mod 4) = 3 Then |
|||
r = pow_mod(n, ((p +1) \ 4), p) |
|||
Print "Solution found:"; r; " and"; p - r |
|||
Print |
|||
Continue For |
|||
End If |
|||
z = 1 |
|||
Do |
|||
z += 1 |
|||
Loop Until legendre_symbol(z, p) = -1 |
|||
c = pow_mod(z, q, p) |
|||
r = pow_mod(n, (q +1) \ 2, p) |
|||
t = pow_mod(n, q, p) |
|||
m = s |
|||
Do |
|||
i = 0 |
|||
If (t Mod p) = 1 Then |
|||
Print "Solution found:"; r; " and"; p - r |
|||
Print |
|||
Continue For |
|||
End If |
|||
Do |
|||
i += 1 |
|||
If i >= m Then Continue For |
|||
Loop Until pow_mod(t, 2 ^ i, p) = 1 |
|||
b = pow_mod(c, (2 ^ (m - i -1)), p) |
|||
r = mul_mod(r, b, p) |
|||
c = mul_mod(b, b, p) |
|||
t = mul_mod(t, c, p)' t = t * b ^ 2 |
|||
m = i |
|||
Loop |
|||
Next |
|||
Data 10, 13, 56, 101, 1030, 10009, 1032, 10009, 44402, 100049 |
|||
Data 665820697, 1000000009, 881398088036, 1000000000039 |
|||
' empty keyboard buffer |
|||
While InKey <> "" : Wend |
|||
Print : Print "hit any key to end program" |
|||
Sleep |
|||
End</lang> |
|||
{{out}} |
|||
<pre>Find solution for n = 10 and p = 13 |
|||
Solution found: 7 and 6 |
|||
Find solution for n = 56 and p = 101 |
|||
Solution found: 37 and 64 |
|||
Find solution for n = 1030 and p = 10009 |
|||
Solution found: 1632 and 8377 |
|||
Find solution for n = 1032 and p = 10009 |
|||
1032 is not a quadratic residue |
|||
Find solution for n = 44402 and p = 100049 |
|||
Solution found: 30468 and 69581 |
|||
Find solution for n = 665820697 and p = 1000000009 |
|||
Solution found: 378633312 and 621366697 |
|||
Find solution for n = 881398088036 and p = 1000000000039 |
|||
Solution found: 791399408049 and 208600591990</pre> |
|||
=={{header|Go}}== |
=={{header|Go}}== |