寻找最近的较小素数的快速算法
Posted
技术标签:
【中文标题】寻找最近的较小素数的快速算法【英文标题】：Fast algorithm for finding the nearest smaller prime number 【发布时间】：2013-05-06 17:34:04 【问题描述】：例如：- 如果给定的数字是 10，我们必须返回 7（因为它是最接近的较小素数）
我能想到的方法是这样的：- Mainloop：测试给定的数字是否为素数（通过应用素数测试）， 如果它是素数，则返回数字，否则将数字减 1 并转到 Mainloop。
但是我必须在 long long int 范围上工作，这需要很多时间。
有没有更好的方法，如果我应该采用上述方式，那么我应该使用哪种素性测试？谢谢:)
【问题讨论】：
只检查奇数可以减少一半的测试。 试试Sieve of Eratosthenes。 那么您需要更快的主要测试。如果它这么慢，你可能正在使用试用部门？例如，使用 Baillie-Pomerance-Selfridge-Wagstaff 测试。 Algorithm to find largest prime number smaller than x 的可能重复项，请参阅 Find a largest prime number less than n 和 Prime number just below a number， 【参考方案1】：如果输入的大小是有界的，那么在预先计算的素数表中查找可能是最快的。
【讨论】：
“但我必须在 long long int 范围上工作” +1 查找表应该是解决方案的一部分，几乎不管其他细节是什么。 但现在我三思而后行，因为素数的分布变得越来越稀疏；如果所有输入都很大，则不需要表的下部（人口更密集）；因此可以将其修剪为更易于管理的尺寸【参考方案2】：除上述之外，还要注意Bertrand's postulate 指出在n<p<2n-2 处始终存在至少一个质数p。所以这给了你一个上限。
【讨论】：
如果您将其与素数筛法相结合，您应该找到的最高素数应低于 1.5*n-1（因为如果 n 是上限，则应低于 0.5*n+1是下界，它与 n 的距离为 0.5*n-1)。所以你应该重复初筛方法直到 floor(sqrt(1.5*n-1))。【参考方案3】：这是 Daniel Fischer 在评论中提到的 Baillie-Wagstaff 伪素性测试的伪代码实现。我们从一个简单的 Eratosthenes 筛开始，稍后我们将需要它。
function primes(n)
ps := []
sieve := makeArray(2..n, True)
for p from 2 to n step 1
if sieve(p)
ps.append(p)
for i from p * p to n step p
sieve[i] := False
return ps
powerMod 函数将基数 b 提高到指数 e，所有计算均以 m 为模；它比先求幂，然后取结果的模要快得多，因为中间计算会很大。
function powerMod(b, e, m)
x := 1
while e > 0
if e % 2 == 1
x := (b * x) % m
b := (b * b) % m
e := floor(e / 2)
return x
数论中的jacobi 函数判断a 是否为二次余数模p。
function jacobi(a, p)
a := a % p
t := 1
while a != 0
while a % 2 == 0
a := a / 2
if p % 8 == 3 or p % 8 == 5
t := -t
a, p := p , a # swap
if a % 4 == 3 and p % 4 == 3
t := -t
a := a % p
if p == 1 return t else return 0
Gary Miller 的强伪素检验基于 Pierre de Fermat 的 小定理，该定理指出如果 p 是一个素数，那么对于任何 a != 0, a ^ (p - 1) == 1 (mod p)。 Miller 的检验比 Fermat 的检验要强一些，因为它不会被 Carmichael Numbers 愚弄。
function isStrongPseudoprime(n, a)
d := n - 1; s := 0
while d % 2 == 0
d := d / 2; s := s + 1
t = powerMod(a, d, n)
if t == 1 return ProbablyPrime
while s > 0
if t == n - 1 return ProbablyPrime
t := (t * t) % n; s := s - 1
return Composite
Miller-Rabin 检验执行 k 强伪素检验，其中 k 通常介于 10 到 25 之间。强伪素检验可能会被愚弄，但如果你执行足够的他们，被愚弄的可能性很小。
function isPrime(n) # Miller-Rabin
for i from 1 to k
a := randInt(2 .. n-1)
if not isStrongPseudoprime(n, a)
return Composite
return ProbablyPrime
该素性测试对于大多数用途来说已经足够了，而且速度也足够快。但如果你想要更强大一点、更快一点的东西，可以使用基于 Lucas 链的测试。这是Lucas链的计算。
function chain(n, u, v, u2, v2, d, q, m)
k := q
while m > 0
u2 := (u2 * v2) % n; v2 := (v2 * v2 - 2 * q) % n
q := (q * q) % n
if m % 2 == 1
t1 := u2 * v; t2 := u * v2
t3 := v2 * v; t4 := u2 * u * d
u, v := t1 + t2, t3 + t4
if u % 2 == 1 u := u + n
if v % 2 == 1 v := v + n
u, v, k := (u / 2) % n, (v / 2) % n), (q * k) % n
m := floor(m / 2)
return u, v, k
由于 John Selfridge，通常使用算法初始化 Lucas 链。
function selfridge(n)
d, s := 5, 1; ds := d * s
repeat
if gcd(ds, n) > 1 return ds, 0, 0
if jacobi(ds, n) == 1 return ds, 1, (1 - ds) / 4
d, s := d + 2, s * -1; ds := d * s
然后卢卡斯伪素测试确定一个数字是素数还是可能是合数。像费马测试一样，它有两种风格，标准和强，和费马测试一样，它可以被愚弄，虽然费马测试的错误是合数可能被错误地报告为素数，但卢卡斯测试错误是素数可能被错误地报告为合数。
function isLucasPseudoprime(n) # standard
d, p, q := selfridge(n)
if p == 0 return n == d
u, v, k := chain(n, 0, 2, 1, p, d, q, (n + 1) / 2)
return u == 0
function isLucasPseudoprime(n) # strong
d, p, q := selfridge(n)
if p == 0 return n == d
s, t := 0, n + 1
while t % 2 == 0
s, t := s + 1, t / 2
u, v, k := chain(n, 1, p, 1, p, d, q, t // 2
if u == 0 or v == 0 return Prime
r := 1
while r < s
v := (v * v - 2 * k) % n; k := (K * k) % n
if v == 0 return Prime
return ProbablyComposite
然后，Baillie-Wagstaff 测试很简单。首先检查输入是否小于 2 或者是完美平方（检查平方根是否为整数）。然后用小于 100 的素数进行试除法可以快速找到大多数复合物，最后对基数 2 进行强伪素数检验（为了更加确定，有些人在基数 3 上添加了强伪素数检验），然后是卢卡斯伪素数检验做出最终确定.
function isPrime(n) # Baillie-Wagstaff
if n < 2 or isSquare(n) return False
for p in primes(100)
if n % p == 0 return n == p
return isStrongPseudoprime(n, 2) \
and isLucasPseudoprime(n) # standard or strong
Baillie-Wagstaff 测试没有已知错误。
一旦你有一个好的素数测试，你可以通过从 n 开始倒数找到小于 n 的最大素数，在第一个素数处停止。
如果你对使用素数编程感兴趣，我在我的博客中谦虚地推荐this essay，或者其他许多与素数相关的博客文章，你可以使用博客上的搜索功能找到这些文章。
【讨论】：
【参考方案4】：查看Miller-Rabin primality test。这是概率性的，但如果你这样做几百次，这几乎可以保证long long范围内的精度。
另外，如果您会使用 Java，BigInteger.isProbablePrime 可以提供帮助。 C\C++ 似乎没有用于测试素数的内置函数。
【讨论】：
【参考方案5】：看起来你正在处理this problem。
正如@Ziyao Wei所说，你可以简单地使用Miller-Rabin primality test来解决它。
这是我的解决方案
#include<cstdio>
#include<cstdlib>
#include<ctime>
short T;
unsigned long long n;
inline unsigned long long multi_mod(const unsigned long long &a,unsigned long long b,const unsigned long long &n)
unsigned long long exp(a%n),tmp(0);
while(b)
if(b&1)
tmp+=exp;
if(tmp>n)
tmp-=n;
exp<<=1;
if(exp>n)
exp-=n;
b>>=1;
return tmp;
inline unsigned long long exp_mod(unsigned long long a,unsigned long long b,const unsigned long long &c)
unsigned long long tmp(1);
while(b)
if(b&1)
tmp=multi_mod(tmp,a,c);
a=multi_mod(a,a,c);
b>>=1;
return tmp;
inline bool miller_rabbin(const unsigned long long &n,short T)
if(n==2)
return true;
if(n<2 || !(n&1))
return false;
unsigned long long a,u(n-1),x,y;
short t(0),i;
while(!(u&1))
++t;
u>>=1;
while(T--)
a=rand()%(n-1)+1;
x=exp_mod(a,u,n);
for(i=0;i<t;++i)
y=multi_mod(x,x,n);
if(y==1 && x!=1 && x!=n-1)
return false;
x=y;
if(y!=1)
return false;
return true;
int main()
srand(time(NULL));
scanf("%hd",&T);
while(T--)
for(scanf("%llu",&n);!miller_rabbin(n,20);--n);
printf("%llu\n",n);
return 0;
【讨论】：
这很糟糕。以上是关于寻找最近的较小素数的快速算法的主要内容，如果未能解决你的问题，请参考以下文章