萌新学数论啦（一）

最新推荐文章于 2020-03-25 21:42:18 发布

原创最新推荐文章于 2020-03-25 21:42:18 发布 · 338 阅读

0 ·

CC 4.0 BY-SA版权

数论专栏收录该内容

11 篇文章

订阅专栏

一.基本概念：

整除性：

一个整数a能够被另一个整数d整除，记作 d|a，意味着对于某个整数k，有a = k*d。

如果d|a，且d >= 0，则d为a的一个约数。

对于任意整数x, y，有 d|a 且 d|b 则 d|(ax+by)

每个整数a都可以被其平凡约束1和a整除，a的非平凡约数也称为a的因子。

公约数；

两个不同时为0的整数a和b的最大公约数表示成gcd(a, b)，有gcd(0,0) = 0, gcd(a,0) = |a|;

定理：∀a,b∈N, a,b ≠ 0，则gcd(a,b)是a与b的线性组合集合{ax+by: x,y∈Z}中的最小正元素

如果gcd(a,b) = 1，则a与b称为互质数，有：gcd(a,p) = gcd(b,p) = 1，则gcd(ab,p) = 1。

唯一因子的分解：

定理：对于所有素数p和所有整数a和b，如果p|ab，则p|a 或 p|b（或两者都成立）。

定理：合数a仅能以一种方式，写成如下乘积的形式：

其中pi为素数，p1 < p2 < ... < pr, 且ei为正整数

二.最大公约数:

gcd递归定理：对于任意非负整数a和正整数b，gcd(a,b) = gcd(b, a%b);

Bezout定理：对于任何整数a和b和它们的的最大公约数d，关于未知数x和y的线性丢番图方程：ax+by = m,

有整数解时当且仅当m是d的倍数。

代码模板：

int gcd(int a, int b)
{
	if(b == 0)
		return a;
	else
		return gcd(b, a%b);
}

int exgcd(int a, int b, int &x, int &y)
{
	if(b == 0){
		x = 1, y = 0;
		return a;
	}
	int d = exgcd(b, a%b, y, x);
	y -= a/b * x;
	return d;
}

最小公倍数：lcm(a,b)，有 ab = gcd(a,b) * lcm(a,b); lcm(S/a, S/b) = S/gcd(a, b)

三.素数筛，快速幂：

int prime[N], vis[N], cnt;	//cnt表示[1,N]中素数个数 
void init(){
    memset(vis, 0, sizeof(vis))
    for(int i = 2; i <= N; i ++){
    	if(!vis[i])
    		prime[cnt++] = i;     
        for(int j = i*i; j < N; j += i){
            vis[j] = 1;
        }
    }
}

该算法的时间复杂度为O(nloglogn)

LL pow_mod(LL a, LL b, LL p)
{ 
    LL ret = 1;
    while(b){
        if(b & 1) 
			ret = (ret * a) % p;
        a = (a * a) % p;
        b >>= 1;
    }
    return ret;
}

四.欧拉函数：

欧拉定理:若n,a为正整数且互质，则a^φ(n) ≡ 1 (mod n)，φ(n)为欧拉函数，欧拉函数是求在1到n-1中与n互质的数的数目。

欧拉函数是积性函数，即若m,n互质，有f(mn) = f(m)*f(n).(所有的幂函数都是积性函数)

//求某个数的欧拉函数值
int getphi(int n)
{
	int m = sqrt(n+0.5);
	int ans = n;
	for(int i = 2; i <= m; i++){
		if(n % i == 0){
			ans = ans/i * (i-1);
			while(n%i == 0){
				n /= i;
			}
		}
	}
	if(n > 1)
		ans = ans/n * (n-1);
	return ans;
} 

//计算1~n的所有数字的欧拉函数值，O(n)
const int N = 1e6+10 ;
int phi[N], prime[N];
int tot;//tot计数，表示prime[N]中有多少质数 
void Euler(){
    phi[1] = 1;
    for(int i = 2; i < N; i ++){
        if(!phi[i]){
            phi[i] = i-1;
            prime[tot ++] = i;
        }
        for(int j = 0; j < tot && 1ll*i*prime[j] < N; j ++){
            if(i % prime[j]) phi[i * prime[j]] = phi[i] * (prime[j]-1);
            else{
                phi[i * prime[j] ] = phi[i] * prime[j];
                break;
            }
        }
    }
}

一个类似埃筛的写法，速率没有上面的快（但是好写啊~）

euler[1]=1;    
for(int i = 2;i < Max;i++)    
    euler[i]=i;    
for(int i = 2;i < Max;i++)    
    if(euler[i] == i)    
        for(int j = i;j < Max; j+=i)    
            euler[j] = euler[j]/i*(i-1);//先进行除法是为了防止中间数据的溢出

最近试了一下，发现线筛竟然快了一点？，在2~1000005数据范围下打表，第一种237.7s, 线筛225.8s。。。。

顺带补充一个威尔逊定理：若p为质数，则p能被(p-1)!+1整除。

五.逆元：

对于a*x = 1 (mod p)，我们把x看成a的(数论)倒数，称x为a关于p的逆元。（PS:a和p互质）

a的逆元，我们用inv(a)来表示，那么(a / b) % p = (a * inv(b) ) % p = (a % p * inv(b) % p) % p。

逆元的求法：

①用快速幂：已知a^(p-1) ≡ 1 (mod p)，两边同时除以a，所以a^(p-2) ≡ inv(a) (mod p)，所以inv(a) = a^(p-2) (mod p)。

pow_mod(a, p-2, p);

②用扩展欧几里得算法：已知ax + by = 1，如果a和b互质，则有解，并且这个解的x就是a关于b的逆元，y是b关于a的逆元。

void ex_gcd(LL a, LL b, LL &x, LL &y, LL &d){
    if (!b) {d = a, x = 1, y = 0;}
    else{
        ex_gcd(b, a % b, y, x, d);
        y -= x * (a / b);
    }
}
LL inv(LL t, LL p){//如果不存在，返回-1 
    LL d, x, y;
    ex_gcd(t, p, x, y, d);
    return d == 1 ? (x % p + p) % p : -1;
}

③.当p是个质数的时候有inv(a) = (p - p / a) * inv(p % a) % p. 下面为逆元线性筛的写法

const int mod = 1000000009;  
const int maxn = 10005;  
int inv[maxn];  
inv[1] = 1;  
for(int i = 2; i < 10000; i++)  
    inv[i] = inv[mod % i] * (mod - mod / i) % mod;

附：求阶乘的逆元

inv[maxn] = mod_pow(fac[maxn], mod-2);  
for(ll i = maxn-1;i >= 0; i--)  
    inv[i] = (inv[i+1] * (i+1)) % mod;

六.中国剩余定理：

我们先来看一下这个东西：线性同余方程.对于ax ≡ b mod n来说，可以将其转化为一个线性不定方程，即ax-ny = b.

定理：假设方程ax ≡ b mod n有解(即有d|b，d = gcd(a,b))，x0是该方程的任意一个解，则该方程对模n恰有d个不同的解.

分别为 xi = x0 + i*(n/d) (i = 1,2,..,d-1).

代码实现如下：

int ext_euclid(int a,int b,int &x,int &y)  //求gcd(a,b)=ax+by
{
    int t, d;
    if (b == 0){
		x = 1;
		y = 0;
		return a;
	}
    d = ext_euclid(b,a % b, x, y);
    t = x;
    x = y;
    y = t-a/b*y;
    return d;
}

void modular_equation(int a, int b, int n, int s[])
{
    int e,i,d;
    int x,y;
    d = ext_euclid(a,n,x,y);
    if (b%d>0)
       printf("No answer!\n");
    else{
    	e=(x*(b/d))%n;
        for(i = 0; i < d; i++){
            s[i] = (e + i* (n/d) ) %n; 
        }
}

如果需要对给定的一个单变量模线性方程组（如下式）进行求解，就要利用到中国剩余定理.

中国剩余定理：假设整数m1,m2,.....,m_{n 两两互质，则对任意的整数：a1,a2,...an,上面的方程组有解。（这是由孙子搞出来的）}

_{可通过如下方式得到：}

_{M = m1 * m2 * ... *}m_n,Mi = M / m_i

t_iM_i ≡ 1 mod m_{i，通解形式为 x = a1t1M1 + a2t2M2 + ... + antnMn+KM, K∈Z}

_{在模M的意义下，方程组只有一个解。}

//n个方程：x=a[i](mod m[i]) (0<=i<n)
LL china(int n, LL *a, LL *m){
    LL M = 1, ret = 0;
    for(int i = 0; i < n; i ++) M *= m[i];
    for(int i = 0; i < n; i ++){
        LL w = M / m[i];
        ret = (ret + w * inv(w, m[i]) * a[i]) % M;
    }
    return (ret + M) % M;
}

如果m1,m2,.....,m_{n两两不互质呢？，请看如下代码：}

typedef pair<LL, LL> PLL;
PLL linear(LL A[], LL B[], LL M[], int n) {	//求解A[i]x = B[i] (mod M[i]),总共n个线性方程组 
    LL x = 0, m = 1;
    for(int i = 0; i < n; i ++) {
        LL a = A[i] * m, b = B[i] - A[i]*x, d = gcd(M[i], a);
        if(b % d != 0)  
			return PLL(0, -1);	//答案不存在，返回-1 
        LL t = b/d * inv(a/d, M[i]/d)%(M[i]/d);
        x = x + m*t;
        m *= M[i]/d;
    }
    x = (x % m + m ) % m;
    return 
		PLL(x, m);//返回的x就是答案，m是最后的lcm值 
}

七.大组合数-卢卡斯定理：

我们不妨先来看看一般的组合数怎么写. (P可以忽略~)

有一种是比较Low的O(n²)求法，是利用的杨辉三角的性质来写，这里就不贴代码了。

另外一种是利用逆元这个东西，这是一个O(n)的求法，这个也就先不写了。

于是，就到了我们的卢卡斯定理：C(n, m) % p = C(n / p, m / p) * C(n%p, m%p) % p

证明：其实整个博客的定理我基本上都不会证明只有有智慧的人才可以看到...

LL PowMod(LL a, LL b, LL MOD){    
    LL ret=1;    
    while(b){    
        if(b&1) 
			ret = (ret * a) % MOD;    
        a = (a*a) % MOD;    
        b >>= 1;    
    }    
    return ret;    
}    
LL fac[100005];    
LL Get_Fact(LL p)
{    
    fac[0] = 1;    
    for(LL i = 1;i <= p; i++)    
        fac[i] = (fac[i-1] * i) % p;  //预处理阶乘   
}    
LL Lucas(LL n,LL m,LL p)
{    
    LL ret = 1;    
    while(n && m){    
        LL a = n%p, b = m%p;    
        if(a<b) 
			return 0;    
        ret = (ret * fac[a] * PowMod(fac[b] * fac[a-b] % p, p-2, p)) % p;    
        n /= p;    
        m /= p;    
    }    
    return ret;    
}