一、定义
素数又称质数,是指一个大于(1)的正整数,如果除了(1)和它本身外,没有其他任何约数。偶素数只有一个为(2)。
对于正实数(x),定义(pi(x))为不大于(x)的素数个数,那么(pi(x)approxfrac{x}{ln(x)})
二、素数的判定
对于单个数或数据范围比较小时,我们采用枚举法:
bool is_prime(int x)
{
for(int i=2;i<=sqrt(x);i++)
if(!(x%i))return 0;
return 1;
}
对于比较大的数据范围,并且要求出范围内的所有质数,我们卡可以采用筛选法
1、(Eratosthenes)筛选法
void get_prime(int N)
{
memset(v,0,sizeof(v));
for(int i=2;i<=N;i++)
{
if(v[i])continue ;
prime[++cnt]=i;
for(int j=i;j<=N/i;j++)v[i*j]=1;
}
}
不过我们发现这个筛法的效率并不高,因为它会重复筛选同一个质数。因此,就有了第二种筛法——快速线性筛法,它的复杂度几乎是线性的。
void get_prime(int N)
{
memset(v,0,sizeof(v));
for(int i=2;i<=N;i++)
{
if(!v[i])prime[++cnt]=i;
for(int j=1;j<=cnt&&prime[j]*i<=N;j++)
{
v[i*prime[j]=1;
if(!(i%prime[j]))break ;
}
}
}
我们考虑每一个数至多被筛选一次,对于当前的(i),如果(i)是素数,显然它和之前的素数不重复,必定筛选出的是新的数;如果(i)是合数,由于循环只会循环到比(i)的最小的质因子还要小的数,所以新加入的质数必定小于(i)的任何质因子,感性理解一下由于质因子序列不减所以不会重复。
三、素数相关的定理
1、唯一分解定理
定理:对于一个整数(a)满足(age2),那么(a)一定可以分解为若干素数的乘积且在不计顺序时唯一。
[a=p_1^{r_1}*p_2^{r_2}*p_3^{r_3}*···*p_n^{r_n}(p_i为素数) ]
2、威尔逊定理
定理:若(p)为素数,那么((p-1)!equiv-1(mod p))成立
逆定理:若((p-1)!equiv-1(mod p))成立,那么(p)为素数
证明:1、必要性
若(p)不是质数,我们假设(a)为(p)的质因子
那么显然(amid(p-1)!),并且(a mid (p-1)!+1)
而(pmid (p-1)!+1)可得(amid (p-1)!+1)前后矛盾
2、充分性
对于(p=2)和(p=3)时,显然定理成立。
当(p ge 5)时,令
[M={2,3,4,···,p-2},N={1,2,···,p-1} ]对于(forall a in M),令
[S={a,2a,···,(p-2)a} ]我们很容易知道对于(S)中的数模(p)都不为(0),且(S)中元素模(p)两两值不相等
所以对于(forall ain M),(exist x in N)使得(axequiv1(mod p))
接下来考虑三种情况:
1、若(x=1),那么(ax\%p=a\%p=a),显然不成立
2、若(x=p-1),那么(ax\%p=[(a-1)p+p-a]\%p=p-a),显然也不成立
3、若(x=a),那么(a^2equiv1(mod p)),即((a-1)(a+1)equiv0(mod p)),那么(a=1)或(a=p-1),显然不行。
所以对于(forall a in M,exist x in M),且(x eq a),使得(axequiv1(mod p))
因此((p-2)!equiv 1(mod p))
而(p-1equiv -1(mod p))
所以((p-1)!equiv-1(mod p))
3、费马定理
定理:若(p)为质数,(a)为正整数,且(a)和(p)互质,则:(a^{p-1}equiv 1(mod p))
这个的证明就比较简单了,从剩余系分析即可。而这其实是费马小定理的一种特例。
Miller_Rabin素数测试
用费马小定理可以有多种素数的测试方法,Miller-Rabin就是其中一种:
其算法步骤如下(设(N)为询问素数):
1、计算奇数(M),使得(N=2^r*M+1)
2、选择随机数(A<N)
3、对于任意(i<r),若(A^{2^i*M} mod N=N-1),则(N)通过了随机数(A)的测试
或者(A^Mmod N=1),则(N)通过了随机数(A)的测试
经过(t)次,(N)不是素数的概率为(frac{1}{4^t})
#include<bits/stdc++.h>
using namespace std;
const int cnt=10;
int qpow(int a,int b,int mod)
{
int res=1;
for(;b;b>>=1,a=a*a%mod)
if(b&1)res=res*a%mod;
return res;
}
bool Miller_Rabin(int n)
{
if(n==2)return 1;
for(int i=0;i<cnt;i++)
{
int a=rand()%(n-2)+2;
if(qpow(a,n,n)!=a)return 0;
}
return 1;
}
int main()
{
srand(time(NULL));
int n;
scanf("%d",&n);
if(Miller_Rabin(n))printf("Probably a prime
");
else printf("A composite
");
}
这里还有第二种写法,不过需要一个推导,也是基于费马小定理
二次探测定理:
[对于质数p,方程x^2equiv1(mod p),有且只有两个根xequiv pm1(mod p) ]
所以我们考虑其逆定理,对于有非(x=1或p-1)的解的模数一定非负。
因此在用这种方法判断素数是,我们先求出(a^m),再不断对这个数进行自乘,我们设(x=a^m mod n,y=a^2 mod n),如果(y=1且x!=1或n-1)显然(p)为合数,可以排除。最后乘完后判断(y)是否为(1)即可。
bool Millar_Rabin(ll n)
{
if(n==2)return 1;
if(n<2||!(n&1))return 0;
ll m=n-1,k=0;
while((m&1)==0){k++;m>>=1;}
for(ll i=0;i<times;i++)
{
ll a=rand()%(n-1)+1;
ll x=qpow(a,m,n),y=0;
for(ll j=0;j<k;j++)
{
y=multi(x,x,n);
if(y==1&&x!=1&&x!=n-1)return 0;
x=y;
}
if(y!=1)return 0;
}
return 1;
}
欧拉定理
费马定理用来求解在素数模下,指数的同余性质,当模数为合数时,就需要用到欧拉定理。
欧拉函数(varphi):对于正整数(n),欧拉函数是小于等于(n)中与(n)互质的数的个数。
引理1:
(①)如果(n)为某一个质数(p),则:(varphi(p)=p-1)
(②)如果(n)为某一个质数的幂次(p^a),则:(varphi(p^a)=(p-1)*p^{a-1})
(③)如果(n)为任意两个互质的数(a、b)的积,则:(varphi(a*b)=varphi(a)*varphi(b))
证明:
(①)显然;
(②)因为比(p^a)小的数有(p^a-1)个。其中,被(p)整除的有(p^{a-1}-1),所以(varphi(p^a)=p^a-1-(p^{a-1}-1)=(p-1)*p^{a-1})
(③)对于小于(a*b)的数共有(a*b-1)个,其中,只有既与(a)互质,又与(b)互质的数才会和(a*b)互质,显然满足条件的个数为(varphi(a)*varphi(b))
引理2:
设(n=p_1^{r_1}*p_2^{r_2}*···*p_n^{r_n})为正整数(n)的素数幂乘积表达式,则
[varphi(n)=n*(1-frac{1}{p_1})*(1-frac{1}{p_2})*···*(1-frac{1}{p_n}) ]证明:
由于素数幂之间显然互质,由引理(1②③)可知:(varphi(n)=varphi(p_1^{r_1})*varphi(p_2^{r_2})*···*varphi(p_n^{r_n}))
(=p_1^{r_1}*p_2^{r_2}*···*p_n^{r_n}*(1-frac{1}{p_1})*(1-frac{1}{p_2})*···*(1-frac{1}{p_n}))
(=n*(1-frac{1}{p_1})*(1-frac{1}{p_2})*···*(1-frac{1}{p_n}))
欧拉定理:若(a)与(m)互质,则(a^{varphi(m)}equiv1(mod m))
Pollard Rho 算法求大数因子
在讲(Pollard Rho)之前我们先提一下另一个大整数算法(Fermat),其算法的实现是先将一个数(M)把而的因数都提取出来,使(M=N*2^k),那么显然(N)是一个奇数,对于一个奇数,如果它是质数,我们就无须分解;如果是合数,我们必定可以转化为(N=c*d)的形式,我们假设(a=frac{c+d}{2}),(b=frac{c-d}{2}(c>d)),那么显然(N=a^2-b^2)。由于不等式(a^2+b^2ge 2ab),我们可知(a=frac{c+d}{2}ge sqrt{c*d}),所以我们可以枚举大于(N)的完全平方数(a^2),判断(a^2-N)是否为完全平方数即可。一组(a、b)就能就出一组(c、d)因子。
因此我们考虑一种更加有效的方法。(Pollard Rho)算法实质上是生成了两个随机数(a、b),对于待分解的大整数(N),我们不断计算(p=gcd(|a-b|,N)),如果(p)不为(1)就又得到了一个因子,否则我们依赖前一个(a、b)计算,判断(a、b)是否形成循环,形成的话就退出。这样我们必定可以得到一个因子或得到(N)为质数,所以可以继续递归分解。
但有一点比较麻烦,一般情况下通常使用(b=(a^2+c)mod n)生成随机数据,但是这样判环会比较复杂,我们考虑一种更简便的方法,有一个(Floyd)发明的简便的算法,就是考虑如何使我们知道自己已经走完一圈了,我们只需要令(B)的速度是(A)的两倍,它们同时同向出发后,再次相遇时就走了一圈。
实际在实现时,我们发现如果(x_i=x_0)时退出,但这样太慢了,我们考虑直接对于一段(x_i)求(gcd),我们直接把(s=s*|y-x|),对于这一段的(s)求(gcd),而这个段的长度我们可以考虑倍增。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll times=10;
const ll N=5500;
ll read()
{
ll res=0,w=1;
char ch=getchar();
while(!isdigit(ch)){if(ch=='-')w=-1;ch=getchar();}
while(isdigit(ch))res=(res<<3)+(res<<1)+(ch^48),ch=getchar();
return res*w;
}
ll gcd(ll a,ll b)
{
return !b?a:gcd(b,a%b);
}
ll qpow(ll a,ll b,ll mod)
{
ll res=1;
for(;b;b>>=1,a=a*a%mod)
if(b&1)res=res*a%mod;
return res;
}
ll multi(ll a,ll b,ll mod)
{
ll res=0;
for(;b;b>>=1,a=(a<<1)%mod)
if(b&1)res=(res+a)%mod;
return res;
}
ll fac[N],r[N],cnt;
bool Millar_Rabin(ll n)
{
if(n==2)return 1;
if(n<2||!(n&1))return 0;
ll m=n-1,k=0;
while((m&1)==0){k++;m>>=1;}
for(ll i=0;i<times;i++)
{
ll a=rand()%(n-1)+1;
ll x=qpow(a,m,n),y=0;
for(ll j=0;j<k;j++)
{
y=multi(x,x,n);
if(y==1&&x!=1&&x!=n-1)return 0;
x=y;
}
if(y!=1)return 0;
}
return 1;
}
ll pollard_rho(ll n,ll c)
{
ll i=1,k=2,s=1;
ll x=rand()%(n-1)+1,y=x;
while(1)
{
i++;
x=(multi(x,x,n)+c)%n;
s=multi(s,abs(y-x),n);
if(y==x||!s)return n;
if(i==k)
{
ll g=gcd(s,n);
if(g>1&&g<n)return g;
y=x;k<<=1;
}
}
}
void find(ll n,ll c)
{
if(n==1)return ;
if(Millar_Rabin(n)){fac[++cnt]=n;return ;}
ll p=n,k=c;
while(p>=n)p=pollard_rho(p,c--);
find(p,k);
find(n/p,k);
}
int main()
{
ll n;
while(~scanf("%lld",&n))
{
cnt=0;
find(n,120);
sort(fac+1,fac+cnt+1);
r[1]=1;
ll k=2;
for(ll i=2;i<=cnt;i++)
if(fac[i]==fac[i-1])r[k-1]++;
else r[k]=1,fac[k++]=fac[i];
for(ll i=1;i<k;i++)
printf("%lld %lld
",fac[i],r[i]);
printf("
");
}
}
例1
这是一道(NOIP)的初赛题
阅读下列程序,写出程序结果
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll g(ll x)
{
if(x<=1)return x;
return (2002*g(x-1)+2003*g(x-2))%2005;
}
int main()
{
ll n;
scanf("%lld",&n);
printf("%lld",g(n));
}
我们分析一下题,实际就是给出递推式求解,我们必须推出通项公式。
因为:(g(0)=0,g(1)=1)
所以:(g(n)=[2002*g(n-1)+2003*g(n-2)]mod 2005)
(=[-3*g(n-1)-2*g(n-2)]mod 2005)
因此:(g(n)+g(n-1)=[-2*(g(n-1)+g(n-2))]mod 2005)
(equiv[(-2)^2*(g(n-2)+g(n-3))]mod 2005)
(······)
(equiv[(-2)^{n-1}*(g(0)+g(1))]mod 2005)
(equiv(-2)^{n-1}mod 2005)
又因为:
(g(n)+2*g(n-1)=[-1*(g(n-1)+2*g(n-2))]mod 2005)
(equiv[(-1)^2*(g(n-2)+2*g(n-3))]mod 2005)
(······)
(equiv[(-1)^{n-1}*(g(1)+2*g(0))]mod 2005)
((-1)^{n-1}mod 2005)
由上面两式可得:
所以我们只需要类似快速幂的做法把指数分解做即可。
例2 Visible Lattice Points
从原点看第一象限的点,给出(n)的范围,求范围内能看到的点的个数
首先我们考虑看到的点一定是对称的,所以我们只要求出下三角或上三角的点数即可。
而我们考虑每次遍历(x=n)这条直线,会加入(varphi(n))个点,所以我们要求的就是(sum_{i=1}^{n}varphi(i))
我们用欧拉函数的线性筛法即可。
#include<bits/stdc++.h>
using namespace std;
int phi[1005],prime[1005],cnt;
bool v[1005];
void pre(int N)
{
phi[1]=1;
for(int i=2;i<=N;i++)
{
if(!v[i])prime[++cnt]=i,phi[i]=i-1;
for(int j=1;j<=cnt&&i*prime[j]<=N;j++)
{
v[i*prime[j]]=1;
if(!(i%prime[j])){phi[i*prime[j]]=phi[i]*prime[j];break ;}
else phi[i*prime[j]]=phi[i]*(prime[j]-1);
}
}
}
int main()
{
pre(1000);
int t,cas=0;
scanf("%d",&t);
while(t--)
{
int n,ans=0;
scanf("%d",&n);
for(int i=1;i<=n;i++)
ans+=phi[i];
printf("%d %d %d
",++cas,n,ans*2+1);
}
}