当前位置： 代码迷 >> 综合 >> 素数一套：Miller-Rabin 素性检验算法Pollard-Rho算法线性筛——Upside down primesDivisions
  详细解决方案

素数一套：Miller-Rabin 素性检验算法Pollard-Rho算法线性筛——Upside down primesDivisions

热度：35   发布时间：2023-11-21 17:17:33.0

部分目录

    • Solved 94 / 304 K Gym 100753K Upside down primes
      • 高效判断素数
        • 快速幂取模
      • 继续 Miller-Rabin 素性检验算法
    • Unsolved 61 / 182 F Gym 100753F Divisions
      • Pollard-Rho算法
        • 朴素的算法，试除法
        • 生日悖论
        • 伪随机数序列
          • Floyd判环算法
        • 改进log n
        • 求最大因数
        • 关于template < class T> ，map和vector用法
      • AC代码
      • 不需要以上算法的方法
        • 附：判断立方的方法
      • 附：线性筛

来源题解目录
比赛link

Solved 94 / 304 K Gym 100753K Upside down primes

先粘一个我们队的ac代码：
（998ms，我直呼大爷好）
这个判断素数，是线性的，判断一个素数可简单用此，不用线性筛

#include<bits/stdc++.h>
using namespace std;
int b[20];
bool check(long long a){
    while(a){
    int tmp=a%10;if(tmp==3||tmp==4||tmp==7)return true;a/=10;}return false;
}
long long change(long long a){
     int cnt=0;long long tmp=1;long long ans=0;while(a){
    b[++cnt]=a%10;a/=10;}for(int i=cnt;i>=1;i--){
    if(b[i]==6)ans+=9*tmp;else if(b[i]==9)ans+=6*tmp;elseans+=b[i]*tmp;tmp*=10;}return ans;
}
int main(){
    long long a;cin>>a;if(a==1){
    cout<<"no";return 0;}int len=sqrt(a);//由于sqrt(2)=1.414<2,故不会进入循环，不用管 
// cout<<"len="<<len<<endl;if(check(a)){
    cout<<"no"<<endl;return 0;}for(int i=2;i<=len;i++){
    if(a%i==0){
    cout<<"no";return 0;}}a=change(a);
// cout<<"a="<<a<<endl;for(int i=2;i<=len;i++){
    if(a%i==0){
    cout<<"no";return 0;}}cout<<"yes"<<endl;return 0;
}

题意：
判断一个数以及它旋转180度后，是否都是素数
范围：n:1~10e16，还好不大，用最初级的方法没有超时是想不到的。所以说，有希望就要勇于尝试。
题目思路：
字符串读入，
几种情况：

  • n=1 ——no
  • n中有347之一 ——no
  • n不是素数 ——no
  • n倒过来不是素数 ——no
  • 以上均未no ——yes

接下来的关键是：
介绍高效判断素数的算法

高效判断素数

接下来正式介绍高效判断素数的算法
写一个框架，自顶向下码：

#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<math.h>
//#include<iostream>
//#include<algorithm>
#ifdef LOCAL
FILE*FP=freopen("text.in","r",stdin);
FILE*fp=freopen("text.out","w",stdout);
#endif
using namespace std;
#define mem(a,b) memset(a,b,sizeof(a))
#define ll long long
#define _forplus(i,a,b) for(register int i=(a); i<=(b); i++)
#define _forsub(i,a,b) for(register int i=(a); i>=(b); i--)
char ch[20]={
    0};
//判断是否素数，注意long long
bool isp(ll n){
    ...判断素数}
//得正数
inline int check(int a){
    if(a==0||a==2||a==5||a==8||a==1)return a;if(a==9)return 6;if(a==6)return 9;
}
ll change(int tag,char*ch,int len){
    ll sum=0;if(tag==1){
    _forplus(i,1,len){
    sum*=10;sum+=ch[i]-'0';}	}else if(tag==2){
    _forsub(i,len,1){
    sum*=10;sum+=check(ch[i]-'0');}}return sum;
}
//得倒数int main(){
    scanf("%s",ch+1);int len=strlen(ch+1);//判断1if(len==1&&ch[1]=='1'){
    printf("no\n");return 0;} //判断有无347 _forplus(i,1,len){
    if(ch[i]=='3'||ch[i]=='4'||ch[i]=='7'){
    printf("no\n");return 0;}}//判断正数、倒数，是素数返回1 ，不是返回0if(!isp(change(1,ch,len))||!isp(change(2,ch,len))){
    printf("no\n");return 0;} printf("yes\n");return 0;
}

接下来是判断一个素数如何高效
初版，904ms，和队长想法一样

bool isp(ll n){
    int len=(int)sqrt(n);_forplus(i,2,len){
    if(n%i==0){
    return false;}}return true;
}

卡常版，不看2的倍数，速度快一倍，451ms

bool isp(ll n){
    if(n==2)return true;if(n%2==0)return false;int len=(int)sqrt(n);for(register int i=3;i<=len;i+=2){
    if(n%i==0){
    return false;}}return true;
}

接下来是大招（甚至帮我判断好了1）：
Miller-Rabin 素性检验算法
参考知乎题解

int prime[10]={
    2, 3, 5, 7, 11, 13, 17, 19, 23, 29};
bool isp(ll n) {
    if (n <= 1) return false;if (n == 2) return true;int s = 0, t = n - 1;while (!(t % 2)) ++s, t >>= 1; // 求解 n-1=2^s*tfor (int i = 0; i < 10 && prime[i] < n; ++i) {
    int a = prime[i];int b = 1, m = a, p = t;while (p) {
     //快速幂，求 b=a^tif (p % 2) b = ((long long) b * m) % n;m = ((long long)m * m) % n;p >>= 1;}if (b == 1) continue;for (int j = 1; j <= s; ++j) {
     // 进行 s 次二次检验int k = ((long long)b * b) % n;if(k == 1 && b != n-1) return false;b = k;}if (b != 1) return false;}return true;
}

很遗憾的是，也许是我没看懂，也许听说该算法有一定的错误率，wa在了test12…
去看看cf的test12:

Test: #12, time: 0 ms., memory: 0 KB, exit code: 0, checker exit code: 1, verdict: WRONG_ANSWER
Input
19
Output
no
Answer
yes
Checker Log
wrong answer 1st words differ - expected: 'yes', found: 'no'

。。。有点无语，只是被一个小小的19卡到了啊。
看看原理再说
以下来自以上网站，讲的很好：

  1. 威尔逊定理：对于自然数 p>1，p 是质数当且仅当无语了打不出公式
  • 恶补同余
    两个整数a，b，若它们除以整数m所得的余数相等，则称a，b对于模m同余
    记作 a ≡ b (mod m) 　读作a同余于b模m，或读作a与b关于模m同余。
    比如 26 ≡ 14 (mod 12)
  • 恶补余数
    余数不能超过除数，且余数应为正整数
    所以-1除以p所得余数为p-1
    但是DEVc++上的%结果是-1，可能与简化计算有关。
    math.h库里的
double fmod(double x, double y)
返回 x 除以 y 的余数。
fmod(-1,p)

fmod(-1,p)的结果也是-1(注意返回的是double,要么(int)转化)
只对正数有效，判断若为%结果为负，"结果+除数"才是余数
但是这里并不需要算对余数，
只要利用

a%b=c
则相当于a多了个c的零头
(就算计算机算负了，也是少了个一小段）
则(a-c)%b=0
则((p-1)!+1)%p=0
又有步步取余变化,等价于：
((p-1)!%p+1)%p

写出代码：

bool isp(ll n) {
    ll res=1;for(register ll i=2;i<n;i++){
    res*=i;res%=n;}res=(res+1)%n;return res==0;
}

虽然效率更慢了（TLE于test25），单独这个并不能减少运算次数。

  1. 费马小定理：如果 p 是素数，那么有 费马小定理
  • 即（a^p-a)%p=0
  • 还是((a^(p-1)-1)* a)%p=0
    但是不是(a^(p-1)-1)%p=0，如1*2%2=0，1%2=1.

a为任意整数.
虽然费马小定理的逆命题是不成立的，但是不排除它在绝大多数情况下都是成立的。
因此用 p对若干个a符合（a^p-a)%p=0 ，判断p是素数，仍然是不准的。
如 a=2，对于经典的卡迈克数561，它虽然是合数（561=3×11×17），但是会被这个算法判定为质数。
卡迈尔数参考
此时代码

bool 你是质数吗(int n) {
    if (n <= 1) return false;int t = 1, m = 2, p = n;while(p) {
     // 快速幂取模if (p % 2) t = ((long long)t * m) % n;m = (m * m) % n;p >>= 1;}t = (t - 2) % n;return t == 0;
}

这里用到了新知识点

快速幂取模

写在：快速幂&&矩阵快速幂&&矩阵乘法&&传矩阵

继续 Miller-Rabin 素性检验算法

二次检验定理：对于质数 p，在0~p-1范围内，满足
在这里插入图片描述
的整数只有 1 和 p-1。
简洁证明的博客

  • 解释：a|b 数学 a能把b整除
    在这里插入图片描述

费马二探
如果p是一个质数，而整数a不是p的倍数，则有a^（p-1）≡1（mod p）

  • 若a，b，c为任意3个整数,m为正整数，且(m,c)=1，则当a·c≡b·c(mod m)时，有a≡b(mod m)。
    证明：a·c≡b·c(mod m)可得ac–bc≡0(mod m)可得(a-b)·c≡0(mod m)。因为(m,c)=1即m,c互质，c可以约去，a– b≡0(mod m)可得a≡b(mod m)。
    所以等价的：a ^（p-1）≡1（mod p）， a ^（p）≡a（mod p）
    参考推导
    在这里插入图片描述

我不由得想，网上写的东西都没有按人的思维去解释，看了好久，现在在练学习方法。。。也发现，动笔+举例子，提高效率。我觉得现在效率低在一直不了解上，以后多这么学
不是耗时多就好，不是笨，专注更重要
因为我想要去滑冰去，就专心看几分钟，就看完了。。。
专心

其实是这样的
a(p-1)≡1（mod p）——费马定理
条件：整数a不是p的倍数
则，a(p-1)%p应为1
若p为质数，x2≡1(mod p)，那么小于p的解只有两个，x=1或x=p?1——二次探测
条件：x2%p为1
则，各x%p应为1或p-1

网上都是代码里，从小到大，一点都不好理解
实际上是从大到小：
测试p,
选定a 为底
以下都要符合，
只要以下有一个不行，就是错:
a(p-1)除以p，余1
p-1化为u*2t（p显然是偶数）
所以上式就是
(a^ (u *2t-1))2余1除以p，
a^ (u *2t-1)应为1或p-1，
若为1或p-1，继续，直到a^u
若为别的，就不是素数
倒过来
x2%p==1&&x!=1&&x!=p-1
则不是素数
若有的
x2%p!=1
则不用看x，进入下一轮，x2当成x看

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
int prime[10]={
    2,3,5,7,11,13,17,19,23,29};
int Quick_Multiply(int a,int b,int c)  //快速积（和快速幂差不多） 
{
    long long ans=0,res=a;while(b){
    if(b&1)ans=(ans+res)%c;res=(res+res)%c;b>>=1;}return (int)ans;
}
int Quick_Power(int a,int b,int c)     //快速幂，这里就不赘述了 
{
    int ans=1,res=a;while(b){
    if(b&1)ans=Quick_Multiply(ans,res,c);res=Quick_Multiply(res,res,c);b>>=1;}return ans;
}
bool Miller_Rabin(int x)     //判断素数 
{
    int i,j,k;int s=0,t=x-1;if(x==2)  return true;   //2是素数 if(x<2||!(x&1))  return false;     //如果x是偶数或者是0,1，那它不是素数 while(!(t&1))  //将x分解成(2^s)*t的样子 {
    s++;t>>=1;}for(i=0;i<10&&prime[i]<x;++i)      //随便选一个素数进行测试 {
    int a=prime[i];int b=Quick_Power(a,t,x);      //先算出a^tfor(j=1;j<=s;++j)    //然后进行s次平方 {
    k=Quick_Multiply(b,b,x);   //求b的平方 if(k==1&&b!=1&&b!=x-1)     //用二次探测判断 return false;b=k;}if(b!=1)  return false;   //用费马小定律判断 }return true;   //如果进行多次测试都是对的，那么x就很有可能是素数 
}
int main()
{
    int x;scanf("%d",&x);if(Miller_Rabin(x))  printf("Yes");else  printf("No");return 0;
}

终于看完了
真的挺怀疑自己的，好笨啊
呜呜呜
花了10个小时，搞一个，寒假也是
比常课难多了，为什么这样子。。。
哎，，，
成果：

#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<math.h>
//#include<iostream>
//#include<algorithm>
#ifdef LOCAL
FILE*FP=freopen("text.in","r",stdin);
FILE*fp=freopen("text.out","w",stdout);
#endif
using namespace std;
#define mem(a,b) memset(a,b,sizeof(a))
#define ll long long
#define _forplus(i,a,b) for(register int i=(a); i<=(b); i++)
#define _forsub(i,a,b) for(register int i=(a); i>=(b); i--)
#define N 7
char ch[20]={
    0};
//判断是否素数，注意long long
int prime[N]={
    2,325,9375,28178,450775,9780504,1795265022};
ll Quick_Multiply(ll a,ll b,ll c)  //快速积（和快速幂差不多） 
{
    ll ans=0,res=a;while(b){
    if(b&1)ans=(ans+res)%c;res=(res+res)%c;b>>=1;}return ans;
}
ll Quick_Power(ll a,ll b,ll c)     //快速幂，这里就不赘述了 
{
    ll ans=1,res=a;while(b){
    if(b&1)ans=Quick_Multiply(ans,res,c);res=Quick_Multiply(res,res,c);b>>=1;}return ans;
}
bool isp(ll x)     //Miller_Rabin判断素数 
{
    int i,j;ll k;//注意溢出！（宁愿全ll） ll s=0,t=x-1;if(x==2)return true;if(x<2||!(x&1))return false;while(!(t&1)){
    s++;t>>=1;}for(i=0;i<N&&prime[i]<x;++i){
    ll a=prime[i];ll b=Quick_Power(a,t,x);for(j=1;j<=s;++j){
    k=Quick_Multiply(b,b,x);if(k==1&&b!=1&&b!=x-1)return false;b=k;}if(b!=1)return false;}return true;
}
/*bool isp(ll n){if(n==2)return true;if(n%2==0)return false;int len=(int)sqrt(n);for(register int i=3;i<=len;i+=2){if(n%i==0){return false;}}return true; }*/
//得正数
inline int check(int a){
    if(a==0||a==2||a==5||a==8||a==1)return a;if(a==9)return 6;if(a==6)return 9;
}
ll change(int tag,char*ch,int len){
    ll sum=0;if(tag==1){
    _forplus(i,1,len){
    sum*=10;sum+=ch[i]-'0';}	}else if(tag==2){
    _forsub(i,len,1){
    sum*=10;sum+=check(ch[i]-'0');}}return sum;
}
//得倒数int main(){
    scanf("%s",ch+1);int len=strlen(ch+1);//判断1if(len==1&&ch[1]=='1'){
    printf("no\n");return 0;} //判断有无347 _forplus(i,1,len){
    if(ch[i]=='3'||ch[i]=='4'||ch[i]=='7'){
    printf("no\n");return 0;}}//判断正数、倒数，是素数返回1 ，不是返回0if(!isp(change(1,ch,len))||!isp(change(2,ch,len))){
    printf("no\n");return 0;} printf("yes\n");return 0;
}

Accepted 30 ms

2021/6/15补记：判断素数的是x，而prime数组中只是底数而已，不是素数

Unsolved 61 / 182 F Gym 100753F Divisions

这是一道分解质因数的题，和上题关系比较密切

给一个数，问它能被多少个数整除？
范围：1~1e18

和上题不太一样，不能快速判断素数了，因为要求有多少个。

那我们有一个定理，只要知道可以分解出的各质数为多少个
如
12=2231
则有
1=2030 3=2031
2=2130 6=2131
4=2230 12=2231
即
n=2q13q2…
答案为
累乘(qi+1) , i from 1 to n.
所以这就是算质因数个数方法
但是打素数表，一个一个判断就超时了

先看看例子：

INPUT OUTPUT
Sample Input 1 Sample Output 1
12 6
Sample Input 2 Sample Output 2
999999999999999989 2
Sample Input 3 Sample Output 3
100000007700000049 4

学会了快速判断素数的算法后，我似乎有点思路了
有思路就多想一下

对于输入n
如999999999999999989，可以快速判断为素数
则输出2
如100000007700000049，不是素数
则寻找子素数1000000072
先除一个100000007，再看是不是素数？
但是，还是慢了
如果枚举试除的话，到100000007是超时的

翻题解，不出所料，思路是对的
而需要学一个新算法：
Pollard-Rho算法
——John Pollard发明的， 能快速找到大整数的一个非1、非自身的因子的算法。

Pollard-Rho算法

Pollard_rho算法的大致流程是 ：
先判断当前数是否是素数（Miller_rabin）了，如果是则直接返回。如果不是素数的话，试图找到当前数的一个因子（可以不是质因子）。然后递归对该因子和约去这个因子的另一个因子进行分解
全B站唯一讲解视频——看博客总是呆懵，看看视频试试是不是更合适
确实理解的多一点
知乎的，他有很多很多算法知识
这个文章关键是要一步一步看，不要跳太快了，一句一句慢慢看，循序渐进

目的：快速找到大整数的一个非1、非自身的因子的算法。

朴素的算法，试除法

方法1：
从2到n-1一个一个试
改进方法：
如果是素数，退出，之后就是2~sqrt(n)
复杂度：O(n)
方法2：
搞笑的方法，随机数看是不是

x = randint(2, n - 1); // 生成2和n-1之间的随机数

不是就再来一次
然而它正是Pollard-Rho算法的基础
O(n)
方法3：

x = randint(2, n - 1); // 生成2和n-1之间的随机数
d = gcd(n, x);       // 求gcd

d不为0，返回d
个人改进猜想：
（其实也可以只产生sqrt(n)内的随机数，以后应该也可以只产生sqrt(n)内的随机数而不用GCD）
（后来发现还是要GCD的，GCD不是只把数降到sqrt（n）)

生日悖论

如果我们不断在某个范围内生成随机整数，很快便会生成到重复的数，期望大约在根号级别

应用到原题上，即是对于最坏情形 n=p2，如果我们不断在 [1,n-1]间生成随机数，那么期望在生成大约sqrt( p )=n1/4个数后，可以出现两个在模p下相同的数（注意[0,n-1]间的随机数模p大致是[0,p-1]间的随机数）。那么这两个数的差的绝对值 abs，就一定满足abs整除p，则gcd(abs,n)是p的倍数，是因子之一

但这件事意义并没有那么大，因为这个高概率是在两两比较下才成立的。所以我们需要一些技巧。

伪随机数序列

Pollard使用一种特别的伪随机数生成器来生成 [0,n-1]间的伪随机数序列：设序列第一个数为 x ， f(x)=(x2+c) mod n，则 x,f(x),f(f(x)),…为一个伪随机数序列。

（也解释了随机数生成原理）

因为每个数都是由前一个数决定的，可以生成的数又是有限的，那么迟早会进入循环。当然，这个循环很可能是混循环

Floyd判环算法

我们在这里使用Floyd判环算法
（跟最短路的Floyd算法是同一个人，但不是同一个算法，这个也叫龟兔赛跑算法）
附：这个UP主的最短路算法学习，这个UP主讲的很好
设置两个变量 t，r，每次判断是否有gcd ( | t - r | , n ) > 1，如果没有，就令 t=f ( t ) , r = f ( f ( r ) ) 。因为 r 跑得更快，如果没有找到答案，最终会与 t 在环上相遇,，t==r，这时退出，换一个 c 重新生成伪随机数。

那么，这有什么好处呢？

注意到若 | t - r | 整除 p , 则 | f ( t ) - f ( r ) | 整除 p ，

由此可得，只要环上距离为 d 的两个数满足条件，那么所有距离为 d 的数都满足条件。在Floyd判环的过程中，每次移动都相当于在检查一个新的距离 d ，这样就不需要进行两两比较了
（这里的 d 是指环上的距离，以相邻为1个单位 )

这个算法的复杂度依赖于这个伪随机数生成器的随机程度，还没有被严格证明。如果它是足够随机的，那么期望复杂度显然是 O ( n1/4log n)。
（我也不知道它是怎么显然的，，，害，可能是：
期望在生成大约sqrt( p )=n1/4个数后，可以出现两个在模p下相同的数）

当前已经比较实用的代码：

ll Pollard_Rho(ll N)
{
    if (N == 4) // 特判4，4找不到解，因为不够随机return 2;if (is_prime(N)) // 特判质数return N;while (1){
    ll c = randint(1, N - 1); // 生成随机的cauto f = [=](ll x) {
     return ((lll)x * x + c) % N; }; // lll表示__int128，防溢出ll t = f(0), r = f(f(0));while (t != r){
    ll d = gcd(abs(t - r), N);if (d > 1)return d;t = f(t), r = f(f(r));}}
}

特判4，4找不到解，因为不够随机
解释：

auto f = [=](ll x) {
     return ((lll)x * x + c) % N; };

这相当于一个函数，

auto f(ll x){
    return ((lll)x * x + c) % N;
}

= [=] 我现在还看不懂，，，还没学c++
区别还有最后的 ；就是了。
2021/6/15补记：参考网站

改进log n

原理：
如果 gcd ( d , n ) > 1 ,
则 gcd ( kd mod n , n ) > 1,
感觉是对的，不知道怎么证明
所以我们可以减少求公因数的次数，即先把一些待选数乘起来，再统一与 n 求公因数

令固定距离C=128

ll Pollard_Rho(ll N)
{
    if (N == 4)return 2;if (is_prime(N))return N;while (1){
    ll c = randint(1, N - 1);auto f = [=](ll x) {
     return ((lll)x * x + c) % N; };ll t = 0, r = 0, p = 1, q;do{
    for (int i = 0; i < 128; ++i) // 令固定距离C=128{
    t = f(t), r = f(f(r));if (t == r || (q = (lll)p * abs(t - r) % N) == 0) // 如果发现环，或者积即将为0，退出break;p = q;}ll d = gcd(p, N);if (d > 1)return d;} while (t != r);}
}

注意，如果积在乘的过程中等于0，那么后面都是0了，可以直接退出，而且这时大概率已经找到了解（除非 | t - r | 整除 n )
（比如，4*3%6==0，则gcd( 4 ,6 )就是解）

求最大因数

给了：

ll max_prime_factor(ll x)
{
    ll fac = Pollard_Rho(x);if (fac == x)return x;elsereturn max(max_prime_factor(fac), max_prime_factor(x / fac));
}

真的这样吗？
比如12
答案应为6
若为
第一层，fac=4
return max(max_prime_factor(4), max_prime_factor(3));
哦，这因子是指
最大素因子
答案为3

关于template < class T> ，map和vector用法

这里

AC代码

特判一下1，注意

#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<math.h>
#include<time.h>
#include<iostream>
#include<algorithm>
#include<random>
#include<map>
#ifdef LOCAL
FILE*FP=freopen("text.in","r",stdin);
FILE*fp=freopen("text.out","w",stdout);
#endif
using namespace std;
#define mem(a,b) memset(a,b,sizeof(a))
#define ll long long
#define lll __int128
#define _forplus(i,a,b) for(register int i=(a); i<=(b); i++)
#define _forsub(i,a,b) for(register int i=(a); i>=(b); i--)
#define INF (1e9+7)
#define N 7
map<int,int>p;
//判断是否素数，注意long long
ll prime[N]={
    2,325,9375,28178,450775,9780504,1795265022};
ll Quick_Multiply(ll a,ll b,ll c)  //快速积（和快速幂差不多） 
{
    ll ans=0,res=a;while(b){
    if(b&1)ans=(ans+res)%c;res=(res+res)%c;b>>=1;}return ans;
}
ll Quick_Power(ll a,ll b,ll c)     //快速幂，这里就不赘述了 
{
    ll ans=1,res=a;while(b){
    if(b&1)ans=Quick_Multiply(ans,res,c);res=Quick_Multiply(res,res,c);b>>=1;}return ans;
}
bool is_prime(ll x)     //Miller_Rabin判断素数 
{
    int i,j;ll k;//注意溢出！（宁愿全ll） ll s=0,t=x-1;if(x==2)return true;if(x<2||!(x&1))return false;while(!(t&1)){
    s++;t>>=1;}for(i=0;i<N&&prime[i]<x;++i){
    ll a=prime[i];ll b=Quick_Power(a,t,x);for(j=1;j<=s;++j){
    k=Quick_Multiply(b,b,x);if(k==1&&b!=1&&b!=x-1)return false;b=k;}if(b!=1)return false;}return true;
}
template <class T>
T randint(T l, T r = 0) // 生成随机数建议用<random>里的引擎和分布，而不是rand()模数，那样保证是均匀分布
{
    static mt19937 eng(time(0));if (l > r)swap(l, r);uniform_int_distribution<T> dis(l, r);return dis(eng);
}
ll Pollard_Rho(ll n)
{
    if (n == 4) // 特判4，4找不到解，因为不够随机return 2;if (n == 1)//特判1 ，但是小心1的因子只有1个，其他质数是2个，它不是质数 return 1;if (is_prime(n)) // 特判质数return n;while (1){
    ll c = randint((ll)1, n - 1); // 生成随机的cauto f = [=](ll x) {
     return ((lll)x * x + c) % n; }; // lll表示__int128，防溢出ll t = f(0), r = f(f(0));while (t != r){
    ll d = __gcd(abs(t - r), n);if (d > 1)return d;t = f(t), r = f(f(r));}}
}
void findp(ll n){
    ll k;k=Pollard_Rho(n);if(k==n){
    p[k]++;return;}findp(k);findp(n/k);
}
int main(){
    ll n;scanf("%lld",&n);if(n==1){
    printf("1\n");return 0;} findp(n);ll s=1;map<int,int>::iterator i;map<int,int>::iterator end=p.end();for(i=p.begin();i!=end;i++){
    s*=i->second+1;}printf("%d\n",s);return 0;
}

不需要以上算法的方法

道高一尺魔高一丈，我找到了不需要以上算法的方法
不得不说，这个方法，用了数论里的“分段”，
乘法是可以分段的，1e18可以化为1e6内分解+剩余几种情况讨论
根本不用直接找因子
并且LL也是1e18，这个方法降到1e6+几个判定，一般是足够的
妙了妙了
可以用于看各素数因子各有多少个，剩余情况Pollard-Rho找也可以

大意是，先把1e6内的素数找到，并分解之，
然后如果剩下的还有，那必由2个以内素数因子组成
因为如果有第三个数，必大于1e6，则整体超过1e18，不对
只要讨论，剩下的是
1，1
2，一个素数
3，因子是2个相同质数
4，因子是2个不同质数

#include <iostream>
#include <algorithm>
#include <string>
#include <cmath>
#include <vector>
#include <sstream>
#include <stack>
#include <queue>
#include <map>
#include <climits>
#include <cstdio>
#include <set>
#include <bitset>
using namespace std;
#define db(a) (cout << (#a) << " = " << (a) << endl)
typedef unsigned long long ll;/** Deterministic Miller-Rabin **/#define MAX_P (1000000ull)
vector<ll> primes;
vector<ll> A({
    2, 3, 5, 7, 11, 13, 17, 19, 23}); // if n < 3,825,123,056,546,413,051, it is enough to test a = 2, 3, 5, 7, 11, 13, 17, 19, and 23.
ll N;ll mulmod(ll a, ll b, ll n)
{
    ll erg = 0;ll r = 0;while(b > 0){
    // unsigned long long gives enough room for base 10 operationsll x = ((a%n) * (b%10)) % n;        for(ll i=0;i<r;i++) x = (x*10)%n;erg = (erg + x) % n;r++;b /= 10;}return erg;
}ll fastexp(ll a, ll b, ll n)
{
    if(b == 0) return 1;if(b == 1) return a%n;ll res = 1;while(b > 0){
    if(b % 2 == 1) res = mulmod(a, res, n);a = mulmod(a, a, n);b /= 2;}return res;
}bool mrtest(ll n)
{
    if(n == 1) return false;ll d = n-1;ll s = 0;while(d % 2 == 0){
    s++;d /= 2;}for(ll j=0;j<(ll)A.size();j++){
    if(A[j] > n-1) continue;ll ad = fastexp(A[j], d, n);if(ad % n == 1) continue;bool notcomp = false;for(ll r=0;r<=max(0ull,s-1);r++){
    ll rr = fastexp(2,r,n);ll ard = fastexp(ad, rr, n);			if(ard % n == n-1) {
    notcomp = true; break;}}if(!notcomp){
    return false;}}return true;
}bool sieve[MAX_P+1];void preprimes(ll l)
{
    ll limit = min(MAX_P, l);for(ll i=0;i<=limit;i++) sieve[i] = false;for(ll p=2;p<=limit;p++){
    if(sieve[p]) continue;sieve[p] = true;primes.push_back(p);ll i = 2*p;while(i <= limit){
    sieve[i] = true;i += p;}}
}// how accurate is this? can we generate cases were this fails?
bool sqtest(ll n)
{
    ll ns = (ll) sqrt((long double)n);if(ns*ns==n) return true;cout<<ns*ns;ns++;if(ns*ns==n) return true;cout<<ns*ns;return false;
}int main()
{
      	cin>>N;ll erg = 1;preprimes(N);ll restN = N;for(ll i=0;i<(ll)primes.size();i++){
    	if(N % primes[i] == 0) {
    			ll tN = N;ll p = 0;	while(tN % primes[i] == 0) {
    tN /= primes[i]; restN /= primes[i]; p++;}// p occurences of primes[i]erg *= p+1;				}}//如果剩下的数存在 （必大于MAX_P） if(restN> 1 && sqtest(restN)) erg *= 3;		// 如果这个数能整开平方 ,则该素数项指数为2 else if(restN > 1 && !mrtest(restN)) erg *= 4;// 如果这个数是合数 ,则该素数项指数为1，1 //为什么只分解成了2个素数 ，//因为如果有第三个数，必大于1e6，则整体超过1e18，不对 else if(restN > 1) erg *= 2;//如果它就是质数，则该项指数1 cout << erg << "\n";return 0;
}

6/15补记：这里素数判断方法大同小异，参考这里

附：判断立方的方法

可以通过i=(ll)pow(x, 1.0/3)求得x的立方根
如同这样如法炮制
pow(i,3)
i++
pow(i,3)
看看其中一个是不是
因为浮点数误差，可能本来是，
比如n，算成n-eps，ll成n-1了
这样子就可以防止之

附：线性筛

埃氏筛可以做到这样：

for(i=2;i<=sqrt(n);i++){
    if(p[i]==0){
    for(j=i*i;j<=n;j+=i){
    //之所以可以i*i开始，是因为i*(i更小的)即(i更小的)*i，即已经算过了！p[j]=1;}}
}

线性筛就完全不会重复：
参考博客

> 算法思路：
> 
>  
> 对于每一个数（无论质数合数）x，筛掉所有不大于x最小质因子的质数乘以x的数。比如对于77,它分解质因数是7*11，那么筛掉所有小于7的质数*77，筛掉2*77、3*77、5*77。
> 
>   好吧，是不是听起来太简单了。。。。没事，重点在证明。
> 
> 算法证明：
> 
>   首先我们要明确证明思路。如果要证明它是对的，只要保证两点：没用重复筛、没有漏筛
> 
>   1、没有重复筛。
> 
>   
> 我们假设一个合数分解成p1*p2*p3，并且保证p1<p2<p3。
> 我们知道，筛掉这个合数的机会有：p1和p2*p3，p2和p1*p3，p3和p1*p2。
> 但我们知道，我们选择的那个质数必须小于那个合数的最小质因子。
> 比如p2和p1*p3，因为p2>p1,所以这样是筛不到的。
> 唯一会筛的是第一种：p1和p2*p3。
> 
>   2、没有漏筛。
> 
>  
> 还是假设把这个合数分解质因数a*b*c，保证a<b<c
> 然后我们设s=b*c，s肯定小于刚才那个合数，说明肯定对于它已经筛过，
> 然后a肯定小于s，因为s=b*c，并且a是最小的因子。
> 说明a*s也就是这个合数一定筛过。
> 
> 证明没看懂的直接看代码吧。。挺好背的。

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cstdlib>
#define _forplus(i,a,b) for( register int i=(a); i<=(b); i++)
#define fastio std::ios::sync_with_stdio(false);std::cin.tie(0);
using namespace std;
#define N 10000100
int prime[N],vis[N]={
    0};
int cnt=0,n=10000000;
void pre(){
    //计算1~n的素数_forplus(i,2,n){
    if(!vis[i])  prime[++cnt]=i;//如果没有筛过，记录素数_forplus(j,1,cnt){
    //其中记录数组里的素数保证严格递增if(i*prime[j]>n)  break;//保证小于n，要不然没有意义vis[i*prime[j]]=1;//筛去这个合数if(!(i%prime[j]))  break;//如果>=这个数的最小质因子，那就结束}}
}
int main(){
    fastio//cin>>n;pre();_forplus(i,1,cnt){
    cout<<prime[i]<<' ';}
}
  相关解决方案

代码迷 - 您身边的软件异常，代码异常，编程异常助手（发现，解决，分享）
Copyright © 2009-2013 DAIMAMI.COM. All rights reserved.