在本系列的前几篇文章中,我们已经对素数测试进行了介绍。
- 原始性测试|第一组(介绍和学校方法)
- 原始性测试|套装2(Fermat方法)
- 原始性测试|第三组(米勒-拉宾)
Solovay–Strassen素数检验是一种概率检验,用于确定数字是复数还是素数。
在深入研究代码之前,我们需要了解一些关键术语和概念,以便能够对该算法进行编码。
背景:
勒让德符号:此符号是为一对整数a和p定义的,使得p为质数。它用(a / p)表示并计算为:
= 0 if a%p = 0
(a/p) = 1 if there exists an integer k such that k2 = a(mod p)
= -1 otherwise.
欧拉证明:
(a/p) = a((p-1)/2)%p Condition (i)
雅可比符号(Jacobian Symbol) :此符号是Legendre符号的推广,其中p替换为n,其中n为
n = p1k1 * .. * pnkn
,则雅可比符号定义为:
(a/n) = ((a/p1)k1) * ((a/p2)k2) *.....* ((a/pn)kn)
如果将n用作质数,则Jacobian等于Legendre符号。这些符号具有某些属性–
1)如果gcd(a,n)!= 1,则(a / n)= 0,因此(0 / n)=0。这是因为,如果gcd(a,n)!= 1,则必须有一些素数pi这样pi将a和n相除。在这种情况下(a / pi)= 0 [根据Legendre Symbol的定义]。
2)(ab / n)=(a / n)*(b / n)。它可以很容易地从(ab / p)=(a / p)(b / p)(这里(a / p)是传奇符号)的事实得出。
3)如果a是偶数,则(a / n)=(2 / n)*((a / 2)/ n)。可以证明:
= 1 if n = 1 ( mod 8 ) or n = 7 ( mod 8 )
(2/n) = -1 if n = 3 ( mod 8 ) or n = 5 ( mod 8 )
= 0 otherwise
4) (a/n) = (n/a) * (-1)((a - 1)(n - 1) / 4) if a and n are both odd.
算法:
我们选择一个数字n来测试其素数,并选择一个随机数a,其范围为[2,n-1],并计算其雅可比行列式(a / n),如果n是质数,则雅可比行列式等于到勒让德雷(Legendre),它将满足Euler给出的条件(i)。如果它不满足给定的条件,则n为复合值,程序将停止。像其他所有概率原始性测试一样,其准确性也与迭代次数成正比。因此,我们对测试进行多次迭代以获得更准确的结果。
注意:我们对计算偶数的雅可比矩阵不感兴趣,因为我们已经知道它们除2外不是素数。
伪代码:
Algorithm for Jacobi:
Step 1 //base cases omitted
Step 2 if a>n then
Step 3 return (a mod n)/n
Step 4 else
Step 5 return (-1)((a - 1)/2)((n - 1)/2)(a/n)
Step 6 endif
Algorithm for Solovay-Strassen:
Step 1 Pick a random element a < n
Step 2 if gcd(a, n) > 1 then
Step 3 return COMPOSITE
Step 4 end if
Step 5 Compute a(n - 1)/2 using repeated squaring
and (a/n) using Jacobian algorithm.
Step 6 if (a/n) not equal to a(n - 1)/2 then
Step 7 return composite
Step 8 else
Step 9 return prime
Step 10 endif
执行:
C++
// C++ program to implement Solovay-Strassen
// Primality Test
#include
using namespace std;
// modulo function to perform binary exponentiation
long long modulo(long long base, long long exponent,
long long mod)
{
long long x = 1;
long long y = base;
while (exponent > 0)
{
if (exponent % 2 == 1)
x = (x * y) % mod;
y = (y * y) % mod;
exponent = exponent / 2;
}
return x % mod;
}
// To calculate Jacobian symbol of a given number
int calculateJacobian(long long a, long long n)
{
if (!a)
return 0;// (0/n) = 0
int ans = 1;
if (a < 0)
{
a = -a; // (a/n) = (-a/n)*(-1/n)
if (n % 4 == 3)
ans = -ans; // (-1/n) = -1 if n = 3 (mod 4)
}
if (a == 1)
return ans;// (1/n) = 1
while (a)
{
if (a < 0)
{
a = -a;// (a/n) = (-a/n)*(-1/n)
if (n % 4 == 3)
ans = -ans;// (-1/n) = -1 if n = 3 (mod 4)
}
while (a % 2 == 0)
{
a = a / 2;
if (n % 8 == 3 || n % 8 == 5)
ans = -ans;
}
swap(a, n);
if (a % 4 == 3 && n % 4 == 3)
ans = -ans;
a = a % n;
if (a > n / 2)
a = a - n;
}
if (n == 1)
return ans;
return 0;
}
// To perform the Solovay-Strassen Primality Test
bool solovoyStrassen(long long p, int iterations)
{
if (p < 2)
return false;
if (p != 2 && p % 2 == 0)
return false;
for (int i = 0; i < iterations; i++)
{
// Generate a random number a
long long a = rand() % (p - 1) + 1;
long long jacobian = (p + calculateJacobian(a, p)) % p;
long long mod = modulo(a, (p - 1) / 2, p);
if (!jacobian || mod != jacobian)
return false;
}
return true;
}
// Driver Code
int main()
{
int iterations = 50;
long long num1 = 15;
long long num2 = 13;
if (solovoyStrassen(num1, iterations))
printf("%d is prime\n",num1);
else
printf("%d is composite\n",num1);
if (solovoyStrassen(num2, iterations))
printf("%d is prime\n",num2);
else
printf("%d is composite\n",num2);
return 0;
}
Java
// Java program to implement Solovay-Strassen
// Primality Test
import java.util.Scanner;
import java.util.Random;
class GFG{
// Modulo function to perform
// binary exponentiation
static long modulo(long base,
long exponent,
long mod)
{
long x = 1;
long y = base;
while (exponent > 0)
{
if (exponent % 2 == 1)
x = (x * y) % mod;
y = (y * y) % mod;
exponent = exponent / 2;
}
return x % mod;
}
// To calculate Jacobian symbol of
// a given number
static long calculateJacobian(long a, long n)
{
if (n <= 0 || n % 2 == 0)
return 0;
long ans = 1L;
if (a < 0)
{
a = -a; // (a/n) = (-a/n)*(-1/n)
if (n % 4 == 3)
ans = -ans; // (-1/n) = -1 if n = 3 (mod 4)
}
if (a == 1)
return ans; // (1/n) = 1
while (a != 0)
{
if (a < 0)
{
a = -a; // (a/n) = (-a/n)*(-1/n)
if (n % 4 == 3)
ans = -ans; // (-1/n) = -1 if n = 3 (mod 4)
}
while (a % 2 == 0)
{
a /= 2;
if (n % 8 == 3 || n % 8 == 5)
ans = -ans;
}
long temp = a;
a = n;
n = temp;
if (a % 4 == 3 && n % 4 == 3)
ans = -ans;
a %= n;
if (a > n / 2)
a = a - n;
}
if (n == 1)
return ans;
return 0;
}
// To perform the Solovay-Strassen Primality Test
static boolean solovoyStrassen(long p,
int iteration)
{
if (p < 2)
return false;
if (p != 2 && p % 2 == 0)
return false;
// Create Object for Random Class
Random rand = new Random();
for(int i = 0; i < iteration; i++)
{
// Generate a random number r
long r = Math.abs(rand.nextLong());
long a = r % (p - 1) + 1;
long jacobian = (p + calculateJacobian(a, p)) % p;
long mod = modulo(a, (p - 1) / 2, p);
if (jacobian == 0 || mod != jacobian)
return false;
}
return true;
}
// Driver code
public static void main (String[] args)
{
int iter = 50;
long num1 = 15;
long num2 = 13;
if (solovoyStrassen(num1, iter))
System.out.println(num1 + " is prime");
else
System.out.println(num1 + " is composite");
if (solovoyStrassen(num2,iter))
System.out.println(num2 + " is prime");
else
System.out.println(num2 + " is composite");
}
}
// This code is contributed by Srishtik Dutta
Python3
# Python3 program to implement Solovay-Strassen
# Primality Test
import random
# modulo function to perform binary
# exponentiation
def modulo(base, exponent, mod):
x = 1;
y = base;
while (exponent > 0):
if (exponent % 2 == 1):
x = (x * y) % mod;
y = (y * y) % mod;
exponent = exponent // 2;
return x % mod;
# To calculate Jacobian symbol of a
# given number
def calculateJacobian(a, n):
if (a == 0):
return 0;# (0/n) = 0
ans = 1;
if (a < 0):
# (a/n) = (-a/n)*(-1/n)
a = -a;
if (n % 4 == 3):
# (-1/n) = -1 if n = 3 (mod 4)
ans = -ans;
if (a == 1):
return ans; # (1/n) = 1
while (a):
if (a < 0):
# (a/n) = (-a/n)*(-1/n)
a = -a;
if (n % 4 == 3):
# (-1/n) = -1 if n = 3 (mod 4)
ans = -ans;
while (a % 2 == 0):
a = a // 2;
if (n % 8 == 3 or n % 8 == 5):
ans = -ans;
# swap
a, n = n, a;
if (a % 4 == 3 and n % 4 == 3):
ans = -ans;
a = a % n;
if (a > n // 2):
a = a - n;
if (n == 1):
return ans;
return 0;
# To perform the Solovay- Strassen
# Primality Test
def solovoyStrassen(p, iterations):
if (p < 2):
return False;
if (p != 2 and p % 2 == 0):
return False;
for i in range(iterations):
# Generate a random number a
a = random.randrange(p - 1) + 1;
jacobian = (p + calculateJacobian(a, p)) % p;
mod = modulo(a, (p - 1) / 2, p);
if (jacobian == 0 or mod != jacobian):
return False;
return True;
# Driver Code
iterations = 50;
num1 = 15;
num2 = 13;
if (solovoyStrassen(num1, iterations)):
print(num1, "is prime ");
else:
print(num1, "is composite");
if (solovoyStrassen(num2, iterations)):
print(num2, "is prime");
else:
print(num2, "is composite");
# This code is contributed by mits
C#
/// C# program to implement Solovay-Strassen
// Primality Test
using System;
using System.Collections.Generic;
class GFG {
// Modulo function to perform
// binary exponentiation
static long modulo(long Base, long exponent, long mod)
{
long x = 1;
long y = Base;
while (exponent > 0)
{
if (exponent % 2 == 1)
x = (x * y) % mod;
y = (y * y) % mod;
exponent = exponent / 2;
}
return x % mod;
}
// To calculate Jacobian symbol of
// a given number
static long calculateJacobian(long a, long n)
{
if (n <= 0 || n % 2 == 0)
return 0;
long ans = 1L;
if (a < 0)
{
a = -a; // (a/n) = (-a/n)*(-1/n)
if (n % 4 == 3)
ans = -ans; // (-1/n) = -1 if n = 3 (mod 4)
}
if (a == 1)
return ans; // (1/n) = 1
while (a != 0)
{
if (a < 0)
{
a = -a; // (a/n) = (-a/n)*(-1/n)
if (n % 4 == 3)
ans = -ans; // (-1/n) = -1 if n = 3 (mod 4)
}
while (a % 2 == 0)
{
a /= 2;
if (n % 8 == 3 || n % 8 == 5)
ans = -ans;
}
long temp = a;
a = n;
n = temp;
if (a % 4 == 3 && n % 4 == 3)
ans = -ans;
a %= n;
if (a > n / 2)
a = a - n;
}
if (n == 1)
return ans;
return 0;
}
// To perform the Solovay-Strassen Primality Test
static bool solovoyStrassen(long p, int iteration)
{
if (p < 2)
return false;
if (p != 2 && p % 2 == 0)
return false;
// Create Object for Random Class
Random rand = new Random();
for(int i = 0; i < iteration; i++)
{
// Generate a random number r
long r = Math.Abs(rand.Next());
long a = r % (p - 1) + 1;
long jacobian = (p + calculateJacobian(a, p)) % p;
long mod = modulo(a, (p - 1) / 2, p);
if (jacobian == 0 || mod != jacobian)
return false;
}
return true;
}
// Driver code
static void Main()
{
int iter = 50;
long num1 = 15;
long num2 = 13;
if (solovoyStrassen(num1, iter))
Console.WriteLine(num1 + " is prime");
else
Console.WriteLine(num1 + " is composite");
if (solovoyStrassen(num2,iter))
Console.WriteLine(num2 + " is prime");
else
Console.WriteLine(num2 + " is composite");
}
}
// This code is contributed by divyeshrabadiya07
PHP
0)
{
if ($exponent % 2 == 1)
$x = ($x * $y) % $mod;
$y = ($y * $y) % $mod;
$exponent = $exponent / 2;
}
return $x % $mod;
}
// To calculate Jacobian
// symbol of a given number
function calculateJacobian($a, $n)
{
if (!$a)
return 0;// (0/n) = 0
$ans = 1;
if ($a < 0)
{
// (a/n) = (-a/n)*(-1/n)
$a = -$a;
if ($n % 4 == 3)
// (-1/n) = -1 if n = 3 (mod 4)
$ans = -$ans;
}
if ($a == 1)
return $ans;// (1/n) = 1
while ($a)
{
if ($a < 0)
{
// (a/n) = (-a/n)*(-1/n)
$a = -$a;
if ($n % 4 == 3)
// (-1/n) = -1 if n = 3 (mod 4)
$ans = -$ans;
}
while ($a % 2 == 0)
{
$a = $a / 2;
if ($n % 8 == 3 || $n % 8 == 5)
$ans = -$ans;
}
//swap
list($a, $n) = array($n, $a);
if ($a % 4 == 3 && $n % 4 == 3)
$ans = -$ans;
$a = $a % $n;
if ($a > $n / 2)
$a = $a - $n;
}
if ($n == 1)
return $ans;
return 0;
}
// To perform the Solovay-
// Strassen Primality Test
function solovoyStrassen($p, $iterations)
{
if ($p < 2)
return false;
if ($p != 2 && $p % 2 == 0)
return false;
for ($i = 0; $i < $iterations; $i++)
{
// Generate a random number a
$a = rand() % ($p - 1) + 1;
$jacobian = ($p +
calculateJacobian($a,
$p)) % $p;
$mod = modulo($a, ($p - 1) / 2, $p);
if (!$jacobian || $mod != $jacobian)
return false;
}
return true;
}
// Driver Code
$iterations = 50;
$num1 = 15;
$num2 = 13;
if (solovoyStrassen($num1, $iterations))
echo $num1," is prime ","\n";
else
echo $num1," is composite\n";
if (solovoyStrassen($num2, $iterations))
echo $num2," is prime\n";
else
echo $num2," is composite\n";
// This code is contributed by ajit
?>
输出 :
15 is composite
13 is prime
运行时间:使用快速算法进行模幂运算,该算法的运行时间为O(k·
n),其中k是我们测试的不同值的数量。
准确性:算法可能返回错误的答案。如果输入n确实是素数,那么输出将始终正确地可能是素数。但是,如果输入n是复合的,则输出可能错误地是素数。然后将数字n称为Euler-Jacobi伪素数。
参考:
- https://www.revolvy.com/topic/Solovay%E2%80%93Strassen%20primality%20test&item_type=topic
- https://www.wikiwand.com/cn/Solovay%E2%80%93Strassen_primality_test
- http://www.cmi.ac.in/~ramprasad/lecturenotes/comp_numb_theory/lecture2021.pdf