隐式Dijkstra:在状态集合中用优先队列求前k小

这种技巧是挺久以前接触的了,最近又突然遇到几道新题,于是总结了一下体会。

这种算法适用的前提是,标题所述的“状态集合”大到不可枚举(否则枚举就行了qaq),且\(k\)一般是在\(10^6\)这个数量级以下。

前置技能:Dijkstra算法,及其思想和正确性证明
传送门1:思想和正确性证明
传送门2:优先队列优化dijkstra

先看一个问题:##

\(m\)(\(2 \leq m \leq 10\))个长度为\(n\)(\(n \leq 10^5\))的整数序列,从每个序列选一个数相加,求所得的和中第\(k\)(\(k \leq 10^5\))大的。

(首先显而易见要把每个序列排序,从大到小)

考虑\(m=2\)的情况:###

二分答案\(A\),就可以对序列1中的每一个元素,找到能够使总和在\(A\)以上的、序列2中可以与它配对的元素集合。易知这个集合是序列2的一个前缀,且前缀的长度随\(i\)递增,故对一个\(A\)求出所有这样的前缀只需要线性时间。而总和大于\(A\)的方案数,就是这些前缀的长度之和。

上代码:
#include <bits/stdc++.h>
using namespace std;
#define iinf 1000000000
#define linf 1000000000000000000LL
#define ulinf 10000000000000000000ull
#define MOD1 1000000007LL
#define mpr make_pair
typedef long long LL;
typedef unsigned long long ULL;
typedef unsigned long UL;
typedef unsigned short US;
typedef pair < int , int > pii;
clock_t __stt;
inline void TStart(){__stt=clock();}
inline void TReport(){printf("\nTaken Time : %.3lf sec\n",(double)(clock()-__stt)/CLOCKS_PER_SEC);}
template < typename T > T MIN(T a,T b){return a<b?a:b;}
template < typename T > T MAX(T a,T b){return a>b?a:b;}
template < typename T > T ABS(T a){return a>0?a:(-a);}
template < typename T > void UMIN(T &a,T b){if(b<a) a=b;}
template < typename T > void UMAX(T &a,T b){if(b>a) a=b;}
int n,m,k,a[2][100005];
bool check(int v){
	int i,j,t,cnt;
	for(t=0;t<n && a[0][t]+a[1][0]<v;++t);
	cnt=n-t;
	for(i=1;i<n;++i){
		for(;t<n && a[0][t]+a[1][i]<v;++t);
		cnt+=n-t;
		if(cnt>=k) return 1;
	}
	return cnt>=k;
}
int main(){
	// inputting start
	// 数据结构记得初始化! n,m别写反!
	scanf("%d%d%d",&n,&m,&k); //m=2
	int i,j;
	for(i=0;i<2;++i){
		for(j=0;j<n;++j){
			scanf("%d",&a[i][j]);
		}
	}
	#ifdef LOCAL
		TStart();
	#endif
	// calculation start
	// 数据结构记得初始化! n,m别写反!
	sort(a[0],a[0]+n);
	sort(a[1],a[1]+n);
	reverse(a[1],a[1]+n);
	int low=0,high=iinf,mid;
	while(low<high){
		mid=((low+high+1)>>1);
		if(check(mid))
			low=mid;
		else
			high=mid-1;
	}
	printf("%d\n",low);
	#ifdef LOCAL
		TReport();
	#endif
	return 0;
}

这种方法不是本文的重点,但是很大一部分可以用马上会介绍的隐式Dijkstra解决的问题,都可以用这种二分的方法,以一个略差的性能解决。

有的问题看似可以用隐式Dijkstra,其实却会出大问题。碰到这种情况,可以考虑改用这类二分答案的做法。

推荐下面这道题。

USACO 2016DEC PLATINUM T3
题解

对于一般情况(\(m \leq 10\)),就要换方法了。###

(敲黑板!重点重点!)

我们在脑海中建一张图,每个结点对应着一种选取方案。方案\(A\)\(B\)有一条有向边,当且仅当\(A\)对应的总和大于\(B\)的总和,边权是两者总和的差的绝对值。

记所有序列的最大数组成\(S\)状态,那么对任意\(T\)状态:\(cost(T)=cost(S)-dist(S,T)\)

那么我们要求的,就是到\(S\)状态距离第\(k\)短的点(包括\(S\)自己)。注意到虽然这是个DAG,但因为结点过多无法DP。

所以,考虑Dijkstra。

先上结论:

Dijkstra在正权图上运行时,优先队列每次弹出的结点到\(S\)的距离,一定是递增的。

证明:

若有\(dist(S,P)>dist(S,Q)\)\(P\)先于\(Q\)弹出,则:

\(P\)弹出前的瞬间,因为\(P\)是优先队列中距\(S\)最近的,所以\(Q\)\(S\)比优先队列中任意状态到\(S\)更近,且\(Q\)\(P\)弹出前不在优先队列里,也从未被压入过。

所以\(Q\)不可能由优先队列中的状态经过若干松弛操作而得到。

\(Q\)不可能在\(P\)弹出后被压入队列,也就不可能在\(P\)之后弹出。

产生矛盾,证毕。

推论:

对任意的\(k\),优先队列里最先弹出的\(k\)个结点,一定是到源点最近的\(k\)个点。

我会做啦!啊哈哈哈哈哈!
跑一遍Dijkstra,优先队列弹出的第\(k\)个点就是答案!

且慢,算一下复杂度。

记各结点的平均度数为\(d\),因为做了\(k\)轮松弛,每轮压入了\(d\)个新结点,故总时间复杂度为:

\(O(kd \cdot log\ k)\),约等于\(10^{56}\)。emmmmm……

(敲黑板!重点又来了!)

这种算法,优化的思路之一是优化连边

易知,对于固定的\(T\),任意\(S\)-\(T\)路径的权值和全都相等。(这只是针对本题的性质,其他此类问题不一定满足。)

所以,如果删掉一些边,使得图的连通性不变,那么答案也不变。

换句话说,要删掉一些边,使得\(S\)到每个点仍有至少一条路径。

所以修改连边策略:

\(A\)\(B\)有一条有向边,仅当\(A\)对应的总和大于\(B\)的总和,边权是两者总和的差的绝对值。

此外,必须满足,\(A\)\(B\)对应的方案,只在一个序列里选的数下标不同(其它\(m-1\)个序列里选的都完全相同),而且,所选的两个下标不同的数,也一定是相邻元素。(数组已排序)

要删掉一些边,使得\(S\)到每个点仍有至少一条路径。

满足了没?满足了。

而现在,\(d \leq m\)。所以总复杂度降为\(O(mk \cdot log\ k)\),bingo!

具体实现的hint:

实际上,完全可以枚举得到每一个点的所有邻居,于是就没有必要建立邻接表了。

我不知道这种算法的正式名称(可能并不存在?),所以就在本文中叫它“隐式Dijkstra”了。

上代码:
// 隐式Dijkstra的做法,代码针对m=2的情况,m<=10的就请自行实现啦,差别不大的
#include <bits/stdc++.h>
using namespace std;
#define iinf 1000000000
#define linf 1000000000000000000LL
#define ulinf 10000000000000000000ull
#define MOD1 1000000007LL
#define mpr make_pair
typedef long long LL;
typedef unsigned long long ULL;
typedef unsigned long UL;
typedef unsigned short US;
typedef pair < int , int > pii;
clock_t __stt;
inline void TStart(){__stt=clock();}
inline void TReport(){printf("\nTaken Time : %.3lf sec\n",(double)(clock()-__stt)/CLOCKS_PER_SEC);}
template < typename T > T MIN(T a,T b){return a<b?a:b;}
template < typename T > T MAX(T a,T b){return a>b?a:b;}
template < typename T > T ABS(T a){return a>0?a:(-a);}
template < typename T > void UMIN(T &a,T b){if(b<a) a=b;}
template < typename T > void UMAX(T &a,T b){if(b>a) a=b;}
int n,m,k,a[2][100005];
struct P{
	int x,y;
	int val(){return a[0][x]+a[1][y];}
	bool operator <(P b) const{
		return a[0][x]+a[1][y]<b.val();
	}
};
P make_P(int x,int y){
	P R;
	R.x=x;R.y=y;
	return R;
}
priority_queue < P > pq;
map < pii , int > vis;
int main(){
	// inputting start
	// 数据结构记得初始化! n,m别写反!
	scanf("%d%d%d",&n,&m,&k); //m=2
	int i,j;
	for(i=0;i<2;++i){
		for(j=0;j<n;++j){
			scanf("%d",&a[i][j]);
		}
	}
	#ifdef LOCAL
		TStart();
	#endif
	// calculation start
	// 数据结构记得初始化! n,m别写反!
	sort(a[0],a[0]+n);
	reverse(a[0],a[0]+n);
	sort(a[1],a[1]+n);
	reverse(a[1],a[1]+n);
	pq.push(make_P(0,0));
	while(!pq.empty()){
		P cur=pq.top();
		pq.pop();
		--k;
		if(!k){
			printf("%d\n",cur.val());
			return 0;
		}
		if(cur.x<n-1 && !vis[mpr(cur.x+1,cur.y)]){
			vis[mpr(cur.x+1,cur.y)]=1;
			pq.push(make_P(cur.x+1,cur.y));
		}
		if(cur.y<n-1 && !vis[mpr(cur.x,cur.y+1)]){
			vis[mpr(cur.x,cur.y+1)]=1;
			pq.push(make_P(cur.x,cur.y+1));
		}
	}
	#ifdef LOCAL
		TReport();
	#endif
	return 0;
}

小结1:隐式Dijkstra的适用条件##

1、合理的数据范围

要注意,这种算法适用的前提是,标题所述的“状态集合”大到不可枚举(否则直接枚举就行了qaq),而且k一般是在\(10^6\)这个数量级以下。

这里的“状态集合”指的是构出的图中,所有结点的集合。

2、易于表示、比较的状态

优先队列里的操作是基于比较的。

如果比较大小的复杂度过高会TLE,如果存储状态的空间复杂度过大会MLE。

3、正权图

边权必须都是非负数。

实践出真知:SGU421##

题意:

给长度为\(n\)(\(n \leq 10^4\))的整数(正负均可)数列,选\(m\)(\(m \leq 13\))个数相乘,问第\(k\)(\(k \leq 10^4\))大的乘积。

题解:

先把数列按正负分成两个数组,再分别排序。

以下是隐式Dijkstra的几个要素:

状态含义:

一个状态代表一种选取方案,即一个正数子集和一个负数子集的二元组。

状态间的大小关系:

即乘积的大小关系。先比较符号(负数数量的奇偶性),再比绝对值(高精度)。

起始状态\(S\)

考虑使用多个起始状态。每个\(S\)是负数集合大小为\(x\)(\(1 \leq x \leq n\))的状态中最大的。

连边方案:

如果状态\(u\)\(v\)只有一个选的数不同,而这两个不同的数是相邻元素,则从较大状态向较小状态连边,边权是权值相除的商。
注意到这里Dijkstra算法对距离的定义不再是路径权值和,而是路径权值积。

答案就是优先队列弹出的第\(k\)个状态的值。

正确性证明:

这张图不是弱连通的,但是每个弱连通分量(包含的状态拥有相同的负数子集大小)都恰有一个起始状态,故正确性依然保持。

上代码:
#include <bits/stdc++.h>
using namespace std;
#define iinf 2000000000
#define linf 1000000000000000000LL
#define ulinf 10000000000000000000ull
#define MOD1 1000000007LL
#define mpr make_pair
typedef long long LL;
typedef unsigned long long ULL;
typedef unsigned long UL;
typedef unsigned int US;
typedef pair < int , int > pii;
clock_t __stt;
inline void TStart(){__stt=clock();}
inline void TReport(){printf("\nTaken Time : %.3lf sec\n",(double)(clock()-__stt)/CLOCKS_PER_SEC);}
template < typename T > T MIN(T a,T b){return a<b?a:b;}
template < typename T > T MAX(T a,T b){return a>b?a:b;}
template < typename T > T ABS(T a){return a>0?a:(-a);}
template < typename T > void UMIN(T &a,T b){if(b<a) a=b;}
template < typename T > void UMAX(T &a,T b){if(b>a) a=b;}
int n,m,k;
vector < int > neg,pos;
struct bigint{
	int len,cnt0;
	int d[85];
	void init(){
		memset(d,0,sizeof(d));
		len=1;
		d[0]=1;
		cnt0=0;
	}
	void reduct(){
		while(len>1 && !d[len-1]) --len;
	}
	void multiply(int v){
		if(!v){
			++cnt0;
			return;
		}
		int i;
		for(i=0;i<len;++i) d[i]*=v;
		for(i=0;i<len;++i){
			d[i+1]+=d[i]/10;
			d[i]%=10;
		}
		while(d[len]){
			d[len+1]+=d[len]/10;
			d[len++]%=10;
		}
	}
	void divide(int v){
		if(!v){
			--cnt0;
			return;
		}
		int i,j,c=0;
		for(i=len-1;i>=0;--i){
			c=c*10+d[i];
			d[i]=0;
			if(c>=v){
				d[i]=c/v;
				c%=v;
			}
		}
		for(i=0;i<len;++i){
			d[i+1]+=d[i]/10;
			d[i]%=10;
		}
		while(d[len]){
			d[len+1]+=d[len]/10;
			d[len++]%=10;
		}
		reduct();
	}
	void print(bool sig){
		if(cnt0){
			printf("0\n");
			return;
		}
		int i;
		reduct();
		if(sig && !(len==1&&d[0]==0)) printf("-");
		for(i=len-1;i>=0;--i) printf("%d",d[i]);
		printf("\n");
	}
};
bool operator <(bigint &A,bigint &B){
	if(A.cnt0) return !B.cnt0;
	if(B.cnt0) return 0;
	if(A.len!=B.len) return A.len<B.len;
	int i;
	for(i=A.len-1;i>=0;--i){
		if(A.d[i]!=B.d[i]) return A.d[i]<B.d[i];
	}
	return 0;
}
struct state{
	int vp[15],vn[15],cp,cn;
	bigint val;
	state(){
		cp=cn=0;
		val.init();
	}
	bool sign(){
		return cn&1;
	}
	bool editp(int p,int d){
		if(vp[p]+d<0 || vp[p]+d>=(int)pos.size()) return 0;
		if((p && vp[p-1]==vp[p]+d)||(p<cp-1 && vp[p+1]==vp[p]+d)) return 0;
		val.divide(pos[vp[p]]);
		vp[p]+=d;
		val.multiply(pos[vp[p]]);
		return 1;
	}
	bool editn(int p,int d){
		if(vn[p]+d<0 || vn[p]+d>=(int)neg.size()) return 0;
		if((p && vn[p-1]==vn[p]+d)||(p<cn-1 && vn[p+1]==vn[p]+d)) return 0;
		val.divide(neg[vn[p]]);
		vn[p]+=d;
		val.multiply(neg[vn[p]]);
		return 1;
	}
	int super_cdd(){
		int ret=cp*101+cn,i;
		for(i=0;i<cp;++i){
			ret=ret*101+vp[i];
		}
		for(i=0;i<cn;++i){
			ret=ret*101+vn[i];
		}
		return ret;
	}
};
const bool operator <(state A,state B){
	if(A.sign()!=B.sign()) return A.sign()>B.sign();
	A.val.reduct();B.val.reduct();
	return (A.sign()?(B.val<A.val):(A.val<B.val));
}
priority_queue < state > pq;
map < int , bool > vis;
int main(){
	// inputting start
	// 数据结构记得初始化! n,m别写反!
	scanf("%d%d%d",&n,&m,&k);
	int i,j;
	for(i=0;i<n;++i){
		scanf("%d",&j);
		if(j<0){
			neg.push_back(-j);
		}
		else{
			pos.push_back(j);
		}
	}
	#ifdef LOCAL
		TStart();
	#endif
	// calculation start
	// 数据结构记得初始化! n,m别写反!
	sort(pos.begin(),pos.end());
	sort(neg.begin(),neg.end());
	for(i=0;i<=m;i+=2){
		state tmp;
		for(j=0;j<i && j<(int)neg.size();++j){
			tmp.vn[tmp.cn++]=(int)neg.size()-1-j;
			tmp.val.multiply(neg[(int)neg.size()-1-j]);
		}
		if(j>=i){
			for(j=0;j<m-i && j<(int)pos.size();++j){
				tmp.vp[tmp.cp++]=(int)pos.size()-1-j;
				tmp.val.multiply(pos[(int)pos.size()-1-j]);
			}
			if(j>=m-i){
				vis[tmp.super_cdd()]=1;
				pq.push(tmp);
			}
		}
	}
	for(i=1;i<=m;i+=2){
		state tmp;
		for(j=0;j<i && j<(int)neg.size();++j){
			tmp.vn[tmp.cn++]=j;
			tmp.val.multiply(neg[j]);
		}
		if(j>=i){
			for(j=0;j<m-i && j<(int)pos.size();++j){
				tmp.vp[tmp.cp++]=j;
				tmp.val.multiply(pos[j]);
			}
			if(j>=m-i){
				vis[tmp.super_cdd()]=1;
				pq.push(tmp);
			}
		}
	}
	while(k--){
		state cur=pq.top();
		pq.pop();
		if(!k){
			cur.val.print(cur.sign());
			break;
		}
		for(i=0;i<cur.cp;++i){
			if(cur.editp(i,(cur.sign()?1:-1))){
				if(!vis[cur.super_cdd()]){
					vis[cur.super_cdd()]=1;
					pq.push(cur);
				}
				cur.editp(i,(cur.sign()?-1:1));
			}
		}
		for(i=0;i<cur.cn;++i){
			if(cur.editn(i,(cur.sign()?1:-1))){
				if(!vis[cur.super_cdd()]){
					vis[cur.super_cdd()]=1;
					pq.push(cur);
				}
				cur.editn(i,(cur.sign()?-1:1));
			}
		}
	}
	#ifdef LOCAL
		TReport();
	#endif
	return 0;
}

实践出真知2:

这题是我打算出出来的一个idea……先等我把它放到oj上再说吧qaq……

UPD:题被枪毙了,大家散了吧

还是丢一道想法类似的题吧:
CS Academy Round 79E

题意:

\(n\)(\(n \leq 10^5\))个整数(正负均可),问第\(k\)(\(k \leq 10^5\))小的子集和。

题解:

先把输入的数升序排列记作数组\(A\)

状态意义:

状态\(\{S,C\}\)表示当前选择的子集和是\(S\),最后一个选择的数下标为\(C\)

状态间的大小关系:

即子集和\(S\)的大小关系。

起始状态\(Source\)

使用多个起始状态。\(Source_i=\{ A_0+A_1+...+A_i , i \}\)

连边方案:

\(\{S,C\}\)\(\{S+A[C+1]-A[C],C+1\}\)连边,边权为\(A_{C+1}-A_C\)
(可以看作把最后一个选取的数往后推了一位)

正确性证明:

因为\(A_{C+1}-A_C \geq 0\),所以边权非负。又易见每一个状态都恰由一个\(Source\)可达,故做法正确。

勘误(来自2020/2/2)

上述算法是假的。
正确的算法思路相近,详见这里
核心思路是,允许两种操作:将末尾元素向后推一位、在末尾元素后面再加一个元素。为了保证第二种操作边权非负,需对答案集合的负数子集取补集。

总结:

抱歉我扯不出什么总结来了,东西都在正文讲完了,大家散了吧

posted @ 2018-05-06 07:54  TianyiQ  阅读(967)  评论(0)    收藏  举报