20190815 UPD
删掉了一些自己没有使用过的模板,增加了部分最近整理的模板
————————————————————————————————————————————————————————
数论
1.拉格朗日插值模板
2018牛客多校第一场F Sum of Maximum
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=1e3+5;
const int mod=1e9+7;
int n,res;
ll a[maxn],b[maxn],ans;
namespace polysum {
#define rep(i,a,n) for (int i=a;i<n;i++)
#define per(i,a,n) for (int i=n-1;i>=a;i--)
const int D=2010;
ll a[D],f[D],g[D],p[D],p1[D],p2[D],b[D],h[D][2],C[D];
ll powmod(ll a,ll b){ll res=1;a%=mod;assert(b>=0);for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;}
ll calcn(int d,ll *a,ll n) { // a[0].. a[d] a[n]
if (n<=d) return a[n];
p1[0]=p2[0]=1;
rep(i,0,d+1) {
ll t=(n-i+mod)%mod;
p1[i+1]=p1[i]*t%mod;
}
rep(i,0,d+1) {
ll t=(n-d+i+mod)%mod;
p2[i+1]=p2[i]*t%mod;
}
ll ans=0;
rep(i,0,d+1) {
ll t=g[i]*g[d-i]%mod*p1[i]%mod*p2[d-i]%mod*a[i]%mod;
if ((d-i)&1) ans=(ans-t+mod)%mod;
else ans=(ans+t)%mod;
}
return ans;
}
void init(int M) {
f[0]=f[1]=g[0]=g[1]=1;
rep(i,2,M+5) f[i]=f[i-1]*i%mod;
g[M+4]=powmod(f[M+4],mod-2);
per(i,1,M+4) g[i]=g[i+1]*(i+1)%mod;
}
/*
猜测:
a[0].. a[m] 指参数 m和a ,即取值的点对(j,a[j])
m=多项式的最高次
n-1求和部分为参数n的作用
即根据a[0]..a[m],求对应的多项表达式在[0,n-1]的和
*/
ll polysum(ll m,ll *a,ll n) { // a[0].. a[m] \sum_{i=0}^{n-1} a[i]
ll b[D];
for(int i=0;i<=m;i++) b[i]=a[i];
b[m+1]=calcn(m,b,m+1);
rep(i,1,m+2) b[i]=(b[i-1]+b[i])%mod;
return calcn(m+1,b,n-1);
}
/*
猜测:
求和部分比polysum多乘了一个指数为i的式子
*/
ll qpolysum(ll R,ll n,ll *a,ll m) { // a[0].. a[m] \sum_{i=0}^{n-1} a[i]*R^i
if (R==1) return polysum(n,a,m);
a[m+1]=calcn(m,a,m+1);
ll r=powmod(R,mod-2),p3=0,p4=0,c,ans;
h[0][0]=0;h[0][1]=1;
rep(i,1,m+2) {
h[i][0]=(h[i-1][0]+a[i-1])*r%mod;
h[i][1]=h[i-1][1]*r%mod;
}
rep(i,0,m+2) {
ll t=g[i]*g[m+1-i]%mod;
if (i&1) p3=((p3-h[i][0]*t)%mod+mod)%mod,p4=((p4-h[i][1]*t)%mod+mod)%mod;
else p3=(p3+h[i][0]*t)%mod,p4=(p4+h[i][1]*t)%mod;
}
c=powmod(p4,mod-2)*(mod-p3)%mod;
rep(i,0,m+2) h[i][0]=(h[i][0]+h[i][1]*c)%mod;
rep(i,0,m+2) C[i]=h[i][0];
ans=(calcn(m,C,n)*powmod(R,n)-c)%mod;
if (ans<0) ans+=mod;
return ans;
}
} // polysum::init();
inline int qpow(int a,int b){
int res=1;
while(b){
if(b&1) res=res*1ll*a%mod;
b>>=1;
a=a*1ll*a%mod;
}
return res;
}
int main(){
// freopen("in.txt","r",stdin);
polysum::init(1005);
while(~scanf("%d",&n)){
for(int i=1;i<=n;i++){
scanf("%lld",&a[i]);
}
sort(a+1,a+1+n);
res=1,ans=0;
for(int i=1;i<=n;i++){
if(a[i]==a[i-1]){
res=res*1ll*a[i]%mod;
continue;
}
b[0]=0;
for(int j=1;j<=n-i+1;j++){//因为多项式括号内作了差,所以当前最高次为n-i+1
b[j]=j*1ll*((qpow(j,n-i+1)-qpow(j-1,n-i+1))%mod)%mod;//(j,b[j])
if(b[j]<0) b[j]+=mod;
}
ans=(ans+res*((polysum::polysum(n-i+1,b,a[i]+1)-polysum::polysum(n-i+1,b,a[i-1]+1))%mod))%mod;
res=res*1ll*a[i]%mod;
}
if(ans<0) ans+=mod;
cout<<ans<<"\n";
}
return 0;
}
2.二次剩余模奇质数
(x^2==n(mod p))
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define random(a,b) (rand()%(b-a+1)+a)
ll quick_mod(ll a, ll b, ll c) { ll ans = 1; while (b) { if (b % 2 == 1)ans = (ans*a) % c; b /= 2; a = (a*a) % c; }return ans; }
ll p;
ll w;//二次域的D值
bool ok;//是否有解
struct QuadraticField//二次域
{
ll x, y;
QuadraticField operator*(QuadraticField T)//二次域乘法重载
{
QuadraticField ans;
ans.x = (this->x*T.x%p + this->y*T.y%p*w%p) % p;
ans.y = (this->x*T.y%p + this->y*T.x%p) % p;
return ans;
}
QuadraticField operator^(ll b)//二次域快速幂
{
QuadraticField ans;
QuadraticField a = *this;
ans.x = 1;
ans.y = 0;
while (b)
{
if (b & 1)
{
ans = ans*a;
b--;
}
b /= 2;
a = a*a;
}
return ans;
}
};
ll Legender(ll a)//求勒让德符号
{
ll ans=quick_mod(a, (p - 1) / 2, p);
if (ans + 1 == p)//如果ans的值为-1,%p之后会变成p-1。
return -1;
else
return ans;
}
ll Getw(ll n, ll a)//根据随机出来a的值确定对应w的值
{
return ((a*a - n) % p + p) % p;//防爆处理
}
ll Solve(ll n)
{
ll a;
if (p == 2)//当p为2的时候,n只会是0或1,然后0和1就是对应的解
return n;
if (Legender(n) == -1)//无解
ok = false;
while (1)//随机a的值直到有解
{
a = random(0, p - 1);
w = Getw(n, a);
if (Legender(w) == -1)
break;
}
QuadraticField ans,res;
res.x = a;
res.y = 1;//res的值就是a+根号w
ans = res ^ ((p + 1) / 2);
return ans.x;
}
int main()
{
ll n,ans1,ans2;
while (scanf("%lld%lld",&n,&p)!=EOF)
{
ok = true;
n %= p;
ans1 = Solve(n);
ans2 = p - ans1;//一组解的和是p
if (!ok)
{
printf("No root\n");
continue;
}
if (ans1 == ans2)
printf("%lld\n", ans1);
else
printf("%lld %lld\n", ans1, ans2);
}
}
3.BSGS求离散对数
unordered_map<ll,ll> ump;
ll solve(ll a,ll b,ll c){
a%=c,b%=c;
if(b==1) return 0;//若b是1则x为0
if(a==0){//需要特判a是c的倍数,0^1=0
if(b==0) return 1;
else return -1;
}
ump.clear();
ll m=ceil(sqrt(c+0.5)),tt=b,tt1=qpow(a,m,c),tt2=tt1;//qpow为快速幂,(a^m)%c
for(ll i=0;i<=m;i++,tt=tt*a%c){
ump[tt%c]=i;
}
for(ll i=1;i<=m;i++,tt2=tt2*tt1%c){
if(ump.find(tt2)!=ump.end()) return i*m-ump[tt2];
}
return -1;
}
4.扩展BSGS
unordered_map<ll,ll> ump;
ll ex_bsgs(ll a,ll b,ll c){
//a^x=b(%c) ,gcd(a,c)!=1,最后形如a1*a^(x-k)=b'(%c'),gcd(a,c')=1
a%=c,b%=c;
if(b==1) return 0;//若b是1则x为0
if(a==0){//需要特判a是c的倍数,0^1=0
if(b==0) return 1;
else return -1;
}
ll k=0,a1=1,d;
while(1){
d=__gcd(a,c);
if(d==1) break;
if(b%d!=0) return -1;
b/=d,c/=d;
k++;
a1=a1*(a/d)%c;
if(b==a1)return k;
}
ump.clear();
ll m=ceil(sqrt(c+0.5)),tt=b,tt1=qpow(a,m,c),tt2=tt1*a1%c;//qpow为快速幂,(a^m)%c
for(ll i=0;i<=m;i++,tt=tt*a%c){
ump[tt]=i;
}
for(ll i=1;i<=m;i++,tt2=tt2*tt1%c){
if(ump.find(tt2)!=ump.end()) return i*m-ump[tt2]+k;
}
return -1;
}
5.扩展欧几里得
求解ax=b(mod n)或ax+ny=b方程,若d=gcd(a,n),b%d==0时,x在1到a*b内有d个解,假设求出来的最小解为x0,则有xi=x0+(i-1)*(n/gc),其中i从0到d-1。
void exgcd(ll a,ll b,ll &x,ll &y){
//ax=gcd(a,b)(mod b),即ax+by=gcd(a,b),求出的解是|x|+|y|最小的
if(b==0) x=1,y=0;
else{
exgcd(b,a%b,x,y);
ll t=x;
x=y,y=t-a/b*y;
}
}
ll solve(ll a,ll b,ll c){//解ax+by=c,即ax=c(mod b)问题,求出最小的符合条件的x的非负解
ll gc=__gcd(a,b);
if(c%gc!=0) return -1;//无解
ll x,y,k=c/gc,tt=b/gc;
exgcd(a,b,x,y);
return (k*x%tt+tt)%tt;//保证答案在0~tt-1,即最小的非负解
}
6.中国剩余定理
/*
n个同余方程
求解ans%m[i]=b[i]
其中m[i]两两互质
通过exgcd求M/m[i]的逆元
解为\sum_{i=1}^{n}(b[i]*(M/m[i]*x))%M
其中M/m[i]*x=1(%m[i])
*/
int n;
ll m[15],b[15];
void exgcd(ll a,ll b,ll &x,ll &y){
//ax=gcd(a,b)(mod b),即ax+by=gcd(a,b),求出的解是|x|+|y|最小的
if(b==0) x=1,y=0;
else{
exgcd(b,a%b,x,y);
ll t=x;
x=y,y=t-a/b*y;
}
}
ll crt(int n,ll m[],ll b[]){
ll M=1,tmp,res=0,x,y;
for(int i=1;i<=n;i++) M*=m[i];
for(int i=1;i<=n;i++){
tmp=M/m[i];
exgcd(tmp,m[i],x,y);
res=(res+b[i]*(tmp*x%M))%M;
}
return res;
}
7.扩展中国剩余定理
洛谷P4777 【模板】扩展中国剩余定理(EXCRT)
保证答案在long long 范围内,但是中间过程可能溢出,因此需要__int128
int n;
ll m[maxn],b[maxn];
void exgcd(__int128 a,__int128 b,__int128 &x,__int128 &y){
//ax=gcd(a,b)(mod b),即ax+by=gcd(a,b),求出的解是|x|+|y|最小的
if(b==0) x=1,y=0;
else{
exgcd(b,a%b,x,y);
ll t=x;
x=y,y=t-a/b*y;
}
}
__int128 solve(__int128 a,__int128 b,__int128 c){//解ax+by=c,即ax=c(mod b)问题,求出最小的符合条件的x的非负解
__int128 gc=__gcd(a,b);
if(c%gc!=0) return -1;//无解
__int128 x,y,k=c/gc,tt=b/gc;
exgcd(a,b,x,y);
return (k*x%tt+tt)%tt;//保证答案在0~tt-1,即最小的非负解
}
ll ex_crt(int n,ll m[],ll b[]){
__int128 M=m[1],B=b[1],x;
for(int i=2;i<=n;i++){
x=solve(M,m[i],b[i]-B);
B+=x*M;
M=M/__gcd(M,(__int128)m[i])*m[i];
}
if(B<0) B+=M;
return (ll)B;
}
8.杜教筛
(51Nod 1244 莫比乌斯函数之和)
const int maxn=2000000;
int prime[150000],num;
int vst[maxn+5],miu[maxn+5];
inline void Pre(){
miu[1]=1;
for (int i=2;i<=maxn;i++){
if (!vst[i]) prime[++num]=i,miu[i]=-1;
for (int j=1;j<=num && (ll)i*prime[j]<=maxn;j++){
vst[i*prime[j]]=1;
if (i%prime[j]==0){
miu[i*prime[j]]=0;
break;
}
miu[i*prime[j]]=miu[i]*miu[prime[j]];
}
}
for (int i=1;i<=maxn;i++) miu[i]+=miu[i-1];
}
unordered_map<ll,int> S;
inline int Sum(ll n){
if (n<=maxn) return miu[n];
if (S.find(n)!=S.end()) return S[n];
int tem=1; ll l,r;
for (l=2;l*l<=n;l++) tem-=Sum(n/l);
for (ll t=n/l;l<=n;l=r+1,t--){
r=n/t;
tem-=(r-l+1)*Sum(t);
}
return S[n]=tem;
}
9.欧拉函数
求单个时间复杂度O(sqrt(n))
ll getphi(ll x){
ll rea=x;
for(ll i=2;i*i<=x;i++){
if(x%i==0){
rea-=rea/i;
while(x%i==0)
x/=i;
}
}
if(x>1) rea-=rea/x;
return rea;
}
O(n)线性筛打表
ll phi[N];
int pri[N],tot;
bool vis[N];
void getphi()
{
phi[1]=1;
for (int i=2;i<N;i++){
if(!vis[i]) pri[++tot]=i,phi[i]=i-1;
for (int j=1;j<=tot;j++){
if(1LL*i*pri[j]>N)break;
vis[1LL*i*pri[j]]=1;
if(i%pri[j]==0){
phi[1LL*i*pri[j]]=phi[i]*pri[j];break;
}
else
phi[1LL*i*pri[j]]=phi[i]*(pri[j]-1);
}
}
}
10.miller_rabbin大素数判定
(似乎常数较小,写https://nanti.jisuanke.com/t/25985这道题时,同样的写法,用另一个板子2000ms左右,用这个18ms...PS:后来听说是因为另一个板子计算时溢出了导致的超时)
UPD:据黄老板说这个快速乘有bug,待更新
#define ll long long
#define ldb long double
const int S=8;
ll mul(ll a,ll b,ll p)//O(1)快速乘
{
ll tmp=(a*b-(ll)((ldb)a/p*b+1e-8)*p);
return tmp<0?tmp+p:tmp;
}
ll quick_pow(ll a,ll b,ll p)
{
ll ans=1; a%=p;
for (ll i=b; i; i>>=1,a=mul(a,a,p))
if (i&1) ans=mul(ans,a,p);
return ans;
}
bool check(ll a,ll n,ll r,ll s)
{
ll ans=quick_pow(a,r,n),p=ans;
for(int i=1; i<=s; i++)
{
ans=mul(ans,ans,n);
if(ans==1&&p!=1&&p!=n-1)return 1;
p=ans;
}
if(ans!=1)return 1;
else return 0;
}
bool miller_rabbin(ll x)
{
if(x<=1) return 0;
if(x==2) return 1;
if(x%2==0) return 0;
ll r=x-1,s=0;
while(r%2==0) r/=2,s++;
for(int i=0; i<S; i++)
if(check(rand()%(x-1)+1,x,r,s))
return 0;
return 1;
}
11.线性筛
void init(){
notprime[1]=1;
For(i,2,N-1){
if(!notprime[i]){
prime[cnt++]=i;
}
for(int j=0;j<cnt&&i*prime[j]<N;j++){
notprime[i*prime[j]]=1;
if(i%prime[j]==0){
break;
}
}
}
}
筛后分解质因数
int pcnt;
ll p[105],e[105];
void fenjie(ll x){
memset(p,0,pcnt*sizeof(ll));
memset(e,0,pcnt*sizeof(ll));
pcnt=0;
for(int i=0;i<cnt&&prime[i]<x;i++){//cnt为筛出来的素数个数
if(x%prime[i]==0){
p[pcnt]=prime[i];
while(x%prime[i]==0)
x/=prime[i],e[pcnt]++;
pcnt++;
}
}
if(x>1){
p[pcnt]=x,e[pcnt++]++;
}
}
筛法求莫比乌斯函数和欧拉函数
const int N = 1e6+10;
int phi[N],mu[N],fac[N];
void init()
{
for(int i=1;i<N;++i) fac[i]=i;
phi[1]=mu[1]=1;
for(int i=2;i<N;++i)
{
if(fac[i]==i)
for(int j=i<<1;j<N;j+=i)
fac[j]=i;
if(i/fac[i]%fac[i]) phi[i]=(fac[i]-1)*phi[i/fac[i]],mu[i]=-mu[i/fac[i]];
else phi[i]=fac[i]*phi[i/fac[i]],mu[i]=0;
}
}
互质性质:若gcd(n,a)=1,则有gcd(n+k*a,a)=1(k为正整数)
对于,有
若则
————————————————————————————————————————————————————————
计算几何及其他
1、判断点是否在直线上
2、求只有一个交点的两直线的交点
bool online(int x1,int y1,int x2,int y2,int x,int y){
//判断点(x,y)是否在直线(x1,y1),(x2,y2)上
return (y2-y)*(x1-x)==(y1-y)*(x2-x);
}
void solve(int x1,int y1,int x2,int y2,int x3,int y3,int x4,int y4){
//求不平行且不重合的两条直线(x1,y1),(x2,y2)和(x3,y3),(x4,y4)的交点
db x=((x1-x2)*(y3-y4)*x4-(x3-x4)*(y1-y2)*x2+(x1-x2)*(x4-x3)*(y4-y2))*1.0/((x3-x4)*(y2-y1)-(x1-x2)*(y4-y3));
db y=(-(y4-y3)*(x2-x1)*y2+(y2-y1)*(x4-x3)*y4+(y2-y1)*(y3-y4)*(x4-x2))*1.0/((y4-y3)*(x1-x2)-(y2-y1)*(x3-x4));
}
定长圆覆盖最多点问题,时间复杂度O(n^2logn)(计蒜客2018ACM/ICPC沈阳网络赛E,先二分R再套这个模板)
const int N = 305;
const double eps = 1e-9;
const double pi = acos(-1.0);
#define sqr(x) ((x) * (x))
struct point
{
double x,y;
void read()
{
scanf("%lf%lf",&x,&y);
}
void print()
{
printf("%lf%lf\n",x,y);
}
double friend dis(const point &a,const point &b)
{
return sqrt(sqr(a.x - b.x) + sqr(a.y - b.y));
}
} p[N];
struct alpha//极角排序
{
double v;
bool flag;
bool friend operator <(const alpha &a,const alpha &b)
{
return a.v < b.v;
}
}alp[N*N];
int solve(int n,double R)//半径为R的圆最多能覆盖多少个点
{
int MAX = 0;
for(int i = 0; i < n; i++)//输入也应为0到n-1
{
int t = 0;
for(int j = 0; j < n; j++)
{
if(i == j)
continue;
double theta,phi,D;
D = dis(p[i],p[j]);
if(D > 2.0 * R)
continue;
theta = atan2(p[j].y - p[i].y,p[j].x - p[i].x);
if(theta < 0)
theta += 2 * pi;
phi = acos(D / (2.0 * R));
alp[t].v = theta - phi + 2 * pi;
alp[t].flag = true;
alp[t + 1].v = theta + phi + 2 * pi;
alp[t + 1].flag = false;
t += 2;
}
sort(alp,alp + t);
int sum = 0;
for(int j = 0; j < t; j++)
{
if(alp[j].flag)
sum ++;
else
sum --;
if(sum > MAX)
MAX = sum;
}
}
return MAX+1;
}
最小圆覆盖,时间复杂度O(n)(半径最小的圆覆盖所有的点)
/*最小圆覆盖*/
/*给定n个点, 让求半径最小的圆将n个点全部包围,可以在圆上*/
#include <cstdio>
#include <iostream>
#include <cstring>
#include <cmath>
#include <cstdlib>
#include <algorithm>
#define EPS 1e-8
using namespace std;
const int maxn = 550;
struct point{
double x, y;
};
int sgn(double x)
{
if (fabs(x) < EPS)
return 0;
return x < 0 ? -1 : 1;
}
double get_distance(const point a, const point b)//两点之间的距离
{
return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
point get_circle_center(const point a, const point b, const point c)//得到三角形外接圆的圆心
{
point center;
double a1 = b.x - a.x;
double b1 = b.y - a.y;
double c1 = (a1 * a1 + b1 * b1) / 2.0;
double a2 = c.x - a.x;
double b2 = c.y - a.y;
double c2 = (a2 * a2 + b2 * b2) / 2.0;
double d = a1 * b2 - a2 * b1;
center.x = a.x + (c1 * b2 - c2 * b1) / d;
center.y = a.y + (a1 * c2 - a2 * c1) / d;
return center;
}
//p表示定点, n表示顶点的个数, c代表最小覆盖圆圆心, r是半径
void min_cover_circle(point *p, int n, point &c, double &r)//找最小覆盖圆(这里没有用全局变量p[], 因为是为了封装一个函数便于调用)
{
random_shuffle(p, p + n);//随机函数,使用了之后使程序更快点,也可以不用
c = p[0];
r = 0;
for (int i = 1; i < n; i++)
{
if (sgn(get_distance(p[i], c) - r) > 0)//如果p[i]在当前圆的外面, 那么以当前点为圆心开始找
{
c = p[i];//圆心为当前点
r = 0;//这时候这个圆只包括他自己.所以半径为0
for (int j = 0; j < i; j++)//找它之前的所有点
{
if (sgn(get_distance(p[j], c) - r) > 0)//如果之前的点有不满足的, 那么就是以这两点为直径的圆
{
c.x = (p[i].x + p[j].x) / 2.0;
c.y = (p[i].y + p[j].y) / 2.0;
r = get_distance(p[j], c);
for (int k = 0; k < j; k++)
{
if (sgn(get_distance(p[k], c) - r) > 0)//找新作出来的圆之前的点是否还有不满足的, 如果不满足一定就是三个点都在圆上了
{
c = get_circle_center(p[i], p[j], p[k]);
r = get_distance(p[i], c);
}
}
}
}
}
}
}
int main()
{
int n;
point p[maxn];
point c; double r;
while (~scanf("%d", &n) && n)
{
for (int i = 0; i < n; i++)
scanf("%lf %lf", &p[i].x, &p[i].y);
min_cover_circle(p, n, c, r);
printf("%.2lf %.2lf %.2lf\n", c.x, c.y, r);
}
return 0;
}
int128模板
#include <bits/stdc++.h>
using namespace std;
void scan(__int128 &x)//输入
{
x = 0;
int f = 1;
char ch;
if((ch = getchar()) == '-') f = -f;
else x = x*10 + ch-'0';
while((ch = getchar()) >= '0' && ch <='9')
x = x*10 + ch-'0';
x *= f;
}
void print(__int128 x)//输出
{
if(x < 0)
{
x = -x;
putchar('-');
}
if(x > 9) print(x/10);
putchar(x%10 + '0');
}
int main()
{
__int128 a, b;
scan(a); scan(b);
print(a + b);
puts("");
print(a*b);
return 0;
}
n的环排列中,不含有相邻[i,i+1](mod n)的排列数有
通项为a[n]=(n-3)*a[n-1]+(n-2)*(2*a[n-2]+a[n-3])
FFT
hdu 1402大数乘法
#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define For(i,a,b) for(int i=a;i<=b;i++)
#define INF 0x3f3f3f3f
#define db double
#define ldb long double
#define m_p make_pair
#define p_b push_back
#define fbo friend bool operator <
const int N =200010;//
const db PI=acos(-1.0);
int n,m,t;
struct Complex{
db x,y;
Complex(db _x=0.0,db _y=0.0){
x=_x,y=_y;
}
Complex operator -(const Complex &b) const{
return Complex(x-b.x,y-b.y);
}
Complex operator +(const Complex &b) const{
return Complex(x+b.x,y+b.y);
}
Complex operator *(const Complex &b) const{
return Complex(x*b.x-y*b.y,x*b.y+y*b.x);
}
};
void change(Complex y[],int len){
int i,j,k;
for(i=1,j=len/2;i<len-1;i++){
if(i<j) swap(y[i],y[j]);
k=len/2;
while(j>=k){
j-=k;k/=2;
}
if(j<k) j+=k;
}
}
void fft(Complex y[],int len,int on){
change(y,len);
for(int h=2;h<=len;h<<=1){
Complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
for(int j=0;j<len;j+=h){
Complex w(1,0);
for(int k=j;k<j+h/2;k++){
Complex u=y[k],t=w*y[k+h/2];
y[k]=u+t;
y[k+h/2]=u-t;
w=w*wn;
}
}
}
if(on==-1){
for(int i=0;i<len;i++){
y[i].x/=len;
}
}
}
Complex x1[N],x2[N];
char str1[N],str2[N];
int sum[N];
int main(){
// freopen("1.txt","r",stdin);
while(scanf("%s %s",str1,str2)==2){
int len1=strlen(str1),len2=strlen(str2),len=1;
while(len<len1*2||len<len2*2) len<<=1;
for(int i=0;i<len1;i++){
x1[i]=Complex(str1[len1-1-i]-'0',0);//
}
for(int i=len1;i<len;i++)
x1[i]=Complex(0,0);
for(int i=0;i<len2;i++)
x2[i]=Complex(str2[len2-1-i]-'0',0);
for(int i=len2;i<len;i++)
x2[i]=Complex(0,0);
fft(x1,len,1);
fft(x2,len,1);
for(int i=0;i<len;i++){
x1[i]=x1[i]*x2[i];
}
fft(x1,len,-1);
for(int i=0;i<len;i++)
sum[i]=(int)(x1[i].x+0.5);
for(int i=0;i<len;i++){
sum[i+1]+=sum[i]/10;
sum[i]%=10;
}
len=len1+len2-1;
while(sum[len]<=0&&len>0) len--;
for(int i=len;i>=0;i--){
printf("%c",sum[i]+'0');
}
printf("\n");
}
return 0;
}
矩阵快速幂
const int mod = 1e9+7;
class Matrix{
public:
ll num[10][10];int n,m;
friend Matrix operator *(Matrix a,Matrix b){
Matrix ans;
ans.n=a.n,ans.m=b.m;
memset(ans.num,0,sizeof(ans.num));
for(int i=0;i<ans.n;i++){
for(int k=0;k<a.m;k++){
for(int j=0;j<ans.m;j++){
ans.num[i][j]=(ans.num[i][j]+a.num[i][k]*b.num[k][j])%mod;
}
}
}
return ans;
}
};
Matrix E(int n){
Matrix a;
memset(a.num,0,sizeof(a.num));
a.n=a.m=n;
for(int i=0;i<n;i++){
a.num[i][i]=1;
}
return a;
}
Matrix ksm(Matrix a,int t){
Matrix c=E(a.m);
while(t){
if(t&1) c=c*a;
t>>=1;
a=a*a;
}
return c;
}
最小生成树
const int N = 100010;
const int M = 500005;
int n,m,fa[N];
struct node{
int u,v,w;
node(){}
node(int _u,int _v,int _w){
u=_u,v=_v,w=_w;
}
}e[M];
int cnt;
void addedge(int u,int v,int w){
e[cnt++]=node(u,v,w);
}
bool cmp(node a,node b){
return a.w<b.w;
}
void init(){
For(i,1,n) fa[i]=i;
}
int find(int x){
return fa[x]==x?x:fa[x]=find(fa[x]);
}
int Kruskal(){
init();
sort(e,e+cnt,cmp);
int num=0,ans=0;
For(i,0,cnt-1){
int u=e[i].u,v=e[i].v,w=e[i].w;
int f1=find(u),f2=find(v);
if(f1!=f2) fa[f1]=f2,ans+=w,num++;
if(num==n-1) break;
}
if(num<n-1) return -1;
return ans;
}
输入挂(比后面那个快很多,数据越多优势越明显
struct FastIO//输入挂
{
static const int S=200;
int wpos;
char wbuf[S];
FastIO():wpos(0){}
inline int xchar()
{
static char buf[S];
static int len=0,pos=0;
if(pos==len) pos=0,len=fread(buf,1,S,stdin);
if(pos==len) exit(0);
return buf[pos++];
}
inline int read()
{
int s=1,c=xchar(),x=0;
while(c<=32) c=xchar();
if(c=='-') s=-1,c=xchar();
for(;'0'<=c&&c<='9';c=xchar()) x=x*10+c-'0';
return x*s;
}
~FastIO()
{
if(wpos) fwrite(wbuf,1,wpos,stdout),wpos=0;
}
}io;
template<class T>
inline bool scan_d(T &ret){
char c;
int sgn;
if(c=getchar(),c==EOF) return 0;
while(c!='-'&&(c<'0'||c>'9')) c=getchar();
sgn=(c=='-')?-1:1;
ret=(c=='-')?0:(c-'0');
while(c=getchar(),c>='0'&&c<='9') ret=ret*10+c-'0';
ret*=sgn;
return 1;
}
polya既可以旋转又可以翻折
ll pow_mod(ll a,ll b)//快速幂求a^b
{
ll res=1;
while(b)
{
if(b&1)
res=(res*a)%mod;
a=(a*a)%mod;
b=b>>1;
}
return res;
}
ll polya(ll m,ll n)//polya定理
{
ll i,ans=0;
for(i=1;i<=n;i++)
ans+=pow_mod(m,__gcd(i,n));//旋转的情况,循环数gcd(i,n),其中i为顺时针旋转i格
if(n&1)ans+=n*pow_mod(m,n/2+1);//奇数的翻转情况。共n条对称轴,每条的循环数均为n/2+1
else ans+=n/2*pow_mod(m,n/2)+n/2*pow_mod(m,n/2+1);//偶数的翻转情况。对称轴共n条,n/2条通过2个点,其余n/2条通过两点之间的中心。
ans=ans%mod*pow_mod(2*n,mod-2)%mod;
return ans;
}
线性基
/*
bool add(long long x);
构造线性基,同时可以判断x是否可成为线性基中的元素
*/
ll d[65];
bool add(ll x){
for(int i=60;i>=0;i--){
if(x&(1ll<<i)){//使用1ll,即长整型的1,因为1<<60显然会溢出
if(!d[i]){
d[i]=x;return true;
}
else x^=d[i];
}
}
return false;
}
/*
void work();
形成与原来等价的线性基,且若d[i]不为0,则第i+1个1仅由d[i]提供
可用来求第k小的异或和
所有元素异或即为最大异或和
*/
void work(){
for(int i=1;i<=60;i++) for(int j=0;j<i;j++) if(d[i]&(1ll<<j)) d[i]^=d[j];
}
杨辉三角求组合数(不取模N至多60)奇偶性:如果,那么为奇数,否则为偶数
typedef long long ll;
const int N=60;
const int mod = 1e9+7;
ll c[N][N];
void C()
{
int i,j;
for(i=0;i<N;i++)
c[i][0]=c[i][i]=1;
for(i=1;i<N;i++)
for(j=1;j<=i;j++)
c[i][j]=(c[i-1][j]+c[i-1][j-1])%mod;
}
求组合数
ll mul(ll a,ll b,ll p)//O(1)快速乘
{
ll tmp=(a*b-(ll)((ldb)a/p*b+1e-8)*p);
return tmp<0?tmp+p:tmp;
}
ll mod_pow(ll a,ll b,ll p)
{
ll ans=1; a%=p;
for (ll i=b; i; i>>=1,a=mul(a,a,p))
if (i&1) ans=mul(ans,a,p);
return ans;
}
ll C(ll n, ll m, ll p){//C(n,m)%p
if(m > n) return 0;
ll ret = 1;
m = min(n - m, m);
for(int i = 1; i <= m; i ++){
ll a = (n + i - m) % p;
ll b = i % p;
ret = ret * (a * mod_pow(b, p - 2, p) % p) % p;
}
return ret;
}
Lucas+中国剩余定理求组合数
#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define For(i,a,b) for(int i=a;i<=b;i++)
#define INF 0x3f3f3f3f
const int N = 100;
#define ldb long double
int t,k;
ll n,m,pp[N];
ll mul(ll a,ll b,ll p)//O(1)快速乘
{
ll tmp=(a*b-(ll)((ldb)a/p*b+1e-8)*p);
return tmp<0?tmp+p:tmp;
}
ll mod_pow(ll a,ll b,ll p)
{
ll ans=1; a%=p;
for (ll i=b; i; i>>=1,a=mul(a,a,p))
if (i&1) ans=mul(ans,a,p);
return ans;
}
ll C(ll n, ll m, ll p){
if(m > n) return 0;
ll ret = 1;
m = min(n - m, m);
for(int i = 1; i <= m; i ++){
ll a = (n + i - m) % p;
ll b = i % p;
ret = ret * (a * mod_pow(b, p - 2, p) % p) % p;
}
return ret;
}
ll Lucas(ll n,ll m,ll p){
if(m == 0ll) return 1;
return C(n%p,m%p,p)*(Lucas(n/p,m/p,p)%p)%p;
}
ll ni[N],li[N],Mi[N];
ll P;
ll exgcd(ll a,ll b,ll &x,ll &y){
if(b==0ll){
x=1;y=0;return a;
}
ll d=exgcd(b,a%b,x,y);
ll t=x;x=y;y=t-a/b*y;
return d;
}
ll solve(ll a,ll n){
ll x,y;
exgcd(a,n,x,y);
return x;
}
ll ans;
int main(){
scanf("%d",&t);
while(t--){
P=1,ans=0;
scanf("%lld %lld %d",&n,&m,&k);
For(i,1,k){
scanf("%lld",&pp[i]);P*=pp[i];
ni[i]=Lucas(n,m,pp[i]);
}
For(i,1,k) li[i]=P/pp[i],Mi[i]=solve(li[i],pp[i]),ans=(ans+mul(mul(Mi[i],li[i],P),ni[i],P))%P;//P最大1e18,这里不用快速乘会溢出,wa哭了
printf("%lld\n",(ans%P+P)%P);
}
return 0;
}
超级卡特兰数:s1=1,s2=1
施罗德数,可通过超级卡特兰数计算:s'0=1,s'1=2,s'i(i>=1)=2*si
奇数个顶点的多边形顶点全连起来,形成的区域个数为:(n-1)*(n-2)*(n*n-3*n+12)/24
kmp的next数组用法
#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;
char ch[1000005];
int next[1000005];
int hash[1000005];
void getnext()
{
int len=strlen(ch);
int i=0,j=-1;
next[0]=-1;
while(i<len)
{
if(j==-1 || ch[i]==ch[j])
{
i++;
j++;
next[i]=j;
// cout << "next[" << i << "]" << " "<<next[i]<<endl; //不懂主函数中while循环 就可以打印出来看看有助于直观理解
hash[j]++; //hash表统计长度为j的子串数目,因为在kmp串匹配中长度相同的串所表示的串是相同的
}
else
j=next[j];
}
}
int main()
{
// freopen("1.txt","r",stdin);
// freopen("3.txt","w",stdout);
while(cin>>ch)
{
int len=strlen(ch);
memset(hash,0,sizeof(hash));
getnext();
int i=len;
int flag=0;
hash[next[len]]--;
while(next[i]) //通过回溯的方法判断整个模式串
{
if(hash[next[i]]) //查找向前回溯后的next是否存在一个长度
{
flag=1;
break;
}
else
{
next[i]=next[next[i]];
}
}
if(!flag)
cout<<"Just a legend"<<endl;
else
{
for(int j=0;j<next[i];j++)
cout<<ch[j];
cout<<endl;
}
}
return 0;
}