浅谈数论在oi中应用
数论的水不深,但是你把握不住
1.快速幂
思路:
如果 \(b=c+d\)
那么我们可以将 \(a^b\)转化为\(a^c*a^d\) 小学生都会的公式
在代码中我们可以在for循环中设一个临时变量i,用来记\(a^i\)并判断\(i\)是否在\(b\)中,如果是将答案乘上\(a^i\)即可,同时b减去i
原理:将指数拆成多个指数(过程中转化为二进制进行计算)
代码:
int powt(int x,int y){
int res=1,xt=x;
while(y){
if(xt*res>p) return -1;
if(y&1) res*=xt;
y>>=1;
xt*=xt;
}
return res;
}
2.素数筛
作用:
求出一段区间内所有的素数
埃氏筛
思路:
将一个素数的所有倍数进行标记(即为合数)
时间复杂度约为\(O(n log n)\)
这个思路比较简单,直接上代码
代码:
void get_primes1(int n) {//埃氏筛
memset(prime, false, sizeof prime);
cnt = 0;
for (int i = 2; i <= n; i++) {
if (!prime[i]) {
primes[cnt++] = i;
for (int j = i; j <= n / i; j++)
prime[j * i] = true;
}
}
}
欧拉筛
思路:
这个算法中我们设一个记录素数的变量primes,记录下素数后再标记primes数组中每个素数的倍数(即为合数)。
void get_primes2(int n) {//欧拉筛(不完整代码)
memset(prime, false, sizeof prime);
cnt = 0;
for (int i = 2; i <= n; i++) {
if (!prime[i]) primes[cnt++] = i;
for (int j = 0; primes[j] <= n / i; j++) {
prime[i * primes[j]] = true;
}
}
}
为了对其进行优化,我们在标记合数循环中添加了一个判断语句:
if (i % primes[j] == 0) break;
这个判断语句其实就是判断primes[j]是否在i中如果在i中证明i中已经包含了接下来的数的最小质因子(即为primes[j]),但是继续往下循环primes[j]会不断变大,那么就没有循环的必要了,就要退出循环。
这里我们说了primes[j]为i中的质因子所以i\(\times\)primes[j+k]可以用后期循环和primes[j]来标记
所以这个代码几乎每个数只跑过一遍,时间复杂度为:\(O(n)\)
完整代码:
void get_primes2(int n) {//欧拉筛
memset(prime, false, sizeof prime);
cnt = 0;
for (int i = 2; i <= n; i++) {
if (!prime[i]) primes[cnt++] = i;
for (int j = 0; primes[j] <= n / i; j++) {
prime[i * primes[j]] = true;
if (i % primes[j] == 0) break;
}
}
}
3.质因数分解
思路:
外层循环判断合数n中是否有\(i\)这个质因数,若有就标记它。
但是\(n\)中可能会有多个\(i\),所以我们就要进行一个循环直到\(n\)没法整除i为止,这个思路比较简单,直接上代码
代码:
pair<int, int> factor[N];
int idx;
void divide(int n) {//质因数分解
idx = 0;
for (int i = 2; i <= n / i; i++) {
if (n % i == 0) {
int k = 0;
while (n % i == 0) k++, n /= i;
factor[idx++] = {i, k};
}
}
if (n > 1) factor[idx++] = {n, 1};
}
4.阶乘\(n!\)的质因数分解
n!=1 * 2 * 3...(n-1)n
思路:
我们首先用欧拉筛筛出所有\(n\)范围内的素数记录到primes数组中。
现在,\(n!\)因数中(指的是1到n)有的数可以分解为x个有的只能分解为x-1个,那么我们就每次除一边,来计算这个值,思路较简单,直接上代码
代码:
void divide_fac(int n) {//阶乘n!的质因数分解
get_primes2(n);
idx = 0;
for (int i = 0; i < cnt; i++) {
int p = primes[i];
int s = 0;
for (int j = n; j; j /= p) s += j / p;
factor[idx++] = {p, s};
}
}
5.最大公约数与最小公倍数
辗转相除法求解
此算法较简单直接上代码
最大公约数:
int GCD(int a, int b)
{
int c;
while (b > 0)
{
c = a % b;
a = b;
b = c;
}
return a;
}
最小公倍数:
两数相乘然后除以最大公约数即为最小公倍数
int LCM(int a,int b)
{
int c;
c = a * b / GCD(a, b);
return c;
}
特殊定理:存在A,B gcd(A,B)=x A=ax B=bx lcm(A,B)=xab a与b互质
那么 AB=abx2=gcd(A,B)lcm(A,B)
推导过程:
假设2数为A,B 可以拆分成A=ax ;
B=bx a&b互质 那么,AB的最小公倍数为abx AB 的最大公约数为 x2者乘起来正好等于A*B
(此过程来源于百度)
扩展欧几里得/裴蜀定理
思路:
扩展欧几里得/裴蜀定理说明了对任何整数 a、b和它们的最大公约数 d ,关于未知数 x以及 y 的线性的丢番图方程。
若\(a,b\)是整数,且\(gcd(a,b)=d\),那么对于任意的整数\(x,y,ax+by\)都一定是d的倍数,特别地,一定存在整数\(x\),\(y\),使\(ax+by=d\)成立。根据我们上面的gcd可以推导出在下一次gcd中a变为b,b变为a%b,即x为y,y为x%y
根据常识
a%b=a-a/b*b
可得:
bx+(a-a/b * b)y=d
bx+ay-(a/b * b)y=d
ay+[x-(a/b)*y]b=d
这我们就推出返回的\(x\)就是\(y\),\(y\)为\(x-(a/b)*y\)
代码:
ll gcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1;
y=0;
return a;
}
ll d=gcd(b,a%b,x,y);
ll t=x;
x=y;
y=t-(a/b)*y;
return d;
}
6.逆元
逆元是什么?
如果\(\frac{n}{a}\) \(mod\) \(p\) \(=\) \(n\times x\) \(mod\) \(p\) 且 \(x \in Z\)
那么我们称 \(x_{min}\) 为\(a\)在模\(p\)意义下的逆元(记作\(x_{min}=a^{-1}\))
为什么a和p不互质就不存在逆元:
证:
令c=gcd(a,p)且c>1(即a,p不互质)
ax=1(mod p)
等价于:
ax+py=1
∵c|a,p
∴c|(ax+py)
∴c|1
∵c>1
所以a,p不存在逆元
证毕
做法1(费马小定理):
我们知道\(a\) 和 \(p\)的情况下:
\(a^{p-1}\) \(mod\) \(p\) \(=\) \(1\) \(mod\) \(p\)
费马小定理证明
那么我们可以得出:
\((a\times a^{p-2})\) \(mod\) \(p\) \(=\) \(1\) \(mod\) \(p\)
则b的逆元为\((a^{p-2}\) \(mod\) \(p)\)为\(a\)在模\(p\)意义下的逆元
即可用快速幂求解
代码:
//注:根据不同的题目请自行开long long
int power(int x,int y,int p){
int t=x;
int ans=1,it=1;
while(y){
if(y&it){
y^=it;
ans*=t;
ans%=p;
}
it<<=1;
t*=t;
t%=p;
}
return ans%p;
}
int fpm(int a,int p){//费马小定理求逆元
return power(a,p-2,p)%p;
}
时间复杂度为:\(O(log_{2}n)\)
此做法由于复杂度较大,且常数较大不建议使用。
P3811 【模板】乘法逆元 在此题中最后两个点可以卡死费马小定理(注意开longlong
做法2(扩展欧几里得):
要求\(a\)在模\(p\)意义下的逆元我们可以将:
\(ax=1(mod\) \(p)\)
转化成:
\(ax-py=1\)
没错你没看错这就是上面我们讲的扩展欧几里得/裴蜀定理(扩展欧几里得)
所以就按照扩展欧几里得的模板求出上面式子中的\(x\)即可
而且这个做法还可以判断最大公约数是否为1(是否存在逆元)
代码:
//根据个人需要改longlong
int exgcd(int a,int b,int &x,int &y){
if(!b){
x=1;
y=0;
return a;
}
int d=exgcd(b,a%b,x,y);
int t=x;
x=y;
y=t-(a/b)*y;
return d;
}
int inv(int a,int p){//扩展欧几里得求逆元
int x,y;
int ans=exgcd(a,p,x,y);
if(ans>1) return -1;//a和p不互质
else return (x%p+p)%p;//防止负数的情况
}
这种求法非常高效,而且常数较小,用与求单个数的逆元非常方便。
然后你再交一下这题:P3811 【模板】乘法逆元
你会发现,最后一个点还是T,吸氧也会T,不要慌,看这么多输出就知道事情不简单,此时就要用复杂度(复杂度假了)最高的算法——递推。
做法3 递推法:
令\(p=ki+r;\){\(k=\lfloor\frac{p}{i}\rfloor,r=p\) \(mod\) \(i\)}(\(i<p,k<p,r>i\))
则有\(ki+r≡0(mod\) \(p)\)①
①式乘\(i^{-1}*r^{-1}\)得:
\(k*r^{-1}+i^{-1}≡0(mod\) \(p)\)
移项得:
\(i^{-1}≡-k*r^{-1}(mod\) \(p)\)
代入\(k=\lfloor\frac{p}{i}\rfloor,r=p\) \(mod\) \(i\):
\(i^{-1}≡-\lfloor\frac{p}{i}\rfloor*(p\) \(mod\) \(i)^{-1}(mod\) \(p)\)
由于\((p\) \(mod\) \(i)<i\)
所以,我们已经求出\((p\) \(mod\) \(i)^{-1}\)
用数组inv[i]记录\(i^{-1}\)(\(i\)的逆元)
则inv[i]=\(-\lfloor\frac{p}{i}\rfloor*inv[p\) \(mod\) \(i]\) \(mod\) \(p\)
证毕;
注:为防止出现逆元为负数在上面式子中要加上p
代码:
//AC_664ms
#include<cstdio>
#include<iostream>
#define ll long long
using namespace std;
const int maxn=3e6+5;
ll inv[maxn]={0,1};
int main(){
int n,p;
scanf("%d%d",&n,&p);
printf("1\n");
for(int i=2;i<=n;i++)
inv[i]=(ll)p-(p/i)*inv[p%i]%p,printf("%d\n",inv[i]);
return 0;
}
这样再提交就AC了
这种做法复杂度较高但是可以批量求逆元,在这类题目中非常好用。
做法4 欧拉定理:
首先我们要知道欧拉定理是什
这里我们需要引入一个函数——欧拉函数
至于这里的缩集指的就是小于等于n且与n互质的数的个数
那么欧拉定理就是
特别的当p为质数时欧拉定理将变为费马小定理
所以欧拉定理相对于费马小定理来说最大的优势便是能求出模数非质数的逆元
易证得a在模p意义下的逆元即为\(a^{\varphi(p)-1}\)显然这可以用快速幂求出来
但欧拉定理的核心在于如何求欧拉函数\(\varphi\)
O(logn)级的复杂度求欧拉函数
对于一个数n,其可以分解为形如:
的形式。
对于一个质数p它的欧拉函数显然是p-1(p!=1)。
那么对于\(p^a\)我们可以用分块的思想求解
我们不难发现对于一个模数p如果p为质数欧拉函数的值显然等于p-1,对于合数p可以在素数筛时用以素数的欧拉函数求出合数欧拉函数,具体实现和思想见代码
#include<cstdio>
int phi[100100],prime[30000];//nprime[100100];
int n,len;
int main()
{
int i,j;
n=50;
phi[1]=1;
//nprime[1]=1;
for(i=2;i<=n;i++)
{
//if(!nprime[i])
if(!phi[i])
{
prime[++len]=i;
phi[i]=i-1;
}
//printf("%d: ",i);
for(j=1;j<=len&&i*prime[j]<=n;j++)
{
//nprime[i*prime[j]]=1;
//printf("%d ",i*prime[j]);
if(i%prime[j]==0)//情况1
{
phi[i*prime[j]]=phi[i]*prime[j];
break;
}
else
phi[i*prime[j]]=phi[i]*(prime[j]-1);
}
//puts("");
}
return 0;
}
然而这个做法复杂度显而易见O(n)显然用O(n)的复杂度求逆元会寄掉(然而这种n线性的复杂度可以解决一些需要预处理的问题
7.中国剩余定理(CRT)
直接看题
中国剩余定理又称孙子定理(所以曹冲和猪为什么被迫害了)
题意就是给了一堆同余方程,让求x的值:
最佳阅读效果(这里引用了OI Wiki的证明)
中国剩余定理 (Chinese Remainder Theorem, CRT) 可求解如下形式的一元线性同余方程组(其中 \(n_1, n_2, \cdots, n_k\) 两两互质):
上面的「物不知数」问题就是一元线性同余方程组的一个实例。
算法流程
- 计算所有模数的积 \(n\);
- 对于第 \(i\) 个方程:
- 计算 \(m_i=\frac{n}{n_i}\);
- 计算 \(m_i\) 在模 \(n_i\) 意义下的 逆元 \(m_i^{-1}\);
- 计算 \(c_i=m_im_i^{-1}\)(不要对 \(n_i\) 取模)。
- 方程组在模 \(n\) 意义下的唯一解为:\(x=\sum_{i=1}^k a_ic_i \pmod n\)。
代码:
// C++ Version
LL CRT(int k, LL* a, LL* r) {
LL n = 1, ans = 0;
for (int i = 1; i <= k; i++) n = n * r[i];
for (int i = 1; i <= k; i++) {
LL m = n / r[i], b, y;
exgcd(m, r[i], b, y); // b * m mod r[i] = 1
ans = (ans + a[i] * m * b % mod) % mod;
}
return (ans % mod + mod) % mod;
}
算法证明:
我们需要证明上面算法计算所得的 \(x\) 对于任意 \(i=1,2,\cdots,k\) 满足 \(x\equiv a_i \pmod {n_i}\)。
当 \(i\neq j\) 时,有 \(m_j \equiv 0 \pmod {n_i}\),故 \(c_j \equiv m_j \equiv 0 \pmod {n_i}\)。又有 \(c_i \equiv m_i \cdot (m_i^{-1} \bmod {n_i}) \equiv 1 \pmod {n_i}\),所以我们有:
即对于任意 \(i=1,2,\cdots,k\),上面算法得到的 \(x\) 总是满足 \(x\equiv a_i \pmod{n_i}\),即证明了解同余方程组的算法的正确性。
因为我们没有对输入的 \(a_i\) 作特殊限制,所以任何一组输入 \(\{a_i\}\) 都对应一个解 \(x\)。
另外,若 \(x\neq y\),则总存在 \(i\) 使得 \(x\) 和 \(y\) 在模 \(n_i\) 下不同余。
**故系数列表 \(\{a_i\}\) 与解 \(x\) 之间是一一映射关系,方程组总是有唯一解。
扩展中国剩余定理(ex_CRT)
与 CRT 的不同
CRT 由于要有求逆元的操作,所以有模数互质的限制。
当遇到模数不互质的情况,普通的办法求CRT就不可行了,所以要用到ex_CRT。 (但ex_CRT基本不考)
CRT和ex_CRT看似功能一样,但实际上算法完全不一样,CRT是基于构造实现的,而exCRT是基于exgcd实现的。
Solution
具体的,对于
可以变成
其中 \(k_1,k_2\in\mathbb{Z}\),考虑什么样的 \(x\) 满足这俩式子,两个式子相减可得
不难发现,这个式子长的和exgcd很像,也就说这两个式子有解的条件是 \(\gcd(m_1,m_2)|a_2-a_1\)
考虑如何求出 \(k_1,k_2\) ,等式两边同时除 \(\gcd(m_1,m_2)\),得到
显然这时 \(k_1,(-k_2)\) 的系数是互质的,于是乎我们可以套扩展欧几里得算法求出一组 \(t_1,t_2\) 满足。
等式两边同时乘 \(\frac{a_2-a_1}{\gcd(m_1,m_2)}\),就能得到 \(k_1,(-k_2)\),具体的
带入 \(x=a_1+k_1m_1=a_2+k_2m_2\),这样我们就解出了这两个式子的解 \(x_0\)。即
$ x\equiv x_0\pmod{\text{lcm}(m_1,m_2)} $
同理我们可以不断合并两个式子的解以求出
ex_CRT--upd on 22/11/23
数论函数
——返回值和参数均为非负整数的函数
积性函数
积性函数定义
积性函数指对于所有互质的整数a和b有性质f(ab)=f(a)f(b)的数论函数。
常见の积性函数
狄利克雷卷积
定义
设两个数论函数函数\(f(n)\quad g(n)\)
那么他们两个的卷积为
记作
常见的狄利克雷卷积函数公式
1:
证明
2:
3:
狄利克雷卷积还没搞明白的说,所以没写