分治最短路

补很久以前的一个坑。

BZOJ4449 Distance on Triangulation

给定一个凸n边形，以及它的三角剖分。再给定q个询问，每个询问是一对凸多边行上的顶点(a,b)，问点a最少经过多少条边(可以是多边形上的边，也可以是剖分上的边)可以到达点b。

n ≤ 50000，q ≤ 100000

题解

显然正多边形的三角剖分是一个平面图，每一条剖分的边可以将正多边形分成有一条重边的两个独立的新多边形，显然这一个过程是可以用分治来实现的。

我们对于分治过程中的多边形进行重新编号，找到两端点数最平均的边割去，对于点集V，边集E和询问集Q分别开三个vector传入函数中。

如果点集的大小等于了3，我们就可以对单个三角形进行直接计算，否则对于集合进行左右的分治，贡献大小使用两次BFS计算，时间复杂度O(n log ⁡n)。

#include<bits/stdc++.h>
#define co const
#define il inline
template<class T> T read(){
	T x=0,w=1;char c=getchar();
	for(;!isdigit(c);c=getchar())if(c=='-') w=-w;
	for(;isdigit(c);c=getchar()) x=x*10+c-'0';
	return x*w;
}
template<class T>il T read(T&x){
	return x=read<T>();
}
using namespace std;
typedef long long LL;

co int N=100000+1;
int n,m;
int he[N],nx[N<<1],to[N<<1],cnt;
int siz[N],le[N],ri[N];
int vis[N],dx[N],dy[N];
int ans[N];

void add_edge(int x,int y){
	to[++cnt]=y,nx[cnt]=he[x],he[x]=cnt;
	to[++cnt]=x,nx[cnt]=he[y],he[y]=cnt;
}

struct edge {int u,v;};
struct quiz {int a,b,id;};

void bfs(int s,int dis[]){
	dis[s]=0;
	queue<int> q;q.push(s);
	while(q.size()){
		int x=q.front();q.pop();
		for(int i=he[x],y;y=to[i],i;i=nx[i])
			if(vis[y]&&dis[y]==1e9) dis[y]=dis[x]+1,q.push(y);
	}
}
void solve(vector<int> ve,vector<edge> ed,vector<quiz> qr){
	if(!qr.size()) return;
	if(ve.size()==3){
		for(unsigned i=0;i<qr.size();++i) ans[qr[i].id]=qr[i].a!=qr[i].b;
		return;
	}
	vector<int> v1,v2;
	vector<edge> e1,e2;
	vector<quiz> q1,q2;
	int n=ve.size(),m=ed.size();
	siz[ve[0]]=1;
	for(int i=1;i<n;++i) siz[ve[i]]=siz[ve[i-1]]+1;
	int x=0,y=0,misz=1e9;
	for(int i=0;i<m;++i){
		int u=ed[i].u,v=ed[i].v,len=siz[v]-siz[u]-1;
		if(max(len,n-len-2)<misz) misz=max(len,n-len-2),x=u,y=v;
	}
	for(int i=0;i<n;++i){
		if(x<=ve[i]&&ve[i]<=y) le[ve[i]]=1,v1.push_back(ve[i]);
		if(ve[i]<=x||ve[i]>=y) ri[ve[i]]=1,v2.push_back(ve[i]);
	}
	for(int i=0;i<m;++i){
		int u=ed[i].u,v=ed[i].v;
		if(le[u]&&le[v]) e1.push_back(ed[i]);
		if(ri[u]&&ri[v]) e2.push_back(ed[i]);
	}
	for(unsigned i=0;i<qr.size();++i){
		int a=qr[i].a,b=qr[i].b;
		if(le[a]&&le[b]) q1.push_back(qr[i]);
		if(ri[a]&&ri[b]) q2.push_back(qr[i]);
	}
	for(int i=0;i<n;++i) vis[ve[i]]=1,dx[ve[i]]=dy[ve[i]]=1e9;
	bfs(x,dx),bfs(y,dy);
	for(unsigned i=0;i<qr.size();++i){
		int a=qr[i].a,b=qr[i].b;
		ans[qr[i].id]=min(ans[qr[i].id],min(min(dx[a]+dx[b],dy[a]+dy[b]),min(dx[a]+dy[b],dy[a]+dx[b])+1));
	}
	for(int i=0;i<n;++i) vis[ve[i]]=le[ve[i]]=ri[ve[i]]=0;
	solve(v1,e1,q1),solve(v2,e2,q2);
}
int main(){
	read(n);
	for(int i=1;i<n;++i) add_edge(i,i+1);
	add_edge(n,1);
	vector<edge> ed;
	for(int i=1;i<=n-3;++i){
		int x=read<int>(),y=read<int>();
		add_edge(x,y);
		if(x>y) swap(x,y);
		ed.push_back((edge){x,y});
	}
	vector<int> ve;
	for(int i=1;i<=n;++i) ve.push_back(i);
	read(m);
	vector<quiz> qr;
	for(int i=1;i<=m;++i){
		int x=read<int>(),y=read<int>();
		if(x>y) swap(x,y);
		qr.push_back((quiz){x,y,i});
	}
	fill(ans+1,ans+m+1,1e9);
	solve(ve,ed,qr);
	for(int i=1;i<=m;++i) printf("%d
",ans[i]);
	return 0;
}

「ZJOI2016」旅行者

小Y来到了一个新的城市旅行。她发现了这个城市的布局是网格状的，也就是有n条从东到西的道路和m条从南到北的道路，这些道路两两相交形成n×m个路口 (i,j) (1 ≤ i ≤ n,1 ≤ j ≤ m)。

她发现不同的道路路况不同，所以通过不同的路口需要不同的时间。通过调查发现，从路口(i,j)到路口(i,j+1)需要时间 r(i,j)，从路口(i,j)到路口(i+1,j)需要时间c(i,j)。注意这里的道路是双向的。

小Y有q个询问，她想知道从路口(x₁,y₁)到路口(x₂,y₂)最少需要花多少时间。

对于所有的测试数据，n×m ≤ 2×10⁴，q ≤ 10⁵。

题解

我们把矩形顺着比较短的一边从中间切成两半，从这条中线上每个点做一次到这个矩形内每个点的Dijkstra，把两个端点都在这个矩形里的询问用中线上每个点到两个端点的最短路更新，然后分治两半边。

复杂度口胡：我们知道这个算法每次更新答案要做比较小的一边的边长次Dijkstra，所以当是一个正方形的时候分治的复杂度最高，而询问一定是传到底的时候复杂度最高。

因为每次切一半两边都是一样的，所以每个询问只要传到底，不管往哪边传都可以，为了方便分析不妨设每个询问往两边各传一半。

T(n,m,Q)表示短边长为n，长边长为m，矩形里有Q个询问的复杂度，T(n,Q)表示边长为n，里边有Q个询问的正方形的复杂度。

[egin{array}{l}{T(n, m, Q)=n^{2} m log (n m)+n Q+2 * Tleft(frac{m}{2}, n, frac{Q}{2} ight)} \ {T(n, m, Q)=n^{2} m log (n m)+n Q+2 *left(left(frac{m}{2} ight)^{2} n log left(frac{n m}{2} ight)+frac{m Q}{4} ight)+4 * Tleft(frac{n}{2}, frac{n}{2}, frac{Q}{4} ight)} \ {T(n, m, Q)=T(n, Q)=n^{3} log left(n^{2} ight)+n Q+frac{n^{3}}{2} * log left(frac{n^{2}}{2} ight)+frac{n Q}{2}+4 * Tleft(frac{n}{2}, frac{Q}{4} ight)} \ {T(n, Q)=3 n^{3} log (n)-frac{n^{3}}{2}+frac{3}{2} n Q+4 Tleft(frac{n}{2}, frac{Q}{4} ight)}end{array} ]

总复杂度

[=sum_{i=1}^{log left(n^{2} ight)}left(frac{3}{2^{i-1}} n^{3} log left(frac{n}{2^{i-1}} ight)-frac{n^{3}}{2^{i}}+frac{3}{2^{i}} n Q ight) ]

把最后的式子化简一下，面积为S的话，复杂度是O(S×√S×log S)的。

#include<bits/stdc++.h>
#define co const
#define il inline
template<class T> T read(){
	T x=0,w=1;char c=getchar();
	for(;!isdigit(c);c=getchar())if(c=='-') w=-w;
	for(;isdigit(c);c=getchar()) x=x*10+c-'0';
	return x*w;
}
template<class T>il T read(T&x){
	return x=read<T>();
}
using namespace std;
typedef long long LL;

co int N=200000+1,INF=0x3f3f3f3f;
int n,m,q,ans[N],dis[N];

struct edge{
	edge*nx;
	int to,w;
}pool[N<<1],*p=pool,*lnk[N];

struct ask{
	int x,y,xx,yy,id;
}seq[N],bin[N];

il int zip(int i,int j){
	return (i-1)*m+j;
}
il void unzip(int v,int&i,int&j){
	i=(v-1)/m+1,j=(v-1)%m+1;
}
il void add_edge(int u,int v,int w){
	*++p=(edge){lnk[u],v,w},lnk[u]=p;
}
il void chkmin(int&x,int y){
	if(y<x) x=y;
}

void spfa(int s,int xl,int xr,int yl,int yr,int fl){
	int t=dis[s];
	for(int i=xl;i<=xr;++i)
		for(int j=yl;j<=yr;++j) dis[zip(i,j)]=fl?dis[zip(i,j)]+t:INF;
	queue<int> q;
	static bool inq[N];
	dis[s]=0,q.push(s),inq[s]=1;
	while(q.size()){
		int u=q.front();q.pop();
		inq[u]=0;
		for(edge*i=lnk[u];i;i=i->nx){
			int xx,yy;
			unzip(i->to,xx,yy);
			if(xx<xl||xx>xr||yy<yl||yy>yr) continue;
			if(dis[u]+i->w<dis[i->to]){
				dis[i->to]=dis[u]+i->w;
				if(!inq[i->to]) q.push(i->to),inq[i->to]=1;
			}
		}
	}
}
void solve(int xl,int xr,int yl,int yr,int ql,int qr){
	if(ql>qr) return;
	if(xr-xl>yr-yl){
		int mid=(xl+xr)>>1;
		for(int i=yl;i<=yr;++i){
			spfa(zip(mid,i),xl,xr,yl,yr,i-yl);
			for(int j=ql;j<=qr;++j)
				chkmin(ans[seq[j].id],dis[zip(seq[j].x,seq[j].y)]+dis[zip(seq[j].xx,seq[j].yy)]);
		}
		int l=ql-1,r=qr+1;
		for(int i=ql;i<=qr;++i){
			if(seq[i].xx<mid&&seq[i].x<mid) bin[++l]=seq[i];
			if(seq[i].x>mid&&seq[i].xx>mid) bin[--r]=seq[i];
		}
		copy(bin+ql,bin+qr+1,seq+ql);
		solve(xl,mid-1,yl,yr,ql,l);
		solve(mid+1,xr,yl,yr,r,qr);
	}
	else{
		int mid=(yl+yr)>>1;
		for(int i=xl;i<=xr;++i){
			spfa(zip(i,mid),xl,xr,yl,yr,i-xl);
			for(int j=ql;j<=qr;++j)
				chkmin(ans[seq[j].id],dis[zip(seq[j].x,seq[j].y)]+dis[zip(seq[j].xx,seq[j].yy)]);
		}
		int l=ql-1,r=qr+1;
		for(int i=ql;i<=qr;++i){
			if(seq[i].yy<mid&&seq[i].y<mid) bin[++l]=seq[i];
			if(seq[i].y>mid&&seq[i].yy>mid) bin[--r]=seq[i];
		}
		copy(bin+ql,bin+qr+1,seq+ql);
		solve(xl,xr,yl,mid-1,ql,l);
		solve(xl,xr,mid+1,yr,r,qr);
	}
}
int main(){
	memset(ans,0x3f,sizeof ans);
	read(n),read(m);
	for(int i=1;i<=n;++i)
		for(int j=1;j<m;++j){
			int val=read<int>();
			add_edge(zip(i,j),zip(i,j+1),val);
			add_edge(zip(i,j+1),zip(i,j),val);
		}
	for(int i=1;i<n;++i)
		for(int j=1;j<=m;++j){
			int val=read<int>();
			add_edge(zip(i,j),zip(i+1,j),val);
			add_edge(zip(i+1,j),zip(i,j),val);
		}
	read(q);
	for(int i=1;i<=q;++i)
		read(seq[i].x),read(seq[i].y),read(seq[i].xx),read(seq[i].yy),seq[i].id=i;
	solve(1,n,1,m,1,q);
	for(int i=1;i<=q;++i) printf("%d
",ans[i]);
	return 0;
}

CodeChef-QGRID Querying on a Grid

大厨在一家教学机构上班。一个学生跑来问大厨一道题，题目如下：

考虑 M 行 N 列，共 M × N 个格点（而非格子）的网格，第 i 行第 j 列的格点为 (i, j)。每个格点均有点权，初始时均为 0。

图中有横竖两类边，所有边均为无向边：

• 对于 i < M，有一条权值为 down(i, j) 的边连接 (i, j) 和 (i + 1, j)；
• 对于 j < N，有一条权值为 right(i, j) 的边连接 (i, j) 和 (i, j + 1)。

我们定义路径的长度为其经过的边的权值和。两个格点 (i₁, j₁) 和 (i₂, j₂) 间的最短路径即为长度最短的路径。当然，格点的权值与路径无关。

题目中有两类操作：

• 1 i1 j1 i2 j2 c：对于 (i₁, j₁) 和 (i_2, j₂) 间的最短路上的所有格点，令其权值加 c；
• 2 i j：求格点 (i, j) 的权值。

1 ≤ M ≤ 3, 1 ≤ N ≤ 100000, 1 ≤ Q ≤ 100000

题解

我们考虑两个点之间的最短路，若其中经过了某个点(x,y)，那么这两点之间的最短路一定在以(x,y)为根的最短路树上。

网格图有一个经典套路就是分治，那么我们对于当前分治区间[l,r]，分别建出纵坐标在这个区间内以mid上的m个点的最短路树。

对于一个操作，我们可以找到所有包含j₁,j₂的最短路树，然后找到最小的一条路径，在树上打标记并用数据结构维护。
具体来说由于最短路树上的标记一定是从根节点到后代节点的路径，那么我们只需要求出树的DFS序，这样一个点的权值就是它后代标记和，用树状数组维护即可。注意若询问点是最短路树根，则答案需要除以2。

复杂度O( (n+Q) log² n )​。实测空间：1289728kB。

#include<bits/stdc++.h>
#define co const
#define il inline
template<class T>T read(){
	T x=0,w=1;char c=getchar();
	for(;!isdigit(c);c=getchar())if(c=='-') w=-w;
	for(;isdigit(c);c=getchar()) x=x*10+c-'0';
	return x*w;
}
template<class T>T read(T&x){
	return x=read<T>();
}
using namespace std;
typedef long long LL;

co int N=300000+1;
int n,m,q;
int id[N][3],tot;
vector<pair<int,LL> > e[N];

struct SPT{
	LL dis[N];
	int fa[N],pos[N],lst[N],dfn;
	vector<int> son[N];
	
	void dfs(int x){
		pos[x]=++dfn;
		for(unsigned i=0;i<son[x].size();++i) dfs(son[x][i]);
		lst[x]=dfn;
	}
	void spfa(int l,int r,int s){
		priority_queue<pair<LL,int>,vector<pair<LL,int> >,greater<pair<LL,int> > > pq;
		static bool vis[N];
		for(int i=l;i<=r;++i) dis[i]=1e18,vis[i]=0;
		dis[s]=0,pq.push(make_pair(dis[s],s));
		while(pq.size()){
			int x=pq.top().second;
			pq.pop();
			if(vis[x]) continue;
			vis[x]=1;
			for(unsigned i=0;i<e[x].size();++i){
				int y=e[x][i].first;
				LL w=e[x][i].second;
				if(y<l||y>r||dis[y]<=dis[x]+w) continue;
				dis[y]=dis[x]+w,fa[y]=x;
				pq.push(make_pair(dis[y],y));
			}
		}
		for(int i=l;i<=r;++i)
			if(i!=s) son[fa[i]].push_back(i);
		dfs(s);
	}
	
	LL c[N];
	il int lowbit(int i) {return i&-i;}
	void modify(int x,LL v){
		for(int p=pos[x];p<=dfn;p+=lowbit(p)) c[p]+=v;
	}
	LL query(int x){
		LL ans=0;
		for(int p=lst[x];p;p-=lowbit(p)) ans+=c[p];
		for(int p=pos[x]-1;p;p-=lowbit(p)) ans-=c[p];
		return ans;
	}
}T[18][3];

void build(int l,int r,int dep){
	int mid=(l+r)>>1;
	for(int i=0;i<n;++i)
		T[dep][i].spfa(id[l][0],id[r][n-1],id[mid][i]);
	if(l<mid) build(l,mid-1,dep+1);
	if(r>mid) build(mid+1,r,dep+1);
}
LL mids;
void calc(int l,int r,int dep,int x,int y){
	int mid=(l+r)>>1;
	for(int i=0;i<n;++i) mids=min(mids,T[dep][i].dis[x]+T[dep][i].dis[y]);
	if(x<=id[mid][n-1]&&y>=id[mid][0]) return;
	if(y<id[mid][0]) calc(l,mid-1,dep+1,x,y);
	else calc(mid+1,r,dep+1,x,y);
}
void modify(int l,int r,int dep,int x,int y,LL w){
	int mid=(l+r)>>1;
	for(int i=0;i<n;++i)
		if(T[dep][i].dis[x]+T[dep][i].dis[y]==mids){
			T[dep][i].modify(x,w),T[dep][i].modify(y,w);
			return;
		}
	if(y<id[mid][0]) modify(l,mid-1,dep+1,x,y,w);
	else modify(mid+1,r,dep+1,x,y,w);
}
LL query(int l,int r,int dep,int x){
	int mid=(l+r)>>1;
	LL ans=0;
	for(int i=0;i<n;++i){
		LL sum=T[dep][i].query(x);
		if(id[mid][i]==x) sum>>=1;
		ans+=sum;
	}
	if(x<id[mid][0]) ans+=query(l,mid-1,dep+1,x);
	if(x>id[mid][n-1]) ans+=query(mid+1,r,dep+1,x);
	return ans;
}

int main(){
	read(n),read(m),read(q);
	for(int i=0;i<m;++i)
		for(int j=0;j<n;++j) id[i][j]=++tot;
	for(int i=0;i<n-1;++i)
		for(int j=0;j<m;++j){
			int x=id[j][i],y=id[j][i+1];
			LL w=read<LL>();
			e[x].push_back(make_pair(y,w)),e[y].push_back(make_pair(x,w));
		}
	for(int i=0;i<n;++i)
		for(int j=0;j<m-1;++j){
			int x=id[j][i],y=id[j+1][i];
			LL w=read<LL>();
			e[x].push_back(make_pair(y,w)),e[y].push_back(make_pair(x,w));
		}
	build(0,m-1,0);
	while(q--){
		if(read<int>()==1){
			int x1=read<int>()-1,y1=read<int>()-1,x2=read<int>()-1,y2=read<int>()-1;
			LL w=read<LL>();
			int x=id[y1][x1],y=id[y2][x2];
			if(x>y) swap(x,y);
			mids=1e18;
			calc(0,m-1,0,x,y);
			modify(0,m-1,0,x,y,w);
		}
		else{
			int x=read<int>()-1,y=read<int>()-1;
			printf("%lld
",query(0,m-1,0,id[y][x]));
		}
	}
	return 0;
}