数论
最大公因数GCD
推荐阅读
辗转相除法
略。
Stein
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
| //Stein
gcd(a,b)
{
if(a==b)
{
gcd=a;
return;
}
if(a为偶数b为奇数)//一定没有2因子
gcd=gcd(a/2,b);
if(a为奇数b为偶数)//同上
gcd=gcd(a,b/2);
if(a,b均为偶数)
gcd=2*gcd(a/2,b/2);//一定有2因子
if(a,b均为奇数)
{
if(a<b) gcd=gcd(b,a);
gcd=gcd(a-b,b);
}
}
|
斐波那契数列最大公因数的性质
gcd(fn,fm)=f(gcd(n,m))
逆元
强烈推荐:http://www.cnblogs.com/linyujun/p/5194184.html
概念
定义
对于互质的a,m∈N,称满足下式的最小x∈N为a在模m意义下的逆元(或数论倒数)。
$$
ax\equiv1\pmod{m}
$$
在代码中,由于逆元的英文为Inverse element,常常用inv(a)表示a的逆元(模数m一般定义为全局变量)
性质
集合A{a|a=1,2,3…p},B{b|b=inv(a),a∈A},则A=B。
换句话说:1~p所有数的逆元都不重复。
求法
费马小定理
费马小定理:p为质数时,$a^{p-1}\equiv1\pmod{p}$
则:$a^{p-2}\equiv inv(a)\pmod{p}$
而$a^{p-2}$可以用快速幂(边算边模)求出。
扩展欧几里得算法
不定方程ax+by=1可以用扩展欧几里得算法求出一组解,而有解的条件是gcd(a,b)=1,即a和b互质。将等式左右两边同模b,可得:
$$
ax\equiv1\pmod{b}
$$
因此求得的不定方程解x也同时是a在模b意义下的逆元。
同理y是b在模a意义下的逆元。
递推
设$x=p \% a,y=⌊\frac{p}{a}⌋$
有:
$$
\begin{aligned}
x+ya&=p\\
(x+ya)&\equiv0\pmod{p}\\
x&\equiv -ya \pmod{p}\\
x*\mathrm{inv}[a]&≡(-y) \pmod{p}\\
\mathrm{inv}[a]&≡(p-y)*\mathrm{inv}[x] \pmod{p}\\
\mathrm{inv}[a]&≡(p-⌊\frac{p}{a}⌋)*\mathrm{inv}[p \% a] \pmod{p}
\end{aligned}
$$
由于p % a< a,所以在代码中就可以这样递推求出1~p所有数的逆元:
1
2
3
4
5
6
| void initInv()
{
inv[1]=1;
for(int i=2;i<=p;++i)
inv[i]=1ll*(p-p/i)*inv(p%i)%p;
}
|
其中等号左边乘上1ll是由于在乘法计算过程中有可能会超出int范围,但是结果一般不会超出,在左边乘上1ll可以将右边的计算转换为long long 类型。
当然,在程序中全部使用long long 一定是没有问题的。
应用
求ans=(a/b) mod m
解法一:逆元 Ans=a*inv(b) mod b
解法二:若已知b|a,则ans=(a mod mb)/b
证明:
$$
\begin{aligned}
设a/b&=km+x\
a&=kmb+bx\
a % bm&=bx\
(a % bm)/b&=x\
\end{aligned}
$$
代码
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
| /*
乘法逆元模板
*/
#include<iostream>
#include<cstring>
#include<cstdlib>
#include<cstdio>
using namespace std;
template<typename T>
void write(T w)
{
if(w<0)putchar('-'),w=-w;
if(w>9)write(w/10);
putchar(w%10+48);
}
const int maxn=int(3*1e6+3);
typedef long long lint;//否则乘的时候会爆
lint p,n;
lint inv[maxn];
int main()
{
scanf("%lld%lld",&n,&p);
inv[1]=1;
puts("1");
for(lint i=2;i<=n;++i)
{
inv[i]=(p-p/i)*inv[p%i]%p;
write(inv[i]);
putchar('\n');
}
return 0;
}
|
裴蜀定理
内容
a,b是整数,且gcd(a,b)=d,∀x,y∈Z,都有d|ax+by。
特别地,一定存在整数x,y,使ax+by=d成立。
推论
gcd(a,b)=1$\iff$存在整数x,y使ax+by=1:。
推广
设a1,a2,a3……an为n个整数,d是它们的最大公约数,存在整数x1……xn使得x1*a1+x2*a2+…xn*an=d。
特别地,如果a1…an互质(不是两两互质),那么存在整数x1……xn使得x1*a1+x2*a2+…xn*an=1。
相关性质
d=gcd(a,b),则ax+by的最小正值为d。(推广后仍成立)
线性筛
朴素线性筛
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
| /*
筛出1-n当中的质数,同时计算欧拉函数、莫比乌斯函数、约数个数
注意其中i*prime[j]<=n,否则可能越界
*/
/*
prime[]数组中的素数是递增的,当i能整除prime[j],
那么i*prime[j+1]这个合数肯定被prime[j]乘以某个数筛掉。
因为i中含有prime[j],prime[j]比prime[j+1]小,即i=k*prime[j],
那么i*prime[j+1]=(k*prime[j])*prime[j+1]=k’*prime[j],
接下去的素数同理。所以不用筛下去了。
即在第一重循环的时候i*prime[j+1]会被k'筛掉
因此,在满足i%prime[j]==0这个条件之前以及第一次
满足改条件时,prime[j]必定是prime[j]*i的最小因子。
*/
flag[maxn],prime[maxn];//flag[i]表示i被筛过,i不是质数
void kick()
{
flag[1]=0;
for(int i=2;i<=n;i++)
{
if(!flag[i])
prime[++k]=i;
for(int j=1;j<=k&&i*prime[j]<=n;j++)
{
flag[i*prime[j]]=1;
if(i%prime[j]==0)break;
}
}
}
|
线性筛欧拉函数
首先请自行了解一些关于积性函数的知识。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
| phi[maxn];
void kick_phi()
{
flag[1]=0;
for(int i=2;i<=n;i++)
{
if(!flag[i])
{
prime[++k]=i;
phi[i]=i-1;
}
for(int j=1;j<=k&&i*prime[j]<=n;j++)
{
flag[i*prime[j]]=1;
if(i%prime[j]==0)
{
phi[i*prime[j]]=phi[i]*prime[j];
/*
解释1:
φ(i)=i(1-1/p1)(1-1/p2)(1-1/p3)(1-1/p4)…..(1-1/pn)
if(i%p==0) φ(i*p)=(i*p) (1-1/p1)(1-1/p2)(1-1/p3)(1-1/p4)…..(1-1/pn)
或者说,p已经是{pi}中的一个
解释2:
i*p相当于p个长度为n的区间
设x (1~i) 且gcd(x,i)=1
那么对i*p中任意x+k*i
有gcd(x+k*i,i)=gcd(i,(x+k*i)%i)=gcd(i,x)=1
所以每个区间
*/
break;
}
else phi[i*prime[j]]=phi[i]*(prime[j]-1);//质数时 积性
}
}
}
|
线性筛莫比乌斯函数
莫比乌斯函数
定义域是N
μ(1)=1
1. 当n存在平方因子时,μ(n)=0
2. 当n是素数或奇数个不同素数之积时,μ(n)=-1
3. 当n是偶数个不同素数之积时,μ(n)=1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
| void kick_mobious()
{
flag[1]=0;
for(int i=2;i<=n;i++)
{
if(!flag[i])
{
prime[++k]=i;
m[i]=-1;//2.
}
for(int j=1;j<=k&&i*prime[j]<=n;j++)
{
flag[i*prime[j]]=1;
if(i%prime[j]==0)
{
m[i*prime[j]]=0;//1.
break;
}
else m[i*prime[j]]=-m[i];//m[i*prime[j]]=m[i]*m[prime[j]]=m[i]*(-1)=-m[i] (积性函数)
}
}
}
|
线性筛约数个数
d[i]表示i的最小因数的次数
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
| void kick_facnum()
{
flag[1]=0;
for(int i=2;i<=n;i++)
{
if(!flag[i])
{
prime[++k]=i;
f[i]=2;
d[i]=1;
}
for(int j=1;j<=k&&i*prime[j]<=n;j++)
{
flag[i*prime[j]]=1;
if(i%prime[j]==0)
{
f[i*prime[j]]=f[i]/(d[i]+1)*(d[i]+2);//at the moment d[i]就是prime[j]的次数
d[i*prime[j]]=d[i]+1;
break;
}
else
{
f[i*prime[j]]=(f[i]<<1);//=f[i]*(1+1)=f[i]*2
d[i*prime[j]]=1;
}
}
}
}
|
题目