数论
Tags:Noip前的大抱佛脚
知识点
Exgcd
\(O(logn)\)求解\(Ax+By=C\)的问题
1、若\(C\%gcd(A,B)!=0\)则无解
2、\(Gcd=gcd(A,B);A/=Gcd,B/=Gcd,C/=Gcd\)
3、代入下面代码求\(Ax+By=1\)
4、\(x*C\),得到一组特解
5、通解为\(\begin{cases}x=x_0+k*B \\y=y_0+k*A\end{cases}\)
void Exgcd(ll a,ll b,ll &x,ll &y)
{
if(!b){x=1;y=0;return;}
Exgcd(b,a%b,y,x);y-=a/b*x;
}
逆元
在\(gcd(A,P)==1\)时,\(A\)在模\(P\)意义下存在逆元(证明可以符合Exgcd有解的证明),其余情况不存在逆元
通常\(P\)为质数时就是\(A^{p-2}\)作为逆元(费马小定理)
欧拉定理(费马小定理) \(x^{-1}\equiv x^{p-2}(mod\ p)\),要求\(p\)为质数
解方程(Exgcd)\(Px+Ay=1,A=\frac{1}{y}(mod\ P)\)
-
线性递推逆元:
令\(a=\lfloor\frac{p}{i}\rfloor,b=p\%i\),则有\(ai+b=p,i\equiv-\frac{b}{a}(mod\ p),\frac{1}{i}\equiv-\frac{a}{b}(mod\ p)\)
由于\(b<i\),所以\(b\)的逆元已经求出,直接可以得到:
inv[i]=p-(p/i*inv[p%i])%p;
gcd
具有一些奇妙的性质,如可合并
- 往数集中加入一个数,要么gcd不变,要么至少变为原来的\(\frac{1}{2}\)(所以可以用来分块了(同样的xor也每次只会增加log次))
欧拉函数\(\varphi(x)\)
表示小于x且与x互质的数的个数
计算公式
#### **欧拉公式**
$$a^{\varphi(p)}\equiv 1\ (mod \ p)\]
降幂公式
- \(x,p\)互质:\(x^k\equiv x^{k\%\varphi(p)}\ (mod\ p)\)
- \(x,p\)不互质:\(x^k\equiv \begin{cases}x^{k}\ (mod\ p),k\le\varphi(p) \\x^{k\%\varphi(p)+\varphi(p)}\ (mod\ p),k>\varphi(p)\end{cases}\)
CRT&EXCRT
\(O(nlogn)\)求解一系列同余方程,如
\]
代码如下(未判无解,无解即exgcd无解)
ll gcd(ll a,ll b) {return !b?a:gcd(b,a%b);}
void Exgcd(ll a,ll b,ll &x,ll &y)
{
if(!b) {x=1;y=0;return;}
Exgcd(b,a%b,y,x);y-=a/b*x;
}
int EXCRT()
{
for(int i=2;i<=n;i++)
{
ll g=gcd(P[i-1],P[i]),C=(A[i]-A[i-1])/g,x,y;
Exgcd(P[i-1]/g,P[i]/g,x,y);
P[i]=P[i-1]*P[i]/g;
A[i]=A[i-1]+P[i-1]*x*C;
A[i]=(A[i]%P[i]+P[i])%P[i];
}
return A[n];
}
方法是每次合并两个方程,手推式子\(P[i]x+A[i]=P[i+1]y+A[i+1]\ \ \ ->\ \ \ P[i]x+P[i+1]y=A[i+1]-A[i]\)
可以发现模数为原来的\(lcm\)(先除以\(gcd\)再乘,很容易爆\(long\ long\)),余数为原来余数加上模数的\(x\)倍
注意经常x算出来是负数所以时刻\(+mod)\%mod !\)
有的时候题目并没有那么裸
- \(x\)前带系数(屠龙勇士)
\(Ax\equiv B(mod\ p)\ \ ->\ \ Ax+pk=B(当且仅当B\%gcd(A,p)==0时有解)\ \ ->\ \ x=x_0+tp\ \ ->\ \ x\equiv x_0(mod\ p)\)
- CRT合并答案
对于一些公式,只适用于模数是质数的情况,而题目要求的模数要是任意正数,所以可以求得
\]
然后用CRT合并答案,求得\(Ans\equiv OUTANS\ (mod\ Mod)\)
但是就所见到的题目来说(任意模数NTT),只有在答案不太大的时候适用,例如该题答案不会超过\(10^9*10^9*10^5=10^{23}\),方法是选取三个乘积大于\(10^{23}\)的NTT模数,得到
\]
令\(M=P_1*P_2\)于是\(Ans=kM+A_1=k_3P_3+A_2\),即\(kM+A_1\equiv A_2(mod\ P_3)\)
由于\(Ans\le10^{23}<P_1P_2P_3=MP_3\),所以\(k<P_3\),根据此同余方程可以求得\(k\)在\(mod\ P_3\)的解也就是\(k\)的实际值
然后就可以带入求得答案\(Ans=kM+A_1\)了,不过会爆\(long\ long\),可以采用这个
ll mul(ll x,ll y,ll p) {return (x*y-(ll)((long double)x/p*y+0.5)*p+p)%p;}
BSGS&EXBSGS
快速(\(O(\sqrt n)\))求解\(A^x=B(mod\ P)\)的\(x\)的解,板子题:[SDOI2013]随机数生成器
普通的BSGS要求P是质数,拓展的可以不用是质数
通过降幂公式,能够知道\(x\le P-1\),令\(M=\sqrt P\),则\(x=iM+j\),可以枚举\(i\)的值,得到\(t=A^{(i+1)M}\),查表看是否存在\(T=A^{M-j}B\)满足条件,得到解就返回,这样就能保证解是最小的了
int BSGS(int A,int B,int P)
{
if(!A%P) return -1;
int M=sqrt(P)+1;Hash.reset();
for(int i=0,t=B;i<M;i++,t=1ll*t*A%P) Hash.Add(t%Mo,t,i);//存入B*A^i的哈希值
for(int i=1,bs=ksm(A,M,P),t=bs;i<=M;i++,t=1ll*t*bs%P)//t=A^{(i+1)M}
if(Hash.Query(t)!=-1) return i*M-Hash.Query(t);
return -1;
}
FFT/NTT/MTT/FWT
板子见最下方
原根相关
定义:
原根类似于FFT的单位复根,使得能够进行取模操作
NTT模数是指有原根且是\(r×2^k+1\)的形式的模数,如\(998244353(3)\)、\(1004535809(3)\)
存在性及判定:
一个数有原根当且仅当它为\(2,4,p,2p,p^r\)(p为奇质数)
\(g^1,g^2...g^{\phi(p)}\)在\(mod\ p\)下各不相同,且\(g^{\phi (p)}=1\),当p为质数时,也就是说\(g^1...g^{p-1}\)分别对应\(1...p-1\)
在\(p\)为质数时,令\(p-1=p_1^{a1}p_2^{a_2}...p_n^{a_n}\),若存在\(g^{\frac{p-1}{p_i}}=1\)则\(g\)不是原根,否则是原根(证明)
注意点
- 两个多项式最高次项分别为\(l1,l2\),那么需要开的长度(包括0)是\(>l1+l2\)的第一个形如\(2^k\)的数
组合公式
-
从\((0,0)\)走到\((n,m)\),不碰到直线\(y=x+b\)的方案数:\(C(n+m,n)-C(n+m,n-b)\)
意义为沿着直线翻折,越过直线的路径对称域从翻着顶点到终点的路径
斯特林数
第一类斯特林数
第一类斯特林数:\(n\)个人坐\(m\)张圆桌的方案数(人不同,圆桌相同,n个元素构成m个圆排列)
\]
第二类斯特林数
第二类斯特林数:\(n\)个球放入\(m\)个盒子的方案数(球不同,盒子相同,n个元素分成m个非空集合)
\]
卡塔兰数
\]
常用数学公式
- 平方和的前缀和 \(\sum_{i=1}^{n}i^2=\frac{n(n+1)(2n+1)}{6}\)
- 立方和的前缀和 \(\sum_{i=1}^{n}i^3=(\frac{n(n+1)}{2})^2\)
技巧经验
容斥
-
从\((0,0)\)到\((n,m)\)且不经过一些禁点的方案数(BZOJ两双手)
\(f[i]=\sum_{j=0}^{j<i}-f[j]*way(j,i),f[0]=-1,way(j,i)\)表示从j点到i点的方案数,按照从\((0,0)\)走到该点的所需步数排序
组合计数
- 可以画成二维平面上的点去移动,辅助理解、推式子
区间筛
-
区间筛质数个数
luogu1835素数密度:筛出长度为\(10^6\)的一段区间内的质数个数,左右端点在int以内
方法:枚举到\(\sqrt r\)然后用\(j=(l-1)/i*i+i\)暴力在这个区间内筛掉合数
-
区间筛\(\varphi()\)之和
luogu3601签到题:筛出长度为\(10^6\)的一段区间的\(\varphi(i)\)之和,左右端点在\(10^{12}\)以内
方法:枚举到\(\sqrt r\)然后\(ans[i]=ans[i]/p*(p-1)\)
-
区间筛\(\mu()\)之和
luogu3653小清新数学题:筛出长度为\(10^6\)的一段区间的\(\mu(i)\)之和,左右端点在\(10^{18}\)以内
方法:枚举到\(10^6\)然后用质数暴力搞,除剩下的数如果不是\(1\)那么只有三种情况:
质数的平方:用\(sqrt\)判掉
一个质数
两个质数之积
上述两种可以用Miller—Rabin素性判断判掉,不过更方便的是\(check\):\(2^{p-1}\%p==1\)、\(3^{p-1}\%p==1\)。这题判两次可以过,当然是有可能不准确的。。
博弈
- Nim游戏:异或和为0则先手必败,否则必胜
- 威佐夫博弈:约定\(a<b\),\(a=\frac{\sqrt 5+1}{2}(a-b)\)时是奇异局势,先手必败。
有趣的式子
gcd有关
-
\(\sum_{i=1}^{n}gcd(i,n)\)
\(=\sum_{d=1}^{n}(d\sum_{i=1}^{n}[gcd(i,n)==d])\)
\(=\sum_{d=1}^{n}d\varphi(\frac{n}{d})\)(当然到这步就可以\(O(n\sqrt n)\)做了)
\(=\sum_{d=1}^{n}\frac{n}{d}\varphi(d)\)
\(=\sum_{d=1}^{n}\frac{n}{d}d\prod_{p_i|d}\frac{p_i-1}{p_i}\)
\(=n\sum_{d=1}^{n}\prod_{p_i|d}\frac{p_i-1}{p_i}\)
考虑每个质数选或者不选的生成函数:\((b_i\frac{p_i-1}{p_i}+1)\),其中\(b_i\)表示\(p_i\)的指数,含义为选了这个质数就会有\(b_i\)种指数
所以最后答案就是\(n\prod_{p_i|n}(b_i\frac{p_i-1}{p_i}+1)\)
数论模板库
黑科技
\(long\ long\)相乘取模
ll mul(ll x,ll y,ll P) {return (x*y-(ll)((long double)x/P*y+0.5)*P+P)%P;}
子集枚举
\(0\)~\(2^n\)的子集个数之和是\(3^n\)(一个位置可以有三种情况:没被选到、选到但是没有被子集枚举到、选到并被枚举到)
以下可以快速枚举子集
for(int i=s;;i=(i-1)&s)
{
cout<<i<<endl;
if(!i) break;
}
高维前缀和
统计子集
for(int p=1;p<=1<<20;p<<=1)
for(int i=0;i<1<<21;i++)
if(i&p) f[i]+=f[i^p];
统计超集
for(int p=1;p<=1<<20;p<<=1)
for(int i=0;i<1<<21;i++)
if(!(i&p)) f[i]+=f[i|p];
各种线性筛
用接近线性的方法筛取一些函数,首先要先找到其值与质数的关联才能做
积性函数都能线性筛
线性筛质数
【模板】线性筛素数
void Get_pri()
{
ispri[1]=1;
for(int i=2;i<=N;i++)
{
if(!ispri[i]) pri[++tot]=i;
for(int j=1;j<=tot&&i*pri[j]<=N;j++)
{
ispri[i*pri[j]]=1;
if(i%pri[j]==0) break;
}
}
}
线性筛\(\mu\)
void Get_mu()
{
ispri[1]=mu[1]=1;
for(int i=2;i<=N;i++)
{
if(!ispri[i]) pri[++tot]=i,mu[i]=-1;
for(int j=1;j<=tot&&i*pri[j]<=N;j++)
{
ispri[i*pri[j]]=1;
if(i%pri[j]) mu[i*pri[j]]=-mu[i];
else break;
}
}
}
线性筛\(\varphi\)
void Get_phi()
{
ispri[1]=phi[1]=1;
for(int i=2;i<=N;i++)
{
if(!ispri[i]) pri[++tot]=i,phi[i]=i-1;
for(int j=1;j<=tot&&i*pri[j]<=N;j++)
{
ispri[i*pri[j]]=1;
if(i%pri[j]) phi[i*pri[j]]=phi[i]*(pri[j]-1);
else {phi[i*pri[j]]=phi[i]*pri[j];break;}
}
}
}
高级算法
Exgcd
上文有讲用法,注意传入的时候要求都除以了\(gcd\)
void Exgcd(int a,int b,int &x,int &y)
{
if(!b) {x=1;y=0;return;}
Exgcd(b,a%b,y,x);y-=a/b*x;
}
Lucas
【模板】卢卡斯定理
用于计算模数比较小的时候的组合数,公式看就是\(C(n,k)=C(n/p,k/p)*C(n\%p,k\%p)\)
然后继续递归,边界条件是\(k\)为\(0\)时答案是\(1\)
int C(int n,int k) {return 1ll*jc[n]*inv[k]%P*inv[n-k]%P;}
int Lucas(int n,int k) {return !k?1:1ll*Lucas(n/P,k/P)*C(n%P,k%P)%P;}
int pre()
{
inv[0]=inv[1]=jc[0]=1;
for(int i=2;i<P;i++) inv[i]=(P-1ll*P/i*inv[P%i]%P)%P;
for(int i=1;i<P;i++) jc[i]=1ll*jc[i-1]*i%P,inv[i]=1ll*inv[i]*inv[i-1]%P;
}
EXCRT
【模板】扩展中国剩余定理(EXCRT)
合并同余方程,复杂度瓶颈在于exgcd
ll mul(ll x,ll y,ll p) {return (x*y-(ll)((long double)x/p*y+0.5)*p+p)%p;}
ll gcd(ll a,ll b) {return !b?a:gcd(b,a%b);}
void Exgcd(ll a,ll b,ll &x,ll &y) {!b?(x=1,y=0):(Exgcd(b,a%b,y,x),y-=a/b*x);}
ll n,B1,P1,B2,P2;
int main()
{
cin>>n>>P1>>B1;n--;
while(n--)
{
cin>>P2>>B2;
ll g=gcd(P1,P2),C=(B2-B1)/g,x,y,P;
Exgcd(P1/g,P2/g,x,y); P=P1/g*P2;
B1=(mul(mul(x,C,P),P1,P)+B1+P)%P; P1=P;
}
cout<<B1<<endl;
}
BSGS
[CQOI2018]破解D-H协议
\(g^a \equiv A(mod\ P)\),给定\(A\),求\(a\)。
方法:
- \(a=im-j\)所以\(g^{im}\equiv A*g^j\),其中\(m=\sqrt P\)能够保证最优复杂度
- 然后从\(0\)到\(m\)枚举\(i\),把\(g^{im}\)存进哈希表
- 暴力判断对于每一个\(j\),哈希表里是否有\(A*g^j\)的值
复杂度\(O(\sqrt PlogP)\)
const int Mo=100003;
int g,A,a,b,P,n;
int ksm(int x,int k)
struct HashTable
{
struct Line{int next,v1,v2;}a[Mo];
int head[Mo],cnt;
void reset() {memset(head,0,sizeof(head));cnt=0;}
void Add(int p,int v1,int v2) {a[++cnt]=(Line){head[p],v1,v2};head[p]=cnt;}
int Query(int x)
{
for(int i=head[x%Mo];i;i=a[i].next)
if(a[i].v1==x) return a[i].v2;return -1;
}
}Hash;
int BSGS(int g,int A,int P)
{
int M=sqrt(P)+1;Hash.reset();
for(int i=0,t=A;i<M;i++,t=1ll*t*g%P) Hash.Add(t%Mo,t,i);
for(int i=1,bs=ksm(g,M),t=bs;i<=M;i++,t=1ll*t*bs%P)
if(Hash.Query(t)!=-1) return i*M-Hash.Query(t);return -1;
}
高斯消元
【模板】高斯消元法
\(O(n^3)\)求解线性方程组。用线性基的方式实现
const int N=110;
const double eps=1e-7;
int n;
double F[N][N],A[N];
int main()
{
cin>>n;
for(int w=1;w<=n;w++)
{
for(int i=1;i<=n+1;i++) cin>>A[i];
for(int i=1;i<=n;i++)
{
if(fabs(A[i])<eps) continue;
if(fabs(F[i][i])<eps)
{
for(int j=i;j<=n+1;j++) F[i][j]=A[j]/A[i];
F[i][i]=1;break;
}
else for(int j=i+1;j<=n+1;j++) A[j]-=F[i][j]*A[i];
}
}
for(int i=1;i<=n;i++)
if(fabs(F[i][i])<eps)
return puts("No Solution"),0;
for(int i=n;i>=1;i--)
for(int j=i-1;j>=1;j--)
F[j][n+1]-=F[j][i]*F[i][n+1],F[j][i]=0;
for(int i=1;i<=n;i++) printf("%.2f\n",F[i][n+1]);
}
线性基
【模板】线性基
很方便地求解一个数集中能够异或出来的数
ll n,x,a[60],ans;
int main()
{
for(cin>>n;n;n--)
{
cin>>x;
for(int i=50;i>=0;i--)
if(x&(1ll<<i))
{
if(!a[i]) {a[i]=x;break;}
else x^=a[i];
}
}
for(int i=50;i>=0;i--)
if((ans^a[i])>ans) ans^=a[i];
return cout<<ans<<endl,0;
}
裴蜀定理
【模板】裴蜀定理
问题:给定\(A_1,A_2...A_k\),求\(A_1x_1+A_2x_2+...A_kx_k\)所能取到的最小值
结论:所有数的绝对值的\(gcd\)
int g,x,n;
int main()
{
cin>>n;while(n--) cin>>x,g=__gcd(g,abs(x));
cout<<g<<endl;
}
FFT
【模板】多项式乘法(FFT)
\(O(nlogn)\)求解两多项式卷积
//FFT
for(int i=1;i<l;i++) r[i]=(r[i>>1]>>1)|((i&1)<<tt);
for(int i=0;i<l;i++) w[i].rl=cos(Pi/l*i),w[i].im=sin(Pi/l*i);
void FFT(Complex *P,int op)
{
for(int i=1;i<l;i++) if(r[i]>i) swap(P[i],P[r[i]]);
for(int i=1;i<l;i<<=1)
for(int j=0,p=i<<1;j<l;j+=p)
for(int k=0;k<i;k++)
{
Complex W=w[l/i*k];W.im*=op;
Complex X=P[j+k],Y=P[j+k+i]*W;
P[j+k]=X+Y;P[j+k+i]=X-Y;
}
}
for(int i=0;i<=m+n;i++) printf("%d ",(int)(A[i].rl/l+0.5));
拉格朗日插值
\(n-1\)次多项式带入\(x=k\)求值。
\]
int Lagrange()
{
cin>>n>>k;
for(int i=1;i<=n;i++) cin>>x[i]>>y[i];
for(int i=1;i<=n;i++)
{
int fz=1,fm=1;
for(int j=1;j<=n;j++)
if(i!=j) fz=1ll*fz*(k-x[j])%P,fm=1ll*fm*(x[i]-x[j])%P;
(ans+=1ll*y[i]*fz%P*ksm(fm,P-2)%P)%=P;
}
cout<<(ans+P)%P<<endl;
}
NTT
FWT
这两个一个用NTT代替,一个用高维前缀和。考到了话我直播CS!