数论
前言:OI-WiKi 和 知乎 里都讲得挺好, 如果哪里没理解一定要弄懂,这里只是总结归纳。
1. 唯一分解定理等
约束个数:
约束和:
2. 欧拉函数
\(\varphi(n)\) 表示 \(n\) 以内与 \(n\) 互质的数的个数。
式子
性质
- \(φ(1)=1\).
- 若p是一个素数,则 \(φ(p)=p-1\).
- 若p是一个素数,则 \(φ(p^k)=(p-1)*p^{k-1}\).
- 对于任意两个正整数 \(a,b\), 且 \(gcd(a,b)=1\),则 \(φ(a*b)=φ(a)*φ(b)\).特别的,对于奇数n\(φ(2n)=φ(n)\).
变形
- p为素数,若 \(n\%p=0\), 则 \(φ(n*p)=p*φ(n)\).
- p为素数,若 \(n\%p≠p\), 则 \(φ(n*p)=(p-1)*φ(n)\).
- 与 \(n\) 互质的数都是成对出现的,且每对的和为 \(n\) ,所以大于2的数的 \(φ(n)\) 都为偶数.
代码
int phi_(int n)
{
int ans=n;
for(int i=2;i*i<=n;i++)
{
if(n%i==0)
{
ans=ans-ans/i;
while(n%i==0)
{
n/=i;
}
}
}
if(n>1) ans=ans-ans/n;
return ans;
}
int phi[maxn], cnt, p[maxn], vis[maxn];
void get_phi(int n)//求1到n的所有欧拉函数值
{
phi[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i])
{
p[cnt++]=i;//质数
phi[i]=i-1;//质数i的欧拉函数为i-1
}
for(int j=0;i*p[j]<=n;j++)
{
int m=i*p[j];
vis[m]=1;
if(i%p[j]==0)
{
phi[m]=p[j]*phi[i];
break;
}
else phi[m]=(p[j]-1)*phi[i];
}
}
}
3. 欧拉定理
若 \(gcd(a, m)=1\), 则 \(a^{φ(m)}≡1(mod\) \(m)\)。
4. 费马小定理(要求模数为素数)
若 \(m\) 为素数, 对于任意整数 \(a\),有\(a^m≡a(mod\) \(a)\).
特殊的,若 \(m\) 为素数,\(gcd(a, m)=1\), 则 \(a^{m-1}≡1(mod\) \(m)\).
5. 扩展欧拉定理
\(a^b=\begin{cases}a^{b\, mod\,φ(p)}&gcd(a,p)=1\\a^b&gcd(a,p)≠1,b<φ(p)(mod\,p)\\a^{b\, mod\,φ(p)+φ(p)}&gcd(a,p)≠1,b≥φ(p)\end{cases}\)
代码就不展示了
6. 扩展欧几里得
用于求 \(ax+by=gcd(a,b)\) 或者 \(ax≡b(mod\,p)\) 的一组整数解。
link:https://oi-wiki.org/math/number-theory/linear-equation/
这个真的非常重要。
构造通解
\(\begin{cases}x=x_0+\frac{b}{gcd(a,b)}\\y=y_{0}-\frac{a}{gcd(a,b)}\end{cases}\)
(考虑 \(ax+by=0\) 的构造)
示例
求不定方程 \(ax+by=c\) 的一组整数解。
- 若 \(gcd(a,b)\) 能被 \(c\) 整除,则有整数解。
先用扩展欧几里德求 \(ax+by=gcd(a,b)\) 的解,再乘以 \(c/gcd(a,b)\) ,即得原方程特解 \((x0,y0)\) 。 - 若 \(gcd(a,b)\) 不能被 \(c\) 整除 ,则无整数解。
代码
void exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1;
y=0;
return ;
}
exgcd(b, a%b, x, y);
int t=x;
x=y;
y=t-a/b*y;
}
int exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1;
y=0;
return a;
}
int ans=exgcd(b, a%b, x, y);
int k=x;
x=y;
y=k-a/b*y;
return ans;
}
函数返回值为gcd,
\(x, y\) 初值为0。
7. 乘法逆元
若 \(a\times x \equiv 1(mod\,b)\), 则称 \(x\) 为 \(a\) 在模 \(b\) 意义下的乘法逆元, 记为 \(a^{-1}\).
注意:并非所有的情况下都存在乘法逆元,当 \(gcd(a,b)=1\) 时,存在乘法逆元。
用法
求 \((a/b)\%p\) 等同于求 \(a\times (b\) 的逆元\()\%p\). (证明见老师ppt)
代码
下列p均为模数。
费马小定理求逆元(p为素数)
qpow(a, p-2, p);
线性求逆元
ny[1]=1;
for(int i=2;i<p;i++)
{
ny[i]=(long long)(p-p/i)*ny[p%i]%p;
}
扩展欧几里得求逆元
int exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1;
y=0;
return a;
}
int ans=exgcd(b, a%b, x, y);
int k=x;
x=y;
y=k-a/b*y;
return ans;
}
if(exgcd(b, p, x, y)==1)
cout<<(x%p+p)%p;
else
cout<<-1;
8. 组合数
\(f[i]\) 存 \(i!\) 的值。
\(g[i]\) 存 \(x!\) 的逆元的值。
也都挺好理解的。
代码
int qpow(int a,int b)
{
int res=1;
while(b)
{
if(b&1)
{
res=res*a%mod;
}
a=a*a%mod;
b>>=1;
}
return res;
}
void init()
{
f[0]=g[0]=1;
for(int i=1;i<=max(n, m)*4;i++)
{
f[i]=(f[i-1]*i)%mod;
g[i]=(g[i-1]*qpow(i, mod-2))%mod;
}
}
int getC(int n, int m)
{
if(n<m)
{
return 0;
}
return f[n]*g[m]%mod*g[n-m]%mod;
}
假如时间复杂度很紧,建议将预处理函数换成以下写法:
inline void init(int n)
{
f[0]=g[0]=1;
for(int i=1;i<=n;i++)
{
f[i]=f[i-1]*i%mod;
}
g[n]=qpow(f[n], mod-2);
for(int i=n-1;i>=0;i--)
{
g[i]=g[i+1]*(i+1)%mod;
}
}
掉过的坑:
- 如果无输出,很可能是N的值太大了。
- 如果输出为0,很可能没有使用init()。
- 如果模数很小,要考虑lucas定理。
- 当没有模数的时候,直接用组合的公式即可(\(f[n]/f[m]\%p/f[n-m]\%p\))。
- 如果不明情况TLE或RE,可以考虑换写法。
补充
可以用杨辉三角推出这样一个式子:
\(c[i][j]\) 表示从i个数里选j个的方案数。
就像这样子:
void init()
{
for(int i=0;i<=n;i++) c[i][0]=1;
for(int i=0;i<=n;i++)
{
for(int j=1;j<=i;j++)
{
c[i][j]=c[i-1][j-1]+c[i-1][j];
}
}
}
9. 错排问题
虽然课件里没有(好像有),但还是挺有意义的,记一下吧。
概要:有 \(n\) 个元素和 \(n\) 个位置,要求每一个元素都在错误的位置上的方案数(\(D(N)\))。(比如数列 \(1,2,3\), 错排方案有 ( \(2,1,3\) ), ( \(3,1,2\) ) 两种。)
然后通项公式为 \(D[i]=(i-1)\times(D[i-1]+D[i-2])\).
证明见这篇博客。
代码:直接预处理出d数组即可。
10. lucas 定理
当p是一个不大的质数时,
$ C_n^m\equiv C_{n/p}^{m/p}\times C_{n \bmod p}^{m \bmod p}(\bmod p) $
证明见老师ppt。
代码
int lucas(int n,int m, int p)
{
if(m==0) return 1;
return lucas(n/p, m/p, p)*getC(n%p, m%p, p)%p;
}
11. 二项式反演
记 \(f_n\) 表示恰好使用 \(n\) 个不同元素形成特定结构的方案数,\(g_n\) 表示从 \(n\) 个不同元素中选出 \(i \geq 0\) 个元素形成特定结构的总方案数。
若已知 \(f_n\) 求 \(g_n\),那么显然有:
若已知 \(g_n\) 求 \(f_n\),那么:
上述已知 \(g_n\) 求 \(f_n\) 的过程,就称为 二项式反演。
12. 波动数列
貌似不是数论??
波动数列:字面意思,每个元素大小波动起伏的数列。(其实我也不知道这是个什么数列,如果有明确定义麻烦告诉我)
性质
- 在一个波动数列中,若两个 \(i\) 与 \(i+1\) 相邻,那么我们直接交换这两个数字就可以组成一个新的波动数列。
- 把波动数列中的每个数字 \(a[i]\) 变成 \((n+1)-a[i]\) 会得到另一个波动数列,且新数列的山峰与山谷情况相反。
- 波动序列有对称性。
13. CRT/中国剩余定理
小学奥数---------韩信点兵
规范一点:
crt就是求解如下形式一元线性同余方程组
\(\begin{cases}x≡r_1(mod\,m_1)\\x≡r_2(mod\,m_2)\\......\\x≡r_n(mod\,m_n)\end{cases}\)
求x的最小整数解。其中模数两两互质。
过程+证明见老师ppt。
代码
void exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1;
y=0;
return ;
}
exgcd(b, a%b, x, y);
int t=x;
x=y;
y=t-a/b*y;
}
int crt(int m[], int r[])
{
int M=1, ans=0;
for(int i=1;i<=n;i++) M=lcm(M, m[i]);
for(int i=1;i<=n;i++)
{
int c=M/m[i], x, y;
exgcd(c, m[i], x, y);
ans=(ans+r[i]*c*x%M)%M;
}
return (ans%M+M)%M;
}
那么对于模数不互质时,
14. EXCRT/扩展中国剩余定理
就是在crt的基础上不保证模数两两互质。
int exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1;
y=0;
return a;
}
int ans=exgcd(b, a%b, x, y);
int t=x;
x=y;
y=t-a/b*y;
return ans;
}
int excrt(int m[], int r[])
{
int m1, m2, r1, r2, p, q;
m1=m[1], r1=r[1];
for(int i=2;i<=n;i++)
{
m2=m[i], r2=r[i];
int d=exgcd(m1, m2, p, q);
if((r2-r1)%d) return -1;//无解
p=p*(r2-r1)/d;
p=(p%(m2/d)+m2/d)%(m2/d);
r1=m1*p+r1;
m1=m1*m2/d;
}
return (r1%m1+m1)%m1;
}
关于CRT的用处:将不是质数的模数分解,最后crt合并。注意数组大小以及循环变量的大小。
15. 卡特兰数
upd on. 24.7.27
为了不让它成为遗忘的过去,我一定要再学一遍!
卡特兰数时一个用于计数问题的序列,就像斐波那契数列一样,它的前几项是:\(1,1,2,5,14,42,132,429,1430,4862......\) (遇到不会的计数题可以自己手算一下前几项的值,如果和卡特兰数一样的话直接套公式啦!)
通项公式
- \(H_n=C_{2n}^{n} -C_{2n}^{n-1}\)
- $H_n=\frac{C_{2n}^{n}}{n+1} $
- \(H_n=\frac{4n-2}{n+1}H_{n-1}\)
递推式
\(H_n=H_0*H_{n-1}+H_1*H_{n-2}+...+H_{n-1}*H_0\)
特征
一种操作数不能超过另一种操作数、两种操作不能有交集、给定一个移动范围,通常是卡特兰数。
16. Prufer序列
\(prufer\)序列是将一棵无根树转化为一个长度为 \(n-2\) 的序列。其中序列的值是无根树除叶节点的节点。
具体操作
(1) 每次选择编号最小的叶节点,将它和它与父节点的连边删去,将它的父节点存入序列。
(2) 循环操作(1),直至无根树里仅剩两个节点。
代码
将无根树转换为\(prufer\)序列:
void solve1()
{
for(int i=1;i<n;i++)
{
cin>>f[i];
deg[f[i]]++;
}
int cur=1;//最小叶节点
for(int i=0;i<=n-2;)
{
while(deg[cur])
{
cur++;
}
p[++i]=f[cur];
while(!--deg[p[i]]&&p[i]<cur)
{
p[i+1]=f[p[i]];
i++;
}
cur++;
}
ans=0;
for(int i=1;i<=n-2;i++)
{
ans^=(1ll*i*p[i]);
}
cout<<ans<<endl;
}
将无根树转换为\(prufer\)序列:
void solve2()
{
for(int i=1;i<=n-2;i++)
{
cin>>p[i];
deg[p[i]]++;
}
int cur=1;
for(int i=0;i<=n-2;)
{
while(deg[cur])
{
cur++;
}
f[cur]=p[++i];
while(!--deg[p[i]]&&p[i]<cur)
{
f[p[i]]=p[i+1];
i++;
}
cur++;
}
f[p[n-2]]=n;
ans=0;
for(int i=1;i<=n-1;i++)
{
ans^=(1ll*i*f[i]);
}
cout<<ans<<endl;
}
时间复杂度为 \(O(n)\)。
性质
- \(prufer\)序列与无根树一一对应。
- 度数为 \(d[i]\) 的节点在\(prufer\)序列中出现 \(d[i]-1\) 次。
- 一个 \(n\) 个节点的完全图的生成树个数为 \(n^{n-2}\)。
- 对于给定度数为 \(d_{1∼n}\) 的无根树有\(\frac{(n-2)!}{\coprod_{i=1}^n(d[i]-1)!}\)种情况。
17. Bsgs算法
给定整数 \(a,b,p\), \(a,p\) 互质。
求满足 \(a^x\equiv b(mod\,p)\) 的最小非负整数 \(x\) 。
具体操作
由扩展欧拉定理 \(a^x\equiv a^{x\,mod\,φ(p)}(mod\,p)\) 可知\(a^x\) 模p意义下的最小循环节为 \(φ(p)\) ,因为 \(φ(p)<p\) ,在0~p范围内,一定能找到最小的 x。
令 \(x=im-j\) ,其中 $ m=ceil(\sqrt{p})$ ,\(i\) 在 \(1至m\) 范围内,\(j\) 在 \(0至m-1\) 范围内。
则 \(a^{im-j}\equiv b(mod\,p)\), 即 \((a^m)^i\equiv ba^j\)。
- 先枚举j,把(\(ba^j,j\))丢入哈希表,结果相同,则取较大的j。
- 再枚举i,计算 \((a^m)^i\) ,到哈希表中查找是否有相等值,找到第一个结束。最小的 \(x=im-j\)。
代码
int bsgs(int a,int b,int p)
{
a%=p;
b%=p;
if(b==1) return 0;
int m=ceil(sqrt(p)),t=b;
h[b]=0;
for(int i=1;i<m;i++)
{
t=t*a%p;
h[t]=i;
}
int mi=1;
for(int i=1;i<=m;i++)
{
mi=mi*a%p;
}
t=1;
for(int i=1;i<=m;i++)
{
t=t*mi%p;
if(h.count(t))
return i*m-h[t];
}
return -1;
}
时间复杂度为\(O(\sqrt{p})\)。
18. exBsgs算法
假如 \(a,p\) 不互质,那怎样求 \(a^x\equiv b(mod\,p)\) ?
这个思路挺好想的,代码如下:
代码
int inv(int a,int p)
{
int x=0, y=0;
exgcd(a, p, x, y);
return (x%p+p)%p;
}
int exbsgs(int a,int b,int p)
{
a%=p, b%=p;
if(b==1||p==1) return 0;
if(!a)
{
if(b) return -1;
return 1;
}
int k=0, t=1;
while(__gcd(a, p)!=1)
{
int gcd=__gcd(a, p);
if(b%gcd) return -1;
++k, b/=gcd, p/=gcd, t=t*(a/gcd)%p;
if(t==b) return k;
}
int ans=bsgs(a, b*inv(t, p)%p, p);
if(ans==-1) return -1;
return ans+k;
}