数论
最大公因数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 ;
}
}
}
}
题目
