HDu 2138 How many prime numbers 高效Miller素数測试
题目就是给出一组数。让我们測试当中有多少个是素数。
求素数有測试sqrt(n)个数的方法。有筛子方法,只是对于本题这种题目来说就都不是高效的。
本题使用Miller Rabin素数測试法。效率奇高。对于不是极其大的整数測试都差点儿是常数时间。令人神往的算法啊。
网上有个程序,好像是什么吉林的模板程序,只是我一直没看懂他是什么思路写的,是个AC的程序,只是却是错误的,呵呵,由于程序一直把9当做素数。
于是上网查找了当中原理,自己写了个程序,效率和他的几乎相同一样。通过时间基本无区别,只是我的思路是依照wiki上的写的。
wiki介绍样例:
Example[edit]
Suppose we wish to determine if n = 221 is prime. We write n − 1 = 220 as 22·55, so that we have s = 2 and d = 55. We randomly select a number a such that a < n, say a = 174. We proceed to compute:
- a20·d mod n = 17455 mod 221 = 47 ≠ 1, n − 1
- a21·d mod n = 174110 mod 221 = 220 = n − 1.
Since 220 ≡ −1 mod n, either 221 is prime, or 174 is a strong liar for 221. We try another random a, this time choosing a = 137:
- a20·d mod n = 13755 mod 221 = 188 ≠ 1, n − 1
- a21·d mod n = 137110 mod 221 = 205 ≠ n − 1.
Hence 137 is a witness for the compositeness of 221, and 174 was in fact a strong liar. Note that this tells us nothing about the factors of 221 (which are 13 and 17). However, the example with 341 in the next section shows how these calculations can sometimes produce a factor of n.
出处:http://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test
參照上面样例就能够敲代码了。挺难写的一个程序,运用好数论知识和位运算才干写出高效的程序:
代码乃属原创,转载请注明作者:
作者 靖心 http://blog.csdn.net/kenden23/article/details/30478591
#include <stdio.h> #include <stdlib.h> #include <time.h> #include <cmath> __int64 s, d; void preCalsd(__int64 N) { s = 0; N--; while ((N & 1) == 0) { s++; N >>= 1LL; } d = N; } __int64 fastPowMod(__int64 base, __int64 num, __int64 mod) { __int64 ans = 1LL; while (num) { if (num & 1) ans = ans * base % mod; base = base * base % mod; num >>= 1LL; } return ans % mod; } bool witness(__int64 a, __int64 N) { __int64 num = d, m = -1; for (int i = 0; i < s; i++) { num = d * (1LL << i); m = fastPowMod(a, num, N); if (m == 1 || m == N-1) return false; } return true; } bool miller_rabin(__int64 N) { if(N == 2) return true; if(N == 1 || ((N&1) == 0)) return false; preCalsd(N); for(int i = 0; i < 3; i++) //i < 50 { __int64 a = 2 + rand() % (N-2); if(witness(a, N)) return false; } return true; } int main() { srand(unsigned(time(NULL))); int n,cnt; __int64 a; while(scanf("%d", &n)!=EOF) { cnt = 0; while(n--) { scanf("%I64d", &a); if(miller_rabin(a)) { cnt++; //printf("%d ", a); } } printf("%d\n",cnt); } return 0; }
相关文章
- 技术干货|如何在微服务架构下构建高效的运维管理平台?
- ReadyDATA为科技研发产业基地构建私有云高效服务
- Xsolla推出全新工具,让您体验更高效的帐户管理
- 智能采购管理系统有哪些应用优势?如何高效提升医药制造业采购管理效率?
- SaaS系统平台的应用优势在哪里?如何高效提升食品制造业数智化发展水平
- B2B商城网站助力企业加快分销速度,构建高效智能的B2B网上分销平台
- 阿里云API错误中心助力开发更高效
- 基于阿里的Node全栈之路(七)[源码分享]——打造高效的开发流程
- 高效地显示Bitmap图片 1 - 有效率地加载大尺寸的位图
- 【ybtoj高效进阶 21266】历经磨难(单调队列优化DP)
- 【ybtoj高效进阶 21251】撰写博客(KMP)(单调队列优化DP)
- 【ybt高效进阶5-5-1】种植方案(状压DP)(轮廓线DP)
- 【ybt高效进阶5-1-2】【luogu P5020】货币系统
- 【ybt高效进阶2-5-3】前缀匹配