lp2634 国家集训队 聪聪可可

这一题本质上和那道「模板 点分治1」是一样的,直接一边统计一边取模即可。

#include<iostream>
#include<cstdio>
#define Fv(i,X) for(int i=h[X];i;i=e[i].nxt) 

typedef long long ll;

inline int Max(int A,int B){
	return A>B?A:B;
}
inline ll gcd(ll A,ll B){
	return B?gcd(B,A%B):A;
}

struct ee{
	int v;
	int w;
	int nxt;
}e[40005];
int h[20005],et=0;
inline void Eadd(int U,int V,int W){
	e[++et]=(ee){V,W,h[U]};
	h[U]=et; 
}
inline void add(int U,int V,int W){
	Eadd(U,V,W);
	Eadd(V,U,W);
}

int n,rt=0,s;
int vis[20005],sz[20005],mx[20005];
inline void dfs0(int X,int FA){
	sz[X]=1,mx[X]=0;
	Fv(i,X){
		if(e[i].v==FA||vis[e[i].v]){
			continue;
		}
		dfs0(e[i].v,X);
		sz[X]+=sz[e[i].v];
		mx[X]=Max(mx[X],sz[e[i].v]);
	}
	mx[X]=Max(mx[X],s-sz[X]);
	if(mx[X]<mx[rt]){
		rt=X;
	}
}
int exist[3],dis[20005],st[20005],st2[20005],tp=0,tp2=0;
ll ans=0;
inline void dfs2(int X,int FA){
	st[++tp]=dis[X];
	Fv(i,X){
		if(e[i].v==FA||vis[e[i].v]){
			continue;
		}
		dis[e[i].v]=dis[X]+e[i].w;
		dfs2(e[i].v,X);
	}
}
inline void calc(int X){
	tp2=0;
	Fv(i,X){
		if(vis[e[i].v]){
			continue;
		}
		tp=0;
		dis[e[i].v]=e[i].w;
		dfs2(e[i].v,X);
		for(int j=1;j<=tp;++j){
			ans+=exist[(3-st[j]%3)%3];
		}
		for(int j=1;j<=tp;++j){
			st2[++tp2]=st[j];
			++exist[st[j]%3];
		}
	}
	for(int i=0;i<3;++i){
		exist[i]=0;
	}
}

inline void dfs1(int X){
	exist[0]=vis[X]=1;
	calc(X);
	Fv(i,X){
		if(vis[e[i].v]){
			continue;
		}
		rt=0;
		s=sz[X];
		dfs0(e[i].v,X);
		dfs1(rt);
	}
}

void init(){
	scanf("%d",&n);
	int u,v,w;
	for(int i=1;i<n;++i){
		scanf("%d%d%d",&u,&v,&w);
		add(u,v,w);
	} 
	rt=0,mx[0]=n;
	dfs0(1,0);
	dfs1(1);
	ans=(ans<<1)+n;
	ll g=gcd(ans,n*n);
	printf("%lld/%lld\n",ans/g,(n*n)/g);
}

int main(){
	init();
	return 0;
}

lp3806 【模板】点分治1

点分治是一种常常用于处理树上路径统计问题的算法。
首先我们可以想到这一题的一种简单的分治做法:找到一个点,从它开始搜,求出它的每一种长度。然后搜它的每一个子树,找到这些子树的根,求出它们的每一种长度并统计答案。
下面的问题便是要如何统计答案。容易想到的一种方案是枚举以某个点为根的链长两两相加得到的长度再容斥掉重复经过某一条边的路径,也就是所有经过它的同一个子节点的链对。
然而仔细想想这种做法的复杂度是有问题的。在菊花图上,它甚至会达到\(n^2\)的复杂度。故而我们必须另外考虑其他做法。
仔细观察数据范围,我们发现询问个数不超过100,这启发我们考虑一种\(O(nmlogn)\)的做法。
当我们求出以某个节点为根且经过前\(i\)个子节点的所有链以后,我们可以直接标记这这些长度的存在性,然后搜索经过第\(i+1\)个子节点的所有链,并枚举判断询问的长度减去这个链的长度剩下的长度是否存在。
这就做完了。

另:这题的数据是真的水。我犯了两个错,一个是没重设为局部大小,一个是找重心写挂了,而且用的还是复杂度错的算法,居然还过了…

#include<iostream>
#include<cstdio>

#define Fv(i,X) for(int i=h[X];i;i=e[i].nxt) 

inline int Max(int A,int B){
	return A>B?A:B;
}
struct ee{
	int v;
	int w;
	int nxt;
}e[20005];
int h[10005],et=0;
inline void Eadd(int U,int V,int W){
	e[++et]=(ee){V,W,h[U]};
	h[U]=et;
}
inline void add(int U,int V,int W){
	Eadd(U,V,W);Eadd(V,U,W);
}
int ans[10005]; 

int n,m,rt,s;

int mx[10005],sz[10005],vis[10005];
inline void dfs0(int X,int FA){
	sz[X]=1,mx[X]=0;
	Fv(i,X){
		if(e[i].v==FA||vis[e[i].v]){
			continue;
		}
		dfs0(e[i].v,X);
		sz[X]+=sz[e[i].v];
		mx[X]=Max(mx[X],sz[e[i].v]);
	}
	mx[X]=Max(mx[X],s-sz[X]);
	if(mx[X]<mx[rt]){
		rt=X;
	}
}
int a[10005];
int dis[10005],tp=0;
inline void dfs2(int X,int FA){
	dis[++tp]=a[X];
	Fv(i,X){
		if(e[i].v==FA||vis[e[i].v]){
			continue;
		}
		a[e[i].v]=a[X]+e[i].w;
		dfs2(e[i].v,X);
	}
}
int exist[10000005],lst[10005],cnt=0,qry[10005];
inline void calc(int X){
	cnt=0;
	Fv(i,X){
		if(vis[e[i].v]){
			continue;
		}
		a[e[i].v]=e[i].w;
		tp=0;
		dfs2(e[i].v,X);
		for(int j=1;j<=m;++j){
			for(int k=1;k<=tp;++k){
				(qry[j]-dis[k]>=0)?ans[j]|=exist[qry[j]-dis[k]]:0; 
			}
		}
		for(int j=1;j<=tp;++j){
			lst[++cnt]=dis[j];
			exist[dis[j]]=1;
		}
	}
	for(int i=1;i<=cnt;++i){
		exist[lst[i]]=0;
	}
}

inline void dfs1(int X){
	exist[0]=vis[X]=1;
	calc(X);
	Fv(i,X){
		if(vis[e[i].v]){
			continue;
		}
		s=sz[X],rt=0;//记得重置根。 
		dfs0(e[i].v,X);
		dfs1(rt);
	}
}

void init(){
	scanf("%d%d",&n,&m);
	int u,v,w,x;
	for(int i=1;i<n;++i){
		scanf("%d%d%d",&u,&v,&w);
		add(u,v,w);
	}
	for(int i=1;i<=m;++i){
		scanf("%d",&x);
		qry[i]=x;
	}
	s=n,mx[0]=n,rt=0;
	dfs0(1,0);
	dfs1(rt);
	for(int i=1;i<=m;++i){
		puts(ans[i]?"AYE":"NAY");
	}
} 

int main(){
	init();
	return 0;
}

lp5175 数列

$$f_{i}=f_{i-1}+a_{i}^2$$
$$a_{i}^2=(xa_{i-1}+ya_{i-2})^2=(x^2a_{i-1}^2+2xya_{i-1}a_{i-2}+y^
2a_{i-2}^2)$$
$$a_{i}a_{i-1}=(xa_{i-1}+ya_{i-2})a_{i-1}=xa_{i-1}^2+ya_{i-1}a_{i-2}$$
故而我们就可以用矩阵加速递推来求这个函数的值。
我们弄出一个如下的基本矩阵:
$$ \left[\begin{matrix}a_{1}^2&a_{2}^2&a_{1}a_{2}&f_{2}\end{matrix}\right]$$
我们不妨将它看作一个四项式,则其中每一项的答案都可以从上述公式推得。故而可以得到这样一个矩阵:
$$ \left[\begin{matrix}0&y^2&0&y^2\\1&x^2&x&x^2\\0&2xy&y&2xy\\0&0&0&1\end{matrix}\right]$$
拿来矩阵快速幂一下即可。

#include<iostream>
#include<cstdio>


typedef long long ll;
const ll MOD=1000000007;

struct Mat{
	ll a[4][4];
	inline void init(){
		for(int i=0;i<4;++i){
			for(int j=0;j<4;++j){
				a[i][j]=0;
			}
		}
	}
	inline Mat operator*(const Mat &B)const{
		Mat C;
		C.init();
		for(int i=0;i<4;++i){
			for(int j=0;j<4;++j){
				for(int k=0;k<4;++k){
					C.a[i][j]=(C.a[i][j]+a[i][k]*B.a[k][j])%MOD;
				}
			}
		}
		return C;
	}
};
ll a1,a2,n,x,y;
void init(){
	scanf("%lld%lld%lld%lld%lld",&n,&a1,&a2,&x,&y);
	if(n==1){
		printf("%lld\n",a1*a1%MOD);
		return;
	}
	if(n==2){
		printf("%lld\n",(a2*a2+a1*a1)%MOD);
	}
	Mat BS,T;
	BS.init(),T.init();
	BS.a[0][0]=a1*a1%MOD,BS.a[0][1]=a2*a2%MOD,BS.a[0][2]=a1*a2%MOD,BS.a[0][3]=(BS.a[0][0]+BS.a[0][1])%MOD;
	T.a[0][0]=0,T.a[0][1]=y*y%MOD,T.a[0][2]=0,T.a[0][3]=y*y%MOD;
	T.a[1][0]=1,T.a[1][1]=x*x%MOD,T.a[1][2]=x,T.a[1][3]=x*x%MOD;
	T.a[2][0]=0,T.a[2][1]=2*x*y%MOD,T.a[2][2]=y,T.a[2][3]=2*x*y%MOD;
	T.a[3][0]=0,T.a[3][1]=0,T.a[3][2]=0,T.a[3][3]=1;
	n-=2;
	while(n){
		if(n&1){
			BS=BS*T;
		}
		T=T*T;
		n>>=1;
	}
	printf("%lld\n",BS.a[0][3]);
}

int main(){
	int T;
	scanf("%d",&T);
	while(T--){
		init();
	}
	return 0;
}

lp5236 【模板】静态仙人掌(圆方树)

按照我们之前的经验,仙人掌上问题往往可以通过圆方树转化为树上问题。
我们发现,最短路在树上是一种非常容易解决的问题。只需预处理到根节点的长度然后死命跑LCA就可以了。
但是在一般图上,最短路的实时处理就会变得很困难。
我们可以尝试通过给圆方树上的边赋上特别的边权来处理这个问题。

对于圆-圆边,赋边权为原边权,这是容易理解的。
对于方点到它的父亲圆点的边,赋边权为0,对于圆点到它的父亲方点的边,赋边权为这个圆点到这个方点的父亲圆点的最短距离。

这时候我们会遇到一个问题,就是Tarjan找环的时候找不到返祖边的边权。
解决方案是记录每一个到根节点在dfs树上的距离,然后当我们找到一个环一路找爸爸并统计长度即可。

然后是计算答案,我们发现,如果询问的两个点的最近公共祖先是一个方点,那么它们的答案不能用普通方法计算。
我一开始的想法是尝试直接用某个节点到父亲方点的边权直接计算答案,但这样会导致答案错误,原因是它无法正确地区分两个点位于圆环的同一侧还是不同侧的情况。
故而,我们需要保存原来的距离它的父亲方点的靠某一侧的距离,故而当lca是方点的时候特殊判断计算即可。

 #include<iostream>
#include<cstdio>
#define Fv(H,A,X) for(int A=H[X];A;A=e[A].nxt)

typedef long long ll;

inline ll Min(ll A,ll B){
    return A<B?A:B;
}
inline int Min(int A,int B){
	return A<B?A:B;
}
inline ll Abs(ll A){
	return A>0?A:-A;
}
inline void Swap(int &A,int &B){
    A^=B^=A^=B;
}

struct ee{
    int v;
    ll w;
    int nxt;
}e[200005];
int h0[10005],h[20005],et=0;
inline void Eadd(int *H,int U,int V,ll W){
    e[++et]=(ee){V,W,H[U]};
    H[U]=et;
}
inline void add(int *H,int U,int V,ll W){
    Eadd(H,U,V,W);
    Eadd(H,V,U,W);
}
int dfn[20005],lw[20005],cnt=0,dep[20005],fa[20005][30];
ll dis[20005],sz[20005];
int nm=0;
int st[20005],tp=0;
inline void dfs0(int X){
    dfn[X]=lw[X]=++cnt;
    Fv(h0,i,X){
    	if(e[i].v==fa[X][0]){
    		continue;
		}
        if(!dfn[e[i].v]){
            fa[e[i].v][0]=X;
            dis[e[i].v]=dis[X]+e[i].w;
            dfs0(e[i].v);
            lw[X]=Min(lw[X],lw[e[i].v]);
    	}else{
			lw[X]=Min(lw[X],dfn[e[i].v]);
			if(e[i].v!=fa[X][0]&&dfn[e[i].v]<dfn[X]){
                ++nm;
                add(h,nm,e[i].v,0);
                ll len=e[i].w;
                for(int j=X;j^e[i].v;j=fa[j][0]){
                	len+=dis[j]-dis[fa[j][0]];
                }
                sz[nm]=len;
                ll nw=e[i].w;
                for(int j=X;j^e[i].v;j=fa[j][0]){
                    add(h,nm,j,Min(nw,len-nw));
                    sz[j]=nw; 
                    nw+=dis[j]-dis[fa[j][0]];
                }
            }
        }
        if(lw[e[i].v]>dfn[X]){
        	add(h,X,e[i].v,e[i].w);
		} 
    }
}

inline void dfs1(int X,int FA){
    dep[X]=dep[FA]+1;
    fa[X][0]=FA;
    Fv(h,i,X){
        if(e[i].v!=FA){
            dis[e[i].v]=dis[X]+e[i].w;
            dfs1(e[i].v,X);
        }
    }
}

int n,m,q;

inline ll calc(int X,int Y){
    int XX=X,YY=Y,lca;
    while(dep[XX]<dep[YY]){
        Swap(XX,YY);
    }
    for(int i=20;~i;--i){
        if(dep[XX]-(1<<i)>=dep[YY]){
            XX=fa[XX][i];
        }
    }
    if(XX==YY){
        lca=XX;
    }else{
        for(int i=20;~i;--i){
            if(fa[XX][i]!=fa[YY][i]){
                XX=fa[XX][i];
                YY=fa[YY][i];
            }
        }
        lca=fa[XX][0];
    }
    ll RT=dis[X]+dis[Y]-(dis[lca]<<1);
    if(lca>n){
        ll P=dis[XX]-dis[lca],Q=dis[YY]-dis[lca];
        RT-=(P+Q);
        RT+=Min(sz[lca]-Abs(sz[XX]-sz[YY]),Abs(sz[XX]-sz[YY]));
    }
    return RT;
}



void init(){
    scanf("%d%d%d",&n,&m,&q);
    nm=n;
    int u,v;
	ll w;
    for(int i=1;i<=m;++i){
        scanf("%d%d%lld",&u,&v,&w);
        add(h0,u,v,w);
    }
    fa[1][0]=0;
    dfs0(1);
    for(int i=1;i<=nm;++i){
        dep[i]=0,dis[i]=0;
    }
    dfs1(1,0);
    for(int j=1;j<=20;++j){
        for(int i=1;i<=nm;++i){
            fa[i][j]=fa[fa[i][j-1]][j-1];
        }
    }
    int x,y;
    for(int i=1;i<=q;++i){
        scanf("%d%d",&x,&y);
        printf("%lld\n",calc(x,y));
    }
}
int main(){
    init();
    return 0;
}

lp3320 SDOI2015 寻宝游戏

我们找到了这样一个结论:
「一个树上点集构成的最小生成树中,若两点间有路径,则此两点的DFS序在树上相近。」
这个结论的证明是较为困难的,但是感性理解还是比较容易的。
有了这个结论,我们就可以解决这一题。
每当一个新的点x即将加入点集的时候,我们只需要计算点集中DFS序刚好比它小的点y和DFS序刚好比它大的点z,然后将答案加上如下的式子:
$$dis_{x,y}+dis_{x,z}-dis_{y,z}$$
删去的时候逆向操作即可。
维护点集可以使用set。(善用STL(大雾))计算路径长度可以LCA+容斥。

注意:记得开long long

#include<iostream>
#include<cstdio>
#include<set>

struct ee{
	int v;
	int w;
	int nxt;
}e[200005];
int h[100005],et=0;
inline void add(int U,int V,int W){
	e[++et]=(ee){V,W,h[U]};
	h[U]=et;
}
int fa[100005][30],dfn[100005],loc[100005],dep[100005],cnt=0;
long long val[100005];
//注意将两种深度分开。 
inline void dfs(int X,int FA){
	fa[X][0]=FA,dfn[X]=++cnt,loc[cnt]=X,dep[X]=dep[FA]+1;
	for(int i=h[X];i;i=e[i].nxt){
		if(e[i].v!=FA){
			val[e[i].v]=val[X]+e[i].w;
			dfs(e[i].v,X);
		}
	}
}

inline int lca(int X,int Y){
	if(dep[X]<dep[Y]){
		std::swap(X,Y);
	}
	for(int i=20;i>=0;--i){
		if(dep[X]-(1<<i)>=dep[Y]){
			X=fa[X][i];
		}
	}
	if(X==Y){
		return X;
	}
	for(int i=20;i>=0;--i){
		if(fa[X][i]!=fa[Y][i]){
			X=fa[X][i];
			Y=fa[Y][i];
		}
	}
	return fa[X][0];
}

inline long long dis(int X,int Y){
	return val[X]+val[Y]-(val[lca(X,Y)]<<1);
}

int n,m;
long long ans=0;
bool vis[100005];
std::set<int> s;
void init(){
	scanf("%d%d",&n,&m);
	int u,v,w;
	for(int i=1;i<n;++i){
		scanf("%d%d%d",&u,&v,&w);
		add(u,v,w);
		add(v,u,w);
	}
	dfs(1,0);
	for(int j=1;j<=20;++j){
		for(int i=1;i<=n;++i){
			fa[i][j]=fa[fa[i][j-1]][j-1];
		}
	}
	int x,y,z;
	long long d;
	std::set<int>::iterator it;
	for(int i=1;i<=m;++i){
		scanf("%d",&x);
		x=dfn[x];
		if(!vis[loc[x]]){
			s.insert(x);
		}
		y=loc[(it=s.lower_bound(x))==s.begin()?*--s.end():*--it];
		z=loc[(it=s.upper_bound(x))==s.end()?*s.begin():*it];
		//注意运算符优先级。 
		if(vis[loc[x]]){
			s.erase(x);
		}
		x=loc[x];
		d=dis(x,y)+dis(x,z)-dis(y,z);
		if(!vis[x]){
			vis[x]=1;
			ans+=d;
		}else{
			vis[x]=0;
			ans-=d;
		}
		//注意前后的x不是一个x 
		printf("%lld\n",ans);
	}
}

int main(){
	init();
	return 0;
}

lp4630 APIO2018 Duathlon 铁人两项

圆方树是一种在仙人掌图上常用的数据结构,但是这并不意味着圆方树只在仙人掌图上有用。事实上,在任何一张图上,我们都可以用相似的方法来构造一棵圆方树。
对于一张图,我们预处理出它的所有点双,然后对每一个点双建一个方点,其他处理方法和仙人掌上圆方树几乎相同。
这里的点双要如何预处理呢?我们考虑一个性质:当我们预处理出一张图的DFS树后,任何一条边属于且仅属于一个点双。
那么,我们就可以尝试找到一个点双中深度最浅的节点,然后由它构造出这个方点。
我们发现,当我们搜索完一个节点的子树后,如果子节点中的某一个节点的lw是这个节点,那么这个节点就是它所在的点双中深度最浅的节点。
这是因为,如果一个节点的lw节点是它的父节点,那么它和它的子树就没有任何返祖边能够到达比当前节点更浅的节点,也就说明如果当前节点在点双内,那么它一定是点双内最浅的节点。
有没有可能当前节点不在点双内呢?这是不可能的。如果有返祖边连向当前节点,那么当前节点显然不会成为割点;如果没有,那么当前节点所属的点双就一定只有两个点。

现在把目光投到这道题。
首先考虑点双的基本性质。我们发现,对于至少有三个点的点双,任取三个点,它们总是一对满足题意的三元组。
这是由于点双的定义,点双里没有割点,自然总是可以找到两点间的两条不相交的简单路径。
拓展这个结论,我们发现,如果一条简单路径经过一个点双,那么显然无论固定哪个点,这条路径都始终是一条简单路径。
故而,对于选定的起点和终点,我们可以分两类讨论。
倘若它们位于同一个点双,那么显然它们之间的满足题意的三元组个数恰好是这个点双的点数-2
对于不属于同一个点双的情况,它们之间有可能经过其他的点双,也有可能不经过。每经过一个点双,它们之间的满足题意的三元组个数就会加上这个点双的点的个数。
同时,答案还要加上起点和终点各自处在的点双的点数个数各自-1。
这要怎么统计呢?我们可以让每个方点的权值为点双中的点的数量,每个圆点的权值为-1,然后枚举点对统计路径和即可。
仔细一想觉得有点不对劲儿:这样的复杂度岂不是要n^2logn级别?妥妥地T啊。
正难则反,我们可以统计每个点对答案的贡献,也就是经过它的路径数乘上它本身的权值。而前者可以通过一个普通的树形DP求得。
这就做完了。

注意:
要注意每个点的子树外的节点的数量是相当于连通块大小减去它的子树大小,而非总节点数减去它的子树大小。故而要注意统计连通块大小。

#include<iostream>
#include<cstdio>
#define Fv(H,A,X) for(int A=H[X];A;A=e[A].nxt)

inline int Min(int A,int B){
	return A<B?A:B;
}

struct ee{
	int v;
	int nxt;
}e[1200005];
int h0[100005],h[200005],et=0;
inline void Eadd(int *H,int U,int V){
	e[++et]=(ee){V,H[U]};
	H[U]=et;
}
inline void add(int *H,int U,int V){
	Eadd(H,U,V);
	Eadd(H,V,U);
}

int n,m;
int dfn[100005],lw[100005],cnt=0,nm=0;
int st[100005],tp=0;
int val[200005],nwsz;
long long ans=0;
inline void dfs(int X){
	dfn[X]=lw[X]=++cnt;
	st[++tp]=X;
	++nwsz;
//	要注意每个点的子树外的节点的数量是相当于连通块大小减去它的子树大小,而非总节点数减去它的子树大小。故而要注意统计连通块大小。 
	Fv(h0,i,X){
		if(!dfn[e[i].v]){
			dfs(e[i].v);
			lw[X]=Min(lw[X],lw[e[i].v]);
			if(lw[e[i].v]==dfn[X]){
//				注意这里应当判定的是lw(v)=dfn(u) 
				val[++nm]=1;
				for(int j=0;j!=e[i].v;--tp){
					++val[nm];
					j=st[tp];
					add(h,nm,j);
				}
				add(h,X,nm);
			}
		}else{
			lw[X]=Min(lw[X],dfn[e[i].v]);
		}
	}
}

int vis[200005],sz[200005];

inline void dfs2(int X){
	vis[X]=1;
	sz[X]=(X<=n);
	Fv(h,i,X){
		if(!vis[e[i].v]){
			dfs2(e[i].v);
			ans+=2ll*val[X]*sz[X]*sz[e[i].v];
			sz[X]+=sz[e[i].v];
		}
	}
	ans+=2ll*val[X]*sz[X]*(nwsz-sz[X]);
}

void init(){
	scanf("%d%d",&n,&m);
	nm=n;
	for(int i=1;i<=n;++i){
		val[i]=-1;
	}
	int u,v;
	for(int i=1;i<=m;++i){
		scanf("%d%d",&u,&v);
		add(h0,u,v);
	}
	for(int i=1;i<=n;++i){
		if(!dfn[i]){
			nwsz=0;
			dfs(i);
			dfs2(i);
			--tp;
		}
	}
	printf("%lld\n",ans);
}

int main(){
	init();
	return 0;
}

lp4244 SHOI2008 仙人掌图 II

仙人掌图是每一个点双都只有一个环的连通图。
对于一个仙人掌图,它在许多方面都有树的性质,只是其中的有一些节点被换成了环。
我们不妨依然采用一种伪树形DP的方式来处理这棵树,我们发现,非环节点的转移都是非常容易达成了,问题在于环上的节点的转移。

首先,我们要找到这个环,这可以用Tarjan算法来完成。
具体来说,我们维护两个数组:dfn和lw,其中前者表示的是某一个节点的dfs序,后者表示的是某一个节点至多走一条返祖边或者父亲边能够到达的最小dfs序。
显然,在搜索中,如果下一个节点未被访问过,我们就可以在搜索完它以后将当前节点的lw与它的lw取较小值。
否则,当前节点通向它的边一定是一条返祖边或者父亲边,那么当前节点的lw值应当与它的dfs序取较小值。

对于一个环,我们如果将它视为一个节点,我们想一想要怎么从它的「孩子」们处转移呢?
我们发现,依次枚举它的每一个孩子,答案显然是:
$$max(f_i+f_j+dis_{i,j})$$
这个式子可以使用单调队列优化。所以我们每一次把一个环提出来瞎DP一下就好了。

#include<iostream>
#include<cstdio>
#include<deque>
#define Fv(A,U) for(int A=h[U];A;A=e[A].nxt)

inline int Min(int A,int B){
	return A<B?A:B;
}
inline int Max(int A,int B){
	return A>B?A:B;
}

struct ee{
	int v;
	int nxt;
}e[200005];
int et=0,h[100005];
inline void Eadd(int U,int V){
	e[++et]=(ee){V,h[U]};
	h[U]=et;
}
inline void add(int U,int V){
	Eadd(U,V);
	Eadd(V,U);
}

int n,m;
int dfn[100005],lw[100005],dep[100005],fa[100005],cnt=0;
int ans=0,f[100005],a[100005];
std::deque<int> q;
inline void calc(int X,int Y){
	int tot=dep[Y]-dep[X]+1;
	int nw=tot;
	for(int i=Y;i!=X;i=fa[i]){
		a[nw--]=f[i];
	}
	a[nw]=f[X];
	for(int i=1;i<=tot;++i){
		a[tot+i]=a[i];
	}
	while(!q.empty()){
		q.pop_front();
	}
	q.push_back(1);
	for(int i=2;i<=(tot<<1);++i){
		while(!q.empty()&&i-q.front()>(tot>>1)){
			q.pop_front();
		}
		if(!q.empty()){
			ans=Max(ans,a[i]+a[q.front()]+i-q.front());
		}
		while(!q.empty()&&a[q.back()]-q.back()<a[i]-i){
			q.pop_back();
		}
		q.push_back(i);
	}
	for(int i=2;i<=tot;++i){
		f[X]=Max(f[X],a[i]+Min(i-1,tot-i+1));
	}
}
inline void dfs(int X,int FA){
	dfn[X]=lw[X]=++cnt;
	Fv(i,X){
		if(e[i].v!=FA){//这一句如果不加会无形降低很多点的lw值,从而使环的大小变成0。 
			if(!dfn[e[i].v]){
				fa[e[i].v]=X;
				dep[e[i].v]=dep[X]+1;
				dfs(e[i].v,X);
				lw[X]=Min(lw[X],lw[e[i].v]);
			}else{
				lw[X]=Min(lw[X],dfn[e[i].v]);
			}
			if(lw[e[i].v]>dfn[X]){
				ans=Max(ans,f[X]+f[e[i].v]+1);
				f[X]=Max(f[X],f[e[i].v]+1);
			}
		}
		
	}
	Fv(i,X){
		if(fa[e[i].v]!=X&&dfn[X]<dfn[e[i].v]){
			calc(X,e[i].v);
		}
	}
}

void init(){
	scanf("%d%d",&n,&m);
	int u,v,x;
	for(int i=1;i<=m;++i){
		u=0;
		scanf("%d",&x);
		for(int j=1;j<=x;++j){
			scanf("%d",&v);
			if(u){
				add(u,v);
			}
			u=v;
		}
	}
	dfs(1,0);
	printf("%d\n",ans);
}

int main(){
	init();
	return 0;
}

lp4382 八省联考2018 劈配

扫一眼题目,看到n个点和m个点匹配的问题,你就会有一种预感——是不是要跑二分图匹配?
是的,就是要跑二分图匹配。
第一问很显然,大力让一个导师能够被多个学生匹配,然后跑一个类似于匈牙利的操作。然后如果匹配满了,我们就考虑将这个导师匹配的其他学生尝试匹配同一个优先级的其他节点。
这个复杂度肯定是很松的。
那么关键在于第二问要怎么处理。
观察发现答案具有单调性。也就是说,如果一个选手在较高的优先度会感到沮丧,那么他在一个较低的优先度也会感到沮丧;
因此我们考虑二分答案来检验。
这样的复杂度大概是三方乘对数级的。时间比较紧张。
我们考虑一种优化策略:发现选手提高到某个优先度了以后,他就相当于顶替了本来在这个优先度的选手的位置。
故而,我们可以对于前缀的选手,储存到每个选手时的匹配状态。每次找到某个优先度,就把他的信息加入到对应的匹配状态中,然后跑上述操作。
这样可以优化到平方乘对数级,可能可以通过此题。
这一题的具体操作——也就是那个一点多匹配的操作并不是特别会,以后需要重新研究。
我好菜啊。

#include<iostream>
#include<cstdio>

int c,n,m;
int mx[205],a[205][205],val[205][205][15],exp[205];
int usd[205],usd1[205],usd2[205][205],vis[205];
int prusd[205][205],prusd1[205][205],prusd2[205][205][205];
inline bool dfs(int X,int V){
	int nw,nxt;
	for(int i=1;i<=val[X][V][0];++i){
		nw=val[X][V][i];
		if(vis[nw]){
			continue;
		}
		vis[nw]=1;
		if(usd[nw]<mx[nw]){
			usd1[X]=nw;
			usd2[nw][++usd[nw]]=X;
			return 1;
		}
		for(int j=1;j<=mx[nw];++j){
			nxt=usd2[nw][j];
			if(dfs(nxt,a[nxt][nw])){
				usd1[X]=nw;
				usd2[nw][j]=X;
				return 1;
			}
		}
	}
	return 0;
}
inline void calc1(){
	for(int i=1;i<=n;++i){
		for(int j=1;j<=m;++j){
			vis[j]=0;
		}
		for(int j=1;j<=m;++j){
			if(!val[i][j][0]){
				continue;
			}
			if(dfs(i,j)){
				break;
			}
		}
		for(int j=1;j<=m;++j){
			prusd[i][j]=usd[j];
		}
		for(int j=1;j<=i;++j){
			prusd1[i][j]=usd1[j];
		}
		for(int j=1;j<=m;++j){
			for(int k=1;k<=usd[j];++k){
				prusd2[i][j][k]=usd2[j][k];
			}
		}
	}
}
inline bool chck(int X,int Y){
	for(int i=1;i<=m;++i){
		usd[i]=prusd[Y-1][i];
	}
	for(int i=1;i<=Y-1;++i){
		usd1[i]=prusd1[Y-1][i];
	}
	for(int i=1;i<=m;++i){
		for(int j=1;j<=usd[i];++j){
			usd2[i][j]=prusd2[Y-1][i][j];
		}
	}
	for(int i=1;i<=m;++i){
		vis[i]=0;
	}
	for(int i=1;i<=exp[X];++i){
		if(!val[X][i][0]){
			continue;
		}
		if(dfs(X,i)){
			return 1;
		}
	}
	return 0;
}
void init(){
//	注意清空数组!!!
	scanf("%d%d",&n,&m);
	for(int i=1;i<=m;++i){
		scanf("%d",&mx[i]);
	}
	for(int i=1;i<=200;++i){
		usd1[i]=0;
		usd[i]=0;
		for(int j=1;j<=200;++j){
			for(int k=0;k<=10;++k){
				val[i][j][k]=0;
			}
		}
	}
	for(int i=1;i<=n;++i){
		for(int j=1;j<=m;++j){
			scanf("%d",&a[i][j]);
			if(a[i][j]){
				val[i][a[i][j]][++val[i][a[i][j]][0]]=j;
			}
		}
	}
	for(int i=1;i<=n;++i){
		a[i][0]=m+1;
	}
	for(int i=1;i<=n;++i){
		scanf("%d",&exp[i]);
	}
	calc1();
	for(int i=1;i<=n;++i){
		printf("%d ",a[i][usd1[i]]);
	}
	puts("");
	int l,r,mid,ans;
	for(int i=1;i<=n;++i){
		if(a[i][usd1[i]]<=exp[i]){
			printf("0 ");
			continue;
		}
		l=1,r=i-1,ans=0;
		while(l<=r){
			mid=(l+r)>>1;
			if(chck(i,mid)){
				ans=mid;
				l=mid+1;
			}else{
				r=mid-1;
			} 
		}
		printf("%d ",i-ans);
	}
	puts("");
}

int main(){
	int T;
	scanf("%d%d",&T,&c);
	while(T--){
		init();
	}
	return 0;
}

lp4245 【模板】任意模数NTT

我们首先考虑一个很暴力的玩法,直接找一个很大很大很大的模数,然后用int128上一个朴素NTT,再将结果对题目给的模数取模,似乎就可以了。
然而仔细考虑一下发现这个做法实在太暴力了(才不是因为找不到合适的模数呢。)我们考虑另一种方案。
现在假设我们已经用好几种模数求出最终的各种值了,那么原值是多少呢?很显然可以用中国剩余定理来求出。
我们计算了一下,发现只要用三个质数,值域就足以用中国剩余定理求出最终值了。
然而我们知道,模运算没有交换律。故而我们不能用常规的中国剩余定理合并,应该需要使用一种类似于扩展中国剩余定理的递推合并方法来合并它们的解。
这就完成了任意模数的NTT。
另外要注意数组大小…

#include<iostream>
#include<cstdio>

inline void Swap(int &A,int &B){
	A^=B^=A^=B;
}

inline long long mlt(long long A,long long X,long long P){
	long long BS=A,RT=0;
	while(X){
		if(X&1){
			RT+=BS;
			RT%=P;
		}
		BS+=BS;
		BS%=P;
		X>>=1;
	}
	return RT;
}

inline long long pw(long long A,long long X,long long P){
	long long BS=A,RT=1;
	while(X){
		if(X&1){
			RT*=BS;
			RT%=P;
		}
		BS*=BS;
		BS%=P;
		X>>=1;
	}
	return RT;
}

const long long MOD[3]={469762049,998244353,1004535809};
const long long MM=468937312667959297;
const long long g0=3,gi[3]={156587350,332748118,334845270};

int L=1,inv[3],R[1<<21|1];
inline void prpr(int LEN){
	int B=0;
	while(L<=LEN){
		L<<=1;
		++B;
	}
	for(int i=0;i<3;++i){
		inv[i]=MOD[i]-(MOD[i]-1)/L;
	}
	for(int i=0;i<L;++i){
		R[i]=R[i>>1]>>1|(i&1)<<(B-1);
	}
}

inline void NTT(int *A,int typ,int P){
	for(int i=0;i<L;++i){
		if(R[i]<i){
			Swap(A[i],A[R[i]]);
		}
	}
	int bs,nw,X,Y,M;
	for(int i=2;i<=L;i<<=1){
		M=i>>1;
		bs=pw(~typ?g0:gi[P],(MOD[P]-1)/i,MOD[P]);
		for(int j=0;j<L;j+=i){
			nw=1;
			for(int k=0;k<M;++k,nw=1ll*nw*bs%MOD[P]){
				X=A[j+k],Y=1ll*nw*A[j+k+M]%MOD[P];
				A[j+k]=(X+Y)%MOD[P];
				A[j+k+M]=(X-Y+MOD[P])%MOD[P];
			}
		}
	}
}
int C[1<<21|1],D[1<<21|1],ans[3][1<<21|1];
inline void FNTT(int *A,int *B,int P){
	for(int i=0;i<L;++i){
		C[i]=A[i],D[i]=B[i];
	}
	NTT(C,1,P);
	NTT(D,1,P);
	for(int i=0;i<L;++i){
		C[i]=1ll*C[i]*D[i]%MOD[P];
	}
	NTT(C,-1,P);
	for(int i=0;i<L;++i){
		ans[P][i]=(int)(1ll*(C[i]+MOD[P])*inv[P]%MOD[P]);
	}
}
int a[1<<21|1],b[1<<21|1],n,m,p;
long long t[3];
void init(){
	scanf("%d%d%d",&n,&m,&p);
	for(int i=0;i<=n;++i){
		scanf("%d",a+i);
	}
	for(int i=0;i<=m;++i){
		scanf("%d",b+i);
	}
	prpr(n+m);
	FNTT(a,b,0);FNTT(a,b,1);FNTT(a,b,2);
	t[0]=pw(MOD[1]%MOD[0],MOD[0]-2,MOD[0]);t[1]=pw(MOD[0]%MOD[1],MOD[1]-2,MOD[1]);t[2]=pw(MM%MOD[2],MOD[2]-2,MOD[2]);//分别求出要用到的三个逆元。 
	long long T1,T2;
	for(int i=0;i<=n+m;++i){
		T1=(mlt(1ll*ans[0][i]*MOD[1]%MM,t[0],MM)+mlt(1ll*ans[1][i]*MOD[0]%MM,t[1],MM))%MM;//直接用CRT合并一二式。 
		T2=((ans[2][i]-T1)%MOD[2]+MOD[2])%MOD[2]*t[2]%MOD[2];//用EXCRT将第三式合并之。 
		printf("%d ",(int)(((T2%p)*(MM%p)%p+T1%p)%p));//求得最终解。 
	}
}

int main(){
	init();
	return 0;
}

lp4777 【模板】扩展中国剩余定理(EXCRT)

$$ \left\{ \begin{matrix}x\equiv a_{1} \pmod{p_1} \\x\equiv a_{2}\pmod{p_2}\\\cdots\\x\equiv a_{n}\pmod{p_i}\end{matrix} \right.$$
对于这个方程,当所有模数两两互质的时候,可以使用中国剩余定理来求解。
然而,倘若并不满足这个条件,普通的中国剩余定理就未必能够奏效了。
我们尝试将中国剩余定理的成立条件进行推广,从而得到一种被称为「扩展中国剩余定理」的算法,以解决模数并不保证互质的情况下的线性同余方程组求解。
对于每个方程的解,它构成了一个解集。而我们需要找到的是这个解集中最小的且满足下一个方程的正整数解。
容易证明,一些线性同余方程组的解集一定是关于所有模数的最小公约数同余的一个数。
于是,我们令\(x_k\)为方程\(1\to k\)的最小正整数解:
$$\begin{equation}\begin{split}P,st:P&=&lcm({p_{i},i\in [1,k]})\\x_{k+1}&=&x_{k}+tP\\x_{k+1}&\equiv&a_{k+1}&\pmod{p_{k+1}}\\x_{k}+tP&\equiv&a_{k+1}&\pmod{p_{k+1}}\\tP&\equiv&a_{k+1}-x_{k}&\pmod{p_{k+1}}\end{split}\end{equation}$$
最后的那个式子可以转化形式然后上扩欧。
这就完成了一种递推求解不保证模数互质的线性同余方程组的方案。

#include<iostream>
#include<cstdio>

long long MOD=1;


inline long long mlt(long long A,long long X,long long P){
	long long BS=A,RT=0;
	while(X){
		if(X&1){
			RT+=BS;
			RT%=P;
		}
		BS+=BS;
		BS%=P;
		X>>=1;
	}
	return RT;
}
inline long long gcd(long long A,long long B){
	return B?gcd(B,A%B):A;
}
inline long long exgcd(long long A,long long B,long long &X,long long &Y,long long D=0){
	return B?(D=exgcd(B,A%B,Y,X),Y-=A/B*X,D):(X=1,Y=0,A);
}
long long n,a[100005],b[100005];
void init(){
	scanf("%lld",&n);
	for(int i=1;i<=n;++i){
		scanf("%lld%lld",&b[i],&a[i]);
	}
	long long X,Y,T,C,ans=a[1];
	MOD=b[1];
	for(int i=2;i<=n;++i){
		C=(a[i]-ans%b[i]+b[i])%b[i];
		X=exgcd(MOD,b[i],T,Y);
		if(C%X!=0){
			puts("F**KING HARD!");
			return;
		}
		T=mlt(T,C/X,b[i]/X);
		ans+=T*MOD;
		MOD/=X;
		MOD*=b[i];
		ans=(ans%MOD+MOD)%MOD;
		
	}
	printf("%lld\n",(ans+MOD)%MOD);
}

int main(){
	init();
	return 0;
}

lp70296 回忆京都

处理方法类似于NOIP2016,预处理完前缀和即可。
注意这里处理前缀和时需要判负数。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<queue>
#include<vector>
using namespace std;

const long long MOD=19260817;
long long C[1005][1005];

inline void prpr(){
	for(int i=0;i<=1000;++i){
		C[i][1]=i;
		C[i][i]=1;
	}
	for(int i=2;i<=1000;++i){
		for(int j=2;j<i;++j){
			C[i][j]=(C[i-1][j]+C[i-1][j-1])%MOD;
		}
	}
	for(int i=1;i<=1000;++i){
		for(int j=1;j<=1000;++j){
			C[i][j]+=(C[i-1][j]+C[i][j-1]-C[i-1][j-1]+MOD)%MOD;
			C[i][j]%=MOD;
		}
	}
}

void init(){
	int x,y;
	scanf("%d%d",&x,&y);
	printf("%lld\n",(C[y][x]+1)%MOD);
}
int main(){
	prpr();
	int T;
	scanf("%d",&T);
	while(T--){
		init();
	}
	return 0;
}

lp70295 整数校验器

细节题,要注意单独判断只有一个负号,形如“-”的情形。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<queue>
#include<vector>
using namespace std;

__int128 l,r;
char ch[200005];
int init(){
	cin>>ch+1;
	int n=strlen(ch+1);
	__int128 nw=0;
	bool bo=1;
	if(ch[1]=='-'){
		if(n==1){
			return 1;
		}
		for(int i=2;i<=n;++i){
			if(i==2&&ch[i]=='0'){
				return 1;
			}
			nw*=10;
			nw+=ch[i]-'0';
			if(-nw<l||-nw>r){
				return 2;
			}
		}
		return 0;
	}else{
		for(int i=1;i<=n;++i){
			if(i==1&&n>1&&ch[i]=='0'){
				return 1;
			}
			nw*=10;
			nw+=ch[i]-'0';
			if(nw<l||nw>r){
				return 2;
			}
		}
		return 0;
	}
}
int main(){
	long long L,R;
	int T;
	scanf("%lld%lld%d",&L,&R,&T);
	l=L,r=R;
	while(T--){
		printf("%d\n",init());
	}
	return 0;
}

lp3868 TJOI2009 猜数字

这是一道中国剩余定理(又名孙子定理、CRT)的模板题。
它最早来自于「孙子算经」,形式化的说,它是一种用于求解如下的线性同余方程组的解的算法。
$$ \left\{ \begin{matrix}x\equiv a_{1} \pmod{p_1} \\x\equiv a_{2}\pmod{p_2}\\\cdots\\x\equiv a_{n}\pmod{p_i}\end{matrix} \right.$$
我们令:
$$\begin{equation}\begin{split}P&=&\prod_{i=1}^{n}p_{i}\\t_{i},st: \frac{P}{p_{i}}t_{i}&\equiv &1\pmod{p_i}\\\exists x,st: x&\equiv &\sum_{i=1}^n \frac{P}{p_{i}}a_{i}t_{i}\\x&\in&\{x_0+kP,k\in Z\pmod{P}\}\end{split}\end{equation}$$
然而这是怎么证明的呢?我们有如下推导:
$$\begin{equation}\begin{split}\forall i,\forall j,st:i\neq j, a_{i}t_{i}\frac{P}{p_{i}}&\equiv&0&\pmod{p_i}\\a_{i}t_{i}\frac{P}{p_i}&\equiv&a_{i} &\pmod{p_i}\\\forall j,\sum_{i=1}^{n}a_{i}t_{i}\frac{P}{p_{i}}&\equiv&a_{j}t_{j}\frac{P}{p_{j}} \equiv a_{j}&\pmod{p_{j}}\\\exists x_0&\equiv& \sum_{i=1}^{n}a_{i}t_{i}\frac{P}{p_{i}}&\pmod{p_{j}}\\\end{split}\end{equation}$$
故而我们就可以求得上述方程的最小正整数解。
另外,看到这题的数据范围就可以猜到又是龟速乘登场的时刻了。
这一题还有个坑点就是读入的数可能是负数…

#include<iostream>
#include<cstdio>

long long MOD=1;


inline long long mlt(long long A,long long X){
	long long BS=A,RT=0;
	while(X){
		if(X&1){
			RT+=BS;
			RT%=MOD;
		}
		BS+=BS;
		BS%=MOD;
		X>>=1;
	}
	return RT;
}
inline long long exgcd(long long A,long long B,long long &X,long long &Y,long long D=0){
	return B?(D=exgcd(B,A%B,Y,X),Y-=A/B*X,D):(X=1,Y=0,A);
}
inline long long pw(long long A,long long X){
	long long BS=A,RT=1;
	while(X){
		if(X&1){
			RT=mlt(RT,BS)%MOD;
		}
		BS=mlt(BS,BS)%MOD;
		X>>=1;
	}
	return RT;
}
long long n,a[100005],b[100005];
void init(){
	scanf("%lld",&n);
	for(int i=1;i<=n;++i){
		scanf("%lld",&a[i]);
	}
	for(int i=1;i<=n;++i){
		scanf("%lld",&b[i]);
		a[i]=(a[i]%b[i]+b[i])%b[i];//注意负数!! 
		MOD*=b[i];
	}
	long long X,Y,T,ans=0;
	for(int i=1;i<=n;++i){
		T=MOD/b[i];
        exgcd(T,b[i],X,Y);
        X=(X%b[i]+b[i])%b[i];;
        ans=(ans+mlt(mlt(T,X),a[i]))%MOD;
	}
	printf("%lld\n",(ans+MOD)%MOD);
}

int main(){
	init();
	return 0;
}

lp2044 NOI2012 随机数生成器

我们不妨将原式递归展开,可以得到形如此的式子:
$$f_{n}=a^nf_{0}+\frac{(a^{n}-1)c}{a-1}$$
对于这个式子的值,我们看似可以上一个快速幂直接求得。
然而仔细观察上式的情况,发现如果出现了奇妙的模数时,\(a-1\)是不一定有逆元的!
这就很令人难受了。
我们找到原来的问题,我们不得不求下述式子的值:
$$f_{n}=a^nf_{0}+\sum_{i=0}^{n-1}a^ic=a^nf_{0}+(\sum_{i=0}^{n-1}a^i)c$$
求值的最大障碍是这个式子:
$$\sum_{i=0}^{n-1}a^i$$
将这个式子分奇偶讨论后递归求解。
我们令:
$$u_{x}=\sum_{i=0}^{x-1}a^i$$
若是偶数,则有:
$$u_x=\sum_{i=0}^{\frac{x-1}{2}}a^i+a^{\frac{x+1}{2}} (\sum_{i=0}^{\frac{x-1}{2}}a^i)$$
即,
$$u_x=(a^{\frac{x+1}{2}}+1)u_{\frac{x-1}{2}}$$
若是奇数,则只需将第一项单独计算。
这样就可以在对数复杂度内对式子求解了。
另外,考虑到模数的数据范围,我们需要使用一种被称为「龟速乘」,也就是「快速幂」在乘法意义上的等同的操作,使得不会出现乘法,也就不会爆范围。

#include<iostream>
#include<cstdio>

long long MOD,a,c,x0,n,g;
//必须注意,不能重载一个全部是标准型的运算符。 
inline long long mlt(long long A,long long X){
	long long BS=A,RT=0,CNT=X;
	while(CNT){
		if(CNT&1){
			RT+=BS;
			RT%=MOD;
		}
		CNT>>=1;
		BS+=BS;
		BS%=MOD;
	}
	return RT;
}
inline long long pw(long long A,long long X){
	long long BS=A,RT=1,CNT=X;
	while(CNT){
		if(CNT&1){
			RT=mlt(RT,BS)%MOD;
		}
		CNT>>=1;
		BS=mlt(BS,BS)%MOD;
	}
	return RT;
}
inline long long calc(long long X){
	if(X==2){
		return (a+1)%MOD;
	}
	if(X==1){
		return 1;
	}
	if(X&1){
		return (pw(a,X-1)+calc(X-1))%MOD;
	}
	return mlt((pw(a,(X+1)>>1)+1),calc(X>>1));
}
void init(){
	scanf("%lld%lld%lld%lld%lld%lld",&MOD,&a,&c,&x0,&n,&g);
	printf("%lld",((mlt(pw(a,n),x0)+mlt(calc(n),c))%MOD)%g);
}

int main(){
	init();
	return 0;
}

lp3803 【模板】多项式乘法(FFT)(NTT实现)

程序运行效率。(评测系统为duck.ac)

我们现在已经有一种快速完成多项式卷积的方案了。然而这种方案还存在两个缺陷。
首先,是浮点型运算带来的精度误差。这在多次乘积之后会不断扩大,最终到一个无法解决的地步。
第二,就是结构体运算带来的大常数。这也使得这道题的时间限制变得比较紧张。
于是,我们考虑另辟蹊径,使用某一种大模数来完成整数多项式卷积。
快速数论变换(Fast Number Transform)就是一种这样的算法。

回到FFT的本质,FFT是一种精确选取取值以后代值优化系数表达和点值表达相互转化的方案。我们尝试在模意义下找到一些满足与单位复根性质类似的数。
一种在数论上被称为原根的数可以满足这一族性质。
首先我们有原根的定义。
形式化的说,若\(g\)模\(p\)的阶等于\(\phi(p)\),则称\(g\)为\(p\)的一个原根。
具体来说,我们认为,若\(g^x \equiv 1\ (mod\ p) \)总是有解,那么我们称\(g\)是质数\(p\)的一个原根。
对于\(\forall g\in Z,p\in P\),总是有\(g^{p-1}\equiv 1\ (mod\ p) \)
那么,对于一个\(p=t2^{k}+1\)我们有,\(g^{\frac{p-1}{2^i}}\ (0\le i\le k)\)在模$latexp$意义下始终有正整数解。这一点和单位复根类似。
通过选取满足上述条件的质数,即\(p=t2^{k}+1\),我们还有这样的性质:
消去引理:
$$g^{\frac{dk(p-1)}{dn}}\equiv g^{\frac{k(p-1)}{n}}$$
折半引理:
$$(g^{\frac{(k+\frac{n}{2})(p-1)}{n}})^2\equiv (g^{\frac{k(p-1)}{n}}g^{\frac{p-1}{2}})^2\equiv (g^{\frac{k(p-1)}{n}}(p-1))^2\equiv (g^{\frac{k(p-1)}{n}})^2$$
求和引理:
$$\sum_{j=0}^{n-1}(g^{\frac{k(p-1)}{n}})^j\equiv \frac{(g^{\frac{k(p-1)}{n}})^n-1}{g^{\frac{k(p-1)}{n}}-1}\equiv \frac{g^{k(p-1)}-1}{g^{\frac{k(p-1)}{n}}-1}\equiv \frac{1^k-1}{g^{\frac{k(p-1)}{n}}}\equiv 0$$
折半引理的一个推论:
$$\sum_{i=0}^{n-1}g^{\frac{i(p-1)}{n}}\equiv 0$$
我们惊讶地发现原根具有单位复根具有的三个定理和一个推论。这说明,我们可以将快速傅里叶变换的操作以原根的形式拓展到模意义下。
我们令:
$$\omega_{n}^{k}=g^{\frac{k(p-1)}{n}}$$
那么我们就可以简单地将快速傅里叶变换中的单位复根全部替换为原根,这样就实现了快速数论变换。
对于这里的模数\(p\),我们可以选取998244353,它是\(2^{23}*119+1\)。

#include<iostream>
#include<cstdio>
#define Swap(A,B) (A^=B^=A^=B)
const long long P=998244353;
const long long g0=3,gi=332748118;
//*****非常重要:gi是332748118!!!!!!! 
int L=1,R[1<<21|1];
long long invL;
inline int pw(int A,int X){
	long long BS=A,RT=1;
	while(X){
		if(X&1){
			RT=RT*BS%P;
		}
		BS=BS*BS%P;
		X>>=1;
	}
	return RT;
}
inline void prpr(int LEN){
	int B=0;
	while(L<=LEN){
		L<<=1;
		++B;
	}
	invL=P-(P-1)/L;
	for(int i=0;i<L;++i){
		R[i]=R[i>>1]>>1|(i&1)<<(B-1);
	}
}
inline void FNTT(int *A,int typ){
	for(int i=0;i<L;++i){
		if(R[i]<i){
			Swap(A[R[i]],A[i]);
		}
	}
	int gn,g,X,Y,M;
	for(int i=2;i<=L;i<<=1){
		M=i>>1;
		gn=pw(~typ?g0:gi,(P-1)/i);
		for(int j=0;j<L;j+=i){
			g=1;
			for(int k=0;k<M;++k,g=1ll*g*gn%P){
				X=A[j+k],Y=1ll*g*A[M+j+k]%P;
				A[j+k]=(X+Y)%P,A[M+j+k]=(X-Y)%P;
			}
		}
	}
}
int n,m,a[1<<21|1],b[1<<21|1];
void init(){
	scanf("%d%d",&n,&m);
	for(int i=0;i<=n;++i){
		scanf("%d",a+i);
	}
	for(int i=0;i<=m;++i){
		scanf("%d",b+i);
	}
	prpr(n+m);
	FNTT(a,1);
	FNTT(b,1);
	for(int i=0;i<L;++i){
		a[i]=1ll*a[i]*b[i]%P;
	}
	FNTT(a,-1);
	for(int i=0;i<=n+m;++i){
		printf("%d ",1ll*(a[i]+P)*invL%P);
	}
}

int main(){
	init();
	return 0;
}

lp3812 【模板】线性基

线性基这个名字似乎是来自于它的数学意义。它的用途一般是在集合异或的题目中简化复杂度。
对于关于一个数集的线性基,它是一个集合。这个集合包含的数大概是数集中的数取对2取对数个。
它满足以下性质:它们中的每个数异或在一起可以得到原数集中任何数的异或和。
它是一个尽可能小的集合。
它是「数集中线性无关的最大子集」。(这一定义会在之后的深入学习中概括)
它的每个数的二进制最高位互不相同。

构造一个线性基可以用贪心的构造法。
对于一个即将插入线性基中的数\(x\),我们必须尽量让它能够以某一种方式被线性基中已有的数构造出来。
我们不妨令线性基中第\(i\)个数表示的是最高位为第\(i\)位的数。
如果线性基中能够构造出这个数,那么显然,从高到低枚举每一个数,每一次都将它异或去最高位的线性基所代表的数,那么最终一定会有两种情况之一:
一:这个数存在,则这个数按照上述方法异或后最终肯定会变成0。
二:这个数不存在,则这个数按照上述方法异或后最终肯定会遇到某一位线性基,使得这个线性基并不存在。那么这个数就可以放到给定线性基上。
上述两个情况是很显然的,因为根本不存在上述两种情况之外的情况。对于任意一种情况它的较高位一定可以被异或掉,而随之诞生的较低位也可以依次按需异或掉。

为什么这样构造出来的东西一定是满足性质一的呢?如果原数集中任意一个数都可以通过线性基里的某些数异或得到,那么任意一些数的异或和也一定可以通过这些数在线性基里的表示法的异或和得到。

回到这一题。这些数中的最大异或和一定是线性基中的最大异或和。
根据最大异或和的贪心构造法(即从高到低尝试异或每一项,如果能变得更大就异或),这样得到的结果必然是一个原数集能够构造出的数。这是有线性基的插入方式决定的。
按照上述的插入方式,一个新的数\(x\)只会在它的某一位作为最高位在线性基中不存在的时候连同它剩余的位被插入。故而任何一个较后面的位的修改都一定是在原数集中可以构造出来的修改。
所以这一题可以用线性基解决。

#include<iostream>
#include<cstdio>

long long bs[105];
int n;
void init(){
	scanf("%d",&n);
	long long x;
	for(int i=1;i<=n;++i){
		scanf("%lld",&x);
		for(int j=50;j>=0;--j){
			if(x&(1ll<<j)){
				if(!bs[j]){
					bs[j]=x;
				}
				x^=bs[j];
			}
		}
	}
	x=0;
	for(int i=50;i>=0;--i){
		if(x<(x^bs[i])){
			x^=bs[i];
		}
	}
	printf("%lld\n",x);
}

int main(){
	init();
	return 0;
}

lp3648 APIO2014 序列分割

首先证明序列分割的次序和最后的答案毫无关系。
这可以用乘法的基本性质证明。例如分成三段,我们有:
$$a(b+c)+bc=b(a+c)+ac$$

我们可以进行DP。
我们记\(f_{i,k}\)表示当前访问到第\(i\)个点,使用了\(k\)次切割的答案。
那么根据题目中贡献的定义,我们有:
$$f_{i,k}=max_{j=k}^{i-1}{f_{j,k-1}+sm_{j}(sm_{i}-sm_{j})}$$ 令\(k\le j_1(sm_i-sm_{j_1})\le f_{j_2,k-1}+sm_{j_2}(sm_i-sm_{j_2})\)$ 展开,移项,合并。 $$f_{j_1,k-1}+sm_{j_1}sm_i-sm_{j_1}^2\le f_{j_2,k-1}+sm_{j_2}sm_i-sm_{j_2}^2$$ $$f_{j_1,k-1}-f_{j_2,k-1}+sm_{j_2}^2-sm_{j_1}^2\le sm_{j_2}sm_i-sm_{j_1}*sm_i$$
$$\frac{f_{j_1,k-1}-f_{j_2,k-1}+sm_{j_2}^2-sm_{j_1}^2}{sm_{j_2}-sm_{j_1}}\le sm_i$$
这个式子显然是可以斜率优化的。
这样做的时空复杂度是\(O(nk)\)的。
考虑到这一题的空间限制比较紧张,我们可以用滚动数组把\(k\)这一维压掉。

#include<iostream>
#include<cstdio>

int q[100005];
long long f[100005][2],sm[100005];
int ans[100005][205],n,K;
inline double slope(int j1,int j2,int k){
	return (sm[j1]^sm[j2])?(double)(f[j1][(k&1)^1]-f[j2][(k&1)^1]+sm[j2]*sm[j2]-sm[j1]*sm[j1])/(double)(sm[j2]-sm[j1]):-1E18;
} 
void init(){
	scanf("%d%d",&n,&K);
	for(int i=1;i<=n;++i){
		scanf("%lld",&sm[i]);
		sm[i]+=sm[i-1];
	}
	for(int i=1;i<=K;++i){
		int l,r;
		l=r=1;
		for(int j=1;j<=n;++j){
			while(l<r&&slope(q[l],q[l+1],i)<=sm[j]){
				++l;
			}
			f[j][i&1]=f[q[l]][(i&1)^1]+sm[q[l]]*(sm[j]-sm[q[l]]);
			ans[j][i]=q[l];
			while(l<r&&slope(q[r-1],q[r],i)>slope(q[r],j,i)){
				--r;
			}
			q[++r]=j;
		}
	}
	int nw=n;
	printf("%lld\n",f[n][K&1]);
	for(int i=K;i>=1;--i){
		nw=ans[nw][i];
		printf("%d ",nw);
	}
}

int main(){
	init();
	return 0;
}

lp3803 【模板】多项式乘法(FFT)(非递归版本)

快速傅里叶变换(FFT)的一种较高效的实现

*由于快速傅里叶变换(FFT)常常被用在通信领域,故而FFT的实现效率成为了一个需要被考察的对象。

事实上,观察FFT的实现过程,我们发现它没有非递归不可的理由。

具体来说,我们发现,如果我们想要把FFT改造为递推模式,最大的障碍就是每一次合并操作。
在快速离散傅里叶变换中,每一次合并操作都会导致两个距离为\(\frac{n}{2}\)的点的值合并在一起。因此,如果将原序列传入,那么得到的就是一个乱序的队列。
我们考虑将序列以什么顺序传入。
我们发现,每一次这样的合并操作本质上都是将两个部分「翻转」了。
以0~7为例:
0 1 2 3 4 5 6 7
0 4 2 6 1 5 3 7
将它们写成二进制:
000 001 010 011 100 101 110 111
000 100 010 110 001 101 011 111
发现变换之后每个数的二进制位都倒转了。
故而我们可以通过倒转二进制位的方式来预处理原序列,使得FFT之后的序列顺序是正常的。
但是不知道是什么原因,这个代码跑出来的速度竟然也是4500ms…这就很让人挠头了。
于是向PinkRabbit巨神请教了几个经验:
1.不用中间数组,直接在原数组上操作。
2.调整赋值操作,使得修改尽量少。
3.将蝴蝶变换预处理为一个数组。
综合几点,就成功地卡到了2500ms
当然还是有待卡常,毕竟这还是一个比较危险的时间。

#include<iostream>
#include<cstdio>
#include<cmath>
#define PI 3.141592653589793238462

struct Complex{
	double r;
	double i;
	Complex(double ri=0.0,double ii=0.0):r(ri),i(ii){}
	inline Complex operator+(const Complex &B)const{
		return (Complex){r+B.r,i+B.i}; 
	}
	inline Complex operator-(const Complex &B)const{
		return (Complex){r-B.r,i-B.i};
	}
	inline Complex operator*(const Complex &B)const{
		return (Complex){r*B.r-i*B.i,r*B.i+i*B.r};
	}
}a[1<<21|1],b[1<<21|1],s[1<<21|1];
int r[1<<21|1];
inline void FFT(Complex *A,int N,int typ){
	//ButterflyTransform
	for(int i=0;i<N;++i){
		(r[i]>i)?std::swap(A[i],A[r[i]]),0:0;
	}
	Complex bs,nw,x,y;
	for(int i=2,mid;i<=N;i<<=1){
		bs=Complex(cos(2.0*PI/i),typ*sin(2.0*PI/i));
		mid=i>>1;
		for(int j=0;j<N;j+=i){
			nw=Complex(1.0,0.0);
			for(int k=0;k<mid;++k,nw=nw*bs){
				x=A[j+k],y=(A[mid+j+k]*nw);
				A[j+k]=x+y;
				A[mid+j+k]=x-y;
			}
		}
	}
}
int n,m,cnt=1,byt=0;;
void init(){
	scanf("%d%d",&n,&m);
	for(int i=0;i<=n;++i){
		scanf("%lf",&a[i].r);
	}
	for(int i=0;i<=m;++i){
		scanf("%lf",&b[i].r);
	}
	while(cnt<=n+m){
		cnt<<=1;
		++byt;
	}
	for(int i=0;i<cnt;++i){
		r[i]=(r[i>>1]>>1)|((i&1)<<(byt-1));
	}
	FFT(a,cnt,1);
	FFT(b,cnt,1);
	for(int i=0;i<cnt;++i){
		a[i]=a[i]*b[i];
	}
	FFT(a,cnt,-1);
	//注意这里typ应该是-1 
	for(int i=0;i<=n+m;++i){
		printf("%d ",(int)(a[i].r/cnt+0.5));
		//注意输出应该强转格式。 
	}
	puts("");
}

int main(){
	init();
	return 0; 
}

lp2597 ZJOI2012 灾难

这一题的思路还是比较清晰的。虽然可以当作支配树模板,但是事实上可以上一个LCA来代替支配树。
具体来说,就是把有多个食物的消费者连到它所有食物的LCA。然后计算一遍子树大小就好了。
注意连点之前要先拓扑排序。

#include<iostream>
#include<cstdio>
#include<queue>

#define Swap(A,B) (A^=B^=A^=B)

struct ee{
	int v;
	int nxt;
}e[2500005];
//h,树上的节点。g,拓扑排序的节点。to,每个节点的所有父亲。 
int h[70000],g[70000],et=0,dep[70000],fa[70000][18],in[70000],to[70000],sz[70000],loc[70000];
int n,tp=0;
inline void add(int *H,int U,int V){
	if(U==V){
		return;
	}
	e[++et]=(ee){V,H[U]};
	H[U]=et;
}
inline void srt(){
	std::queue<int> q;
	for(int i=1;i<=n;++i){
		if(!in[i]){
			q.push(i);
		}
	}
	int p=0;
	while(!q.empty()){
		p=q.front();
		q.pop();
		loc[++tp]=p;
		for(int i=g[p];i;i=e[i].nxt){
			--in[e[i].v];
			if(!in[e[i].v]){
				q.push(e[i].v);
			}
		}
	}
}
inline int lca(int X,int Y){
	if(dep[X]<dep[Y]){
		Swap(X,Y);
	}
	for(int i=17;i>=0;--i){
		if(dep[X]-(1<<i)>=dep[Y]){
			X=fa[X][i]; 
		}
	}
	if(X==Y){
		return X;
	}
	for(int i=17;i>=0;--i){
		if(fa[X][i]!=fa[Y][i]){
			X=fa[X][i],Y=fa[Y][i];
		}
	}
	return fa[X][0];
}
inline void prpr(){
	int nw=0;
	for(int i=1,j;i<=n;++i){
		j=loc[i];
		nw=e[to[j]].v;
		for(int k=to[j];k;k=e[k].nxt){
			nw=lca(nw,e[k].v);
		}
//		请注意这里的nw可能不存在任何父亲节点。 
		add(h,nw,j);
		fa[j][0]=nw;
		dep[j]=dep[nw]+1;
		for(int k=1;k<=17;++k){
			fa[j][k]=fa[fa[j][k-1]][k-1];
		}
	}
}
inline void dfs(int X){
	for(int i=h[X];i;i=e[i].nxt){
		dfs(e[i].v);
		sz[X]+=sz[e[i].v];
	}
	++sz[X];
}
void init(){
	scanf("%d",&n);
	int x; 
	for(int i=1;i<=n;++i){
		scanf("%d",&x);
		while(x){
			++in[i];
			add(g,x,i);
			add(to,i,x);
			scanf("%d",&x);
		}
	}
	srt();
	prpr();
	dfs(0);
	for(int i=1;i<=n;++i){
		printf("%d\n",sz[i]-1);
	}
}

int main(){
	init();
	return 0;
}

lp3615 如厕计划

观察题面我们可以发现,如果男性比女性多的话,显然是无法符合要求的;否则,女性使用男女通用单间的个数不应该超过男女性的人数差。
我们不妨把部分的男女性人数差个女性视为男性——这是一种贪心的做法。于是就把女性多于男性的情况归约成了男性和女性相同的情况。
在这种情况下,任何时候都不应该有女性进入男女通用单间。这也就意味着,如果令男性为1,女性为-1,那么在任何时候,数列的前缀和都应该保证为大于等于-1。
然而,要保证一个前缀和始终大于某个数是很难维护的——这是因为位于后面的信息我们无法得到。我们考虑维护后缀和,这样就可以预先得到位于后面的信息了。
于是,我们考虑对每一段计算它的贡献。

#include<iostream>
#include<cstdio>
#include<cstring>

inline long long Max(long long A,long long B){
	return A>B?A:B;
}
inline long long Min(long long A,long long B){
	return A<B?A:B;
}

long long n,m;
int len;
char ch[200005];
long long dlt=0,a[200005],mx[200005],l[200005];
void init(){
	scanf("%lld%lld",&n,&m);
	for(int i=1;i<=m;++i){
		std::cin>>ch+1;
		len=strlen(ch+1);
		scanf("%lld",&l[i]);
		for(int j=len;j>=1;--j){
			a[i]+=(ch[j]=='M'?1:-1); 
			mx[i]=Max(mx[i],a[i]);
		}
		dlt+=a[i]*l[i];
	}
	if(dlt>0){
		puts("-1");
		return ;
	}
	long long nw=0,ans=1;
	for(int i=m;i>=1;--i){
		if(a[i]>0){
			ans=Max(ans,(l[i]-1)*a[i]+nw+mx[i]);
		}else{
			ans=Max(ans,nw+mx[i]);
		}
		nw+=l[i]*a[i];
	}
	printf("%lld",ans-1);
}

int main(){
	freopen("queue.in","r",stdin);
	freopen("queue.out","w",stdout);
	init();
	return 0;
}