今天下午LYP神犇给讲了一坨数论= =发现太神了直接没办法愉快玩耍了
整理一下里面所有内容的模板。
这是一个史前巨坑,这点时间肯定不够整理完的。。。以后慢慢填坑orz
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#define MAXN 1000000
using namespace std;
bool not_prime[MAXN];
int prime_number[MAXN];
int nu;
int factor[MAXN];//顺手记录一下最小质因子
void prime()//线筛
{
for (int i=2;i<=MAXN;i++)
{
if (!not_prime[i])
{
prime_number[++nu]=i;
factor[i]=i;
}
for (int j=1;j<=nu&&i*prime_number[j]<=MAXN;j++)
{
not_prime[i*prime_number[j]]=1;
f[i*prime_number[j]]=prime_number[j];
if (i%prime_number[j]==0) break;
}
}
}
int mu[MAXN]={0,1};
void mu()//莫比乌斯函数
{
for (int i=2;i<=MAXN;i++)
{
if (!not_prime[i])
mu[i]=-1;
for (int j=1;j<=nu&&i*prime_number[j]<=MAXN;j++)
{
not_prime[i*prime_number[j]]=1;
if (i%prime_number[j]==0)
{
mu[i*prime_number[j]]=0;
break;
}
else
{
mu[i*prime_number[j]]=-mu[i];
}
}
}
}
int phi[MAXN];
void phi()//欧拉函数
{
for (int i=2;i<=MAXN;i++)
{
if (!not_prime[i])
{
phi[i]=i-1;
}
for (int j=1;j<=nu&&i*prime_number[j]<=MAXN;j++)
{
not_prime[i*prime_number[j]]=1;
if (i%prime_number[j]==0)
{
phi[i*prime_number[j]]=phi[i]*prime_number[j];
break;
}
else
{
phi[i*prime_number[j]]=phi[i]*(prime_number[j]-1);
}
}
}
}
int inv[MAXN];
//以下连续逆元
void inv1()//费马小定理
{
for (int i=1;i<=MAXN;i++)
inv[i]=pow(i,phi[p]-1)%p;//p根据题目要求不同来定
}
void inv2()//基于线筛的求法
{
for (int i=2;i<=MAXN;i++)
{
int a=p/i;
int b=p%i;
inv[i]=(p-a)*inv[b]%p;
}
}
void inv3()//基于阶乘的求法
{
unsigned long long a[MAXN]={1},rev[MAXN];//这个阶乘略大0-0不想写高精度了unsigned凑活看着吧
for (int i=1;i<=MAXN;i++)
a[i]=a[i-1]*i;
rev[MAXN]=pow(a[MAXN],phi[p]-1)%p;
for (int i=MAXN-1;i>=0;i--)
rev[i]=rev[i+1]*(i+1)%p;
for (int i=1;i<=MAXN;i++)
inv[i]=a[i-1]*rev[i];
}
void inv4()//基于a*i+b=p的求法,此做法仅用于p为素数
{
inv[1]=1;
for (int i=2;i<=MAXN;i++)
{
int a=p/i;
int b=p%i;
inv[i]=(p-a)*inv[b]%p;
}
}
int power(int x,int n)//快速幂非递归的
{
int ret=1;
while (n)
{
if (n&1)
ans*=x;
n>>=1;
x*=x;
}
}
int powermod(int x,int n,int p)//x^n mod p模幂运算
{
if (!n) return 1;
int ret=powermod((x*x)%p,n>>1,p);
if (n&1) ret=(ret*x)%p;
return ret;
}
int gcd(int a,int b)//欧几里得算法最大公约数
{
if (!b) return a;
else return (b,a%b);
}
int lcm(int a,int b)//最小公倍数
{
return a*b/gcd(a,b);
}
int ex_gcd(int a,int b,int &x,int &y)//扩展欧几里得 ax+by=gcd(a,b)
{
int temp;
if (a%b==0)
{
x=0;
y=1;
return a;
}
int ret=ex_gcd(b,a%b,x,y);
temp=x;
x=y;
y=t-a/b*y;
return ret;
}
int CRT(int n,int a[],int mod[])//中国剩余定理
{
int temp=1,ret=0,t,y;
for (int i=1;i<=n;i++)
temp*=mod[i];
for (int i=1;i<=n;i++)
{
int w=temp/mod[i];
t=ex_gcd(m[i],w,t,y);
x=(x+y*w*a[i])%temp;
}
while (x<=0) x+=temp;
return x;
}//x=a[1] (mod m[1])
//x=a[2] (mod m[2])
//x=a[3] (mod m[3])
//...
int inv(int a,int p)//配合下文食用
{
int d,x,y;
d=exgcd(a,p,x,y);
if (d==1)
return (x+p)%p;
else
return -1;
}
int C(int n,int m,int p)//普通组合
{
if(n<m)return 0;
return fac[n]*inv(fac[m]*fac[n-m]%p,p)%p;
}
int _C(int n,int m,int p)//二逼组合
{
if(n<m)return 0;
int ans=1,res=1;
for(int i=1;i<=n;i++) ans=ans*i%p;
for(int i=1;i<=n-m;i++) res=res*i%p;
for(int i=1;i<=m;i++) res=res*i%p;
res=inv(res,p);
ans=(ans*res)%p;
return ans;
}
int Lucas(int n,int m,int p)//文艺组合 Lucas定理求组合数C(n,m)n=a*p+b m=a’*p+b’ C(n,m)=C(a,a’)*C(b,b’)
{
int ret=1;
while (n&&m)
{
int a=n%p,b=m%p;
if (a<b) return 0;
ret=(ret*C(a,b,p))%p;
n/=p;
m/=p;
}
return ret;
}
int BSGS(int a,int b,int p)//BSGS算法 a^x=b (mod p) 不想手写hash改成了map其实这套模板是跟ZKY学长学的我还不会BSGS
{
int m=sqrt(p)+.5,temp=inv(power(a,m,p),p),t=1;
map<int,int>hash;
hash[1]=0;
for (int i=1;i<m;i++)
{
t=t*a%p;
hash[t]=i;
}
for (int i=0;i<=m;i++)
{
if (hash.count(b))
return i*m+hash[b];
b=b*temp%p;
}
return -1;
}
int get_ans(int a,int b,int k)//枚举除法取值,Mobius反演常用 在这之前要维护一个特别的前缀和,内容根据题目而定
{
int last=0,ret=0;
a/=k;b/=k;
for (int i=1;i<=min(a,b);i=last+1)
{
last=min(a/(a/i),b/(b/i));
ret+=(prev[last]-prev[i-1])*(a/i)*(b/i);
}
}
/*int main()
{
}*/
应该没错吧QUQ如果一旦有任何沙茶错误请诸位神犇立刻给我指出QUQ毕竟我是如此的沙茶