📜  原始性测试|第4组(索洛维·斯特拉森)

📅  最后修改于: 2021-04-24 22:30:27             🧑  作者: Mango

在本系列的前几篇文章中,我们已经对素数测试进行了介绍。

  • 原始性测试|第一组(介绍和学校方法)
  • 原始性测试|套装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·

log^3

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