Loading

如何科学地判断负环和求有负权图的最短路(几乎完更)

主要内容复读自 EI 的博客。部分符号可能会有不同。

起因是这样的,想做一道 0/1 分数规划+判负环的,结果发现出题人卡了 bfs 判负环,放了 dfs 判负环过。遂到谷群求教 dfs 判负环的有关事宜,结果被群友教育了,于是就来学了qwq。
(虽然说这个算法只能做整数边权的图,做不了上面那道题

目标:\(O(m\sqrt{n}\log|w|)\) 时间内对整数边权的有向图判断负环/求出最短路。

1. 前置知识

首先我们要对原来的图进行一些改造。

1. Pricing

如果你学过 Johnson 全源最短路 的话,那你应该知道这样一个操作:
如果能对每个点 \(i\) 求出点权 \(p(i)\),使得 \(\forall u,v:w_p(u,v):=w(u,v)+p(u)-p(v)\geqslant0\),那么我们可以更新边权为 \(w_p(u,v)\),然后在新的图上跑 dijkstra。

如果你不会这个操作,这里是对该方法的简单解释。 注意到对于任意一条 $s$ 到 $t$ 的路径都有 $w_p(s,v_0)+w_p(v_0,v_1)+\cdots+w_p(v_k,t)=w(s,v_0)+w(v_0,v_1)+\cdots+w(v_k,t)+p(s)-p(t)$,由 $p(s)-p(t)$ 是定值知该操作对最短路径没有影响。又有保证了 $\forall u,v:w_p(u,v)\geqslant0$,于是我们就能使用 dijkstra 算法了。

另外如果进行了一些赋权操作后出现了负环,那么原来的图也必定存在负环(证明方法同上,因为是环路,所有的 \(p\) 值都消掉了)

很不幸的是在 Johnson 全源最短路中,我们的做法是先跑了一遍 SPFA,所以这个求法对我们的问题没有帮助。

2. Scaling

但是如果图的边权是整数,那么如果你学过弱多项式复杂度费用流,你会知道这样的操作:
\(w'(u,v):=\left\lceil\dfrac{w(u,v)}{2}\right\rceil\),则 \(w(u,v)=\begin{cases}2w'(u,v),&w(u,v)\bmod2=0\\2w'(u,v)-1,&w(u,v)\bmod2=1\end{cases}\)
此时我们先对调整后的边权求出权值 \(p'(i)\),那么我们令 \(p(i)=2p'(i)\),则容易发现 \(w_p(u,v)=w(u,v)+p(u)-p(v)\geqslant-1\)
此时我们就将问题转化为了在 \(O(m\sqrt{n})\) 时间内调整 \(p(i)\) 使 \(\forall u,v:w_p(u,v)\geqslant0\)

对上述方法的一些补充 Q1: 你这个操作会不会导致误判负环啊。

A1: 首先如果直接将边权除以二什么都不会影响。对于奇数我们做上取整相当于将边权增加了 \(0.5\),那么只要原来的图里没有负环,新的图里也不会有负环。这样我们就不会出现误判负环的情况。

Q2: 你没说上面算法的递归终点啊。

A2: 事实上上述操作的目的是将边权都变成非负的。注意到 \(\left\lceil-\dfrac{1}{2}\right\rceil=0\),事实上进行足够多次上述操作后负数都会变成 \(0\),这时图中所有的边权都是非负的,令 \(p(i)\equiv0\) 即可。

然后你还需要一些离散数学知识。

3. Dilworth 定理相关

温馨提示:该部分内容较多。建议先粗略阅读,学完算法后再来这里看具体细节。

Dilworth 定理参考自北大的《组合数学》。后面的那个定理是我自己证的。

实际上这部分内容就有些剧透的意味了。

直接上定理:
设有限偏序集 \((X,\leqslant)\) 的宽度为 \(m\),则存在划分 \(X=\bigcup\limits_{i=1}^{m}C_i\),使得每一个 \(C_i\) 都是链。

我们可以用 OI 中的东西来简单理解:
考虑一个有限 DAG。在这里我们定义 DAG 中的链为一个点集,其中任意取出两个点 \(x,y\),必有 \(x\) 能到达 \(y\)\(y\) 能到达 \(x\)(也就是说链可以使不连续的)。相对地,我们定义 DAG 中的反链为一个点集,其中从该点集中的任何一个点经该 DAG 都不能到达点集中的其他任何一个点。特别地,单独的一个点既是链又是反链。我们定义这条反链的长度为这个点集的大小,同时定义一个 DAG 的宽度为这个 DAG 的最长反链的长度。

然后这个定理的意思就是宽度为 \(m\) 的 DAG 一定可以被划分成不相交的 \(m\) 个链。(实际上划分就已经有不相交的意思了)

事实上还有一个类似的定理,不过证明起来简单地多。 定理:我们称 DAG 的高度为 DAG 中最长链的长度,则高度为 $n$ 的 DAG 能被划分为 $n$ 个反链。

证明:我们每次取出 DAG 中所有出度为 \(0\) 的点,这些点显然构成了一个反链。并且因为 DAG 中所有的极大链(即不存在其他的链使这条链包含于它)都一定包含一个出度为 \(0\) 的点,而最长链一定是极大链,故操作后的 DAG 的高度一定为 \(n-1\)。重复上述操作即得结论。

我们来证明一下 Dilworth 定理。不妨还是用 DAG 的语言进行叙述。
考虑使用数学归纳法。

  • 如果 DAG 中只有 \(1\) 个点,显然成立。
  • 假设原结论对有 \(1,2,\cdots,k-1\) 个点的 DAG 都成立,我们来考虑有 \(k\) 个点的 DAG。首先我们取一条极大链(即不存在其他的链使这条链包含于它) \(C\),我们考虑去掉这个链后的 DAG。
    • 若去掉这个链后的 DAG 的宽度是 \(m-1\),那么原结论成立。
    • 否则,去掉链后的 DAG 的宽度一定是 \(m\)。这个时候我们考虑去掉链后的 DAG 的一个反链 \(A=\{a_1,a_2,\cdots,a_m\}\)。然后,我们定义 \(S^-\) 为原 DAG 上能到达 \(A\) 中的点的集合, \(S^+\) 为从 \(A\) 中的点经原 DAG 上的边能到达的点的集合。特别地,\(A\) 中的点同属两个集合。容易发现 \(S^-\cup S^+\) 就是原 DAG,并且 \(S^-\cap S^+=A\)。这个时候我们就可以根据归纳假设,将 \(S^-\)\(S^+\) 都划分为 \(m\) 个链。记包含了 \(a_i\) 的链为 \(S_i^-\)\(S_i^+\),容易发现将 \(S_i^-\)\(S_i^+\) 通过 \(a_i\) 拼成一条新链,原来的 DAG 就被划分为了 \(m\) 个链。
关于证明中的两个“容易发现”

1.1 反证,首先 \(S^-\cup S^+\) 中的点一定在原 DAG 上。若存在原 DAG 上的点 \(x\) 满足 \(x\notin S^-\)\(x\notin S^+\),则 \(A\cup \{x\}\) 是一个反链,有原 DAG 的宽度至少为 \(m+1\),矛盾。故 \(S^-\cup S^+\) 就是原 DAG。

1.2 首先 \(A\) 中的点一定属于 \(S^-\cap S^+\)。另一方面,设 \(x\in S^-\cap S^+\),则存在 \(a_i,a_j\) 满足 \(x\) 能到达 \(a_i\)\(a_j\) 能到达 \(x\),故 \(a_j\) 能到达 \(a_i\)。因为 \(A\) 是反链,故只能有 \(a_i=a_j\),即 \(x\) 能到达 \(a_i\)\(a_i\) 能到达 \(x\)。又因为这是个 DAG,于是必有 \(x=a_i\),即 \(x\in A\)。故 \(S^-\cap S^+=A\)

2 我们来证明 \(a_i\)\(S_i^-\) 中出度为 \(0\)。反证,若存在 \(x\in S_i^-\) 使 \(a_i\) 能到达 \(x\)\(x\) 不为 \(a_i\)。由 \(S^-\) 定义知存在 \(a_j\) 使 \(x\) 能到达 \(a_j\)。故 \(a_i\) 能到达 \(a_j\)\(a_i\neq a_j\),这与 \(A\) 是反链矛盾。同理可证 \(a_i\)\(S_i^+\) 中入度为 \(0\)。由此即得 \(S_i^-\cup S_i^+\) 是链。

放一张图(我自己画的)直观理解:
如左图,这是一个宽度为 \(3\) 的 DAG。我们第一步操作是选出了一条极大链(红色),然后发现去掉极大链后的 DAG 宽度还是 \(3\)。于是我们从去掉极大链后的 DAG 中取出一条反链(绿色),并定义 \(S^-\)(蓝色)和 \(S^+\)(黄色)。
然后看右图,根据归纳假设将两个集合分别划分成 \(3\) 条链,然后这些链可以通过中间的反链连接起来成为一条长链,于是原来的 DAG 就被划分为了 \(3\) 条链。

image

有了 Dilworth 定理,我们就可以证明下面的结论了:
对于一个结点数为 \(mn\) 的 DAG,必然存在长度至少为 \(m\) 的链或长度至少为 \(n\) 的反链。

证明:假设最长链的长度为 \(h\),最长反链的长度为 \(w\)。反证,假设有 \(h<m\)\(w<n\)。根据 Dilworth 定理,我们可以将原 DAG 划分成 \(w\) 个链 \(C_1,C_2,\cdots,C_w\)。计算 DAG 中的结点个数 \(S\),得 \(S=|C_1|+|C_2|+\cdots+|C_w|\leqslant hw<mn\),矛盾。原命题得证。

显然的推论:结点数为 \(n\) 的 DAG 中链和反链的最大长度必然有一个大于等于 \(\sqrt{n}\)

2. 算法

我们只考虑原图中边权是 \(-1\)\(0\) 的部分。严谨地说,我们记原图为 \(G=(V,E)\),定义 \(E_p=\{(u,v)|w_p(u,v)\leqslant0\},\ V_p=\{v|\exist u:(u,v)\in E'\or(v,u)\in E'\},\ G_p=(V_p,E_p)\)

然后我们对 \(G_p\) 进行缩点。容易发现若 \(G_p\) 的一个 SCC 中存在负边,那么一定会存在一个负环。若原图不存在负环,那么每个 SCC 中的边边权都是 \(0\),不会产生任何影响,直接缩掉就行了。这样我们将 \(G_p\) 转化为了较容易处理的 DAG。

我们发现 DAG 有这样的一个性质:
取一个点集 \(S\),使 \(\forall s\in S\)\(s\)\(G_p\) 可达的点都在 \(S\) 中。这样由 \(S\) 连出的边的边权都是正的。
然后我们将 \(S\) 内所有点的 \(p\) 值减 \(1\)。容易发现这个操作有两点影响:

  • 我们成功地将连入 \(S\) 的负边变成了 \(0\) 边。
  • 我们成功地让 S 连出的边权是 1 的边变成了 0 边,于是上面取子图缩点的步骤要从头再来。

好吧至少我们没有产生新的负边。
容易发现如果我们定义 \(S_u\) 为点 \(u\) 先经过一条 \(-1\) 边,再经 \(G_p\) 能到达的所有点(不包含 \(u\) 本身),那么我们就处理掉了一个点连出的所有负边。但这样做复杂度将高达 \(O(nm)\)。于是我们考虑一次处理掉多个点。

对于上面的方法,我们先对这个 DAG 定义一个图 \(K=(V_K,E_K)\)\(V_K\) 由所有含有边权为 \(-1\) 的出边的点组成。\(E_K\) 的定义是这样的:对于 \(u,v\in V_K\),若有 \(u\) 能先经过一条边权为 \(-1\) 的边,再通过 DAG 能够到达 \(v\),则我们由 \(u\)\(v\) 连一条边。那么接下来我们就对 \(K\) 进行操作。

进行具体操作之前我们先来画一下图,明确一下定义。

image

如图是缩完点后的 \(G_p\),这是个 DAG。我们在这个 DAG 上定义了图 \(K\)(红色的点和边)。很明显 \(K\) 也是个 DAG。但是这张图里并没有画出 \(K\) 所有的边,这些红边只是描述了 \(K\) 上的可达性。也就是说在上面的图里,如果能从一个红点经红边到达另一个红点,那么在真正的 \(K\) 中这两个点之间就会连一条边。

然后我们就可以对 \(K\) 进行操作了。我们的方法是对 \(K\) 中的链和反链分别设计算法。

先看简单的反链。
考虑一个反链 \(A\),根据定义,从其中任意取出两点 \(u,v\),有 \(u\notin S_v\)\(v\notin S_u\)。此时直接对 \(S:=\bigcup\limits_{u\in A}S_u\) 进行操作,那么从这条反链连出的所有的负边就都被处理掉了。容易发现复杂度是简单的 \(O(m)\)

再来看链的情况。
考虑链 \(C=\{u_1,u_2,\cdots,u_n\}\)。根据性质我们不妨设 \(\forall 1\leqslant i\leqslant n-1: u_{i+1}\in S_{u_i}\)。我们的方法是,考虑定义这样一个点集 \(T_u\),它由。然后依次对 \(T_{u_n},T_{u_{n-1}},\cdots,T_{u_1}\) 进行操作。
首先我们容易发现一开始有 \(T_{u_n}\subset T_{u_{n-1}}\subset\cdots\subset T_{u_1}\)
但是,值得注意的是,在我们执行了若干个操作后,尚未执行操作的点集会发生变化。但仍然满足上述包含关系。

链和反链上的算法都讲完了,现在让我们把它们合起来。
前面说过结点数为 \(n\) 的 DAG 中链和反链的最大长度必然有一个大于等于 \(\sqrt{n}\)。我们只需要 \(O(m)\) 在 DAG 上做一次 dp 求每个点连出的最长链(这个是非常简单的),然后检查极大点构成的反链。这样我们就能一次处理掉根号个数的点。

事实上如果我们这么做的话,算法在根号次数内就会完成(这个证明我就不会了qaq

蒯来了 EI 的证明 经过 $\sqrt{|K|}$ 轮,集合大小每轮乘以一个不超过 $1 - \frac 1{\sqrt{|K|}}$ 的因子,那么集合变成原来的 $\leq 1/\mathrm{e}$ 倍,那么经过 $$\sum_{i\geq 0} \sqrt{\frac{|K|}{\mathrm{e}^i}} \leq \frac 1{1 - \mathrm {e}^{-1/2}} \sqrt{|K|}$$ 轮之后一定完成。

3. 代 码 实 现

如果前面的理论部分还有细节不清楚,我们最后再讲解一下 EI 的实现。

#include <bits/stdc++.h>

using namespace std;

const int V = 5005, E = 5005;

int N, M;
int u[E], v[E], w[E];

int p[V];
bool mark[V];
vector<pair<int, int>> G[V], G2[V];

int t, N2;
stack<int> stk;
int dfn[V], low[V], belongs[V], topo[V];
void tarjan(int u) {
	dfn[u] = low[u] = ++t; stk.push(u);
	for (const auto &e : G[u]) {
		int v = e.first, w = e.second + p[u] - p[v];
		if (w > 0) continue;
		if (!dfn[v]) {
			tarjan(v);
			low[u] = min(low[u], low[v]);
		} else if (belongs[v] == -1) {
			low[u] = min(low[u], dfn[v]);
		}
	}
	if (dfn[u] == low[u]) {
		topo[N2++] = u;
		int v;
		do {
			v = stk.top(); stk.pop();
			belongs[v] = u;
		} while (v != u);
	}
}

bool SCC() {
	t = N2 = 0;
	memset(dfn, 0, sizeof(dfn));
	memset(mark, 0, sizeof(mark));
	memset(belongs, -1, sizeof(belongs));
	for (int i = 0; i != N; ++i) G2[i].clear();

	for (int i = 0; i != N; ++i) if (!dfn[i]) tarjan(i);
#ifdef ELEGIA_DEBUG
	int marks = 0;
#endif // ELEGIA_DEBUG
	for (int i = 0; i != N; ++i) for (const auto &e : G[i]) {
		int j = e.first, w = e.second + p[i] - p[j];
		if (w < 0) {
			if (belongs[i] == belongs[j]) return true;
			mark[belongs[i]] = true;
#ifdef ELEGIA_DEBUG
			++marks;
#endif // ELEGIA_DEBUG
		}
		G2[belongs[i]].emplace_back(belongs[j], w);
	}
#ifdef ELEGIA_DEBUG
	cerr << "marks: " << marks << endl;
#endif // ELEGIA_DEBUG
	return false;
}

int vis[V], pref[V], sec[V], len[V], chn[V];
vector<int> delay[V];
void chain(int start) {
	int len = 0;
	while (start != -1) {
		chn[++len] = start;
		start = sec[start];
	}
	memset(vis, 0, sizeof(vis));
	queue<int> q;
	for (int i = len; i; --i) {
		int s = chn[i];
		for (const auto &e : G2[s]) {
			int v = e.first, w = e.second;
			if (vis[v] || w >= 0) continue;
			vis[v] = i;
			q.push(v);
		}
		for (int u : delay[i]) if (!vis[u]) {
			vis[u] = i;
			q.push(u);
		}
		delay[i].clear();
		while (!q.empty()) {
			int u = q.front(); q.pop();
			for (const auto &e : G2[u]) {
				int v = e.first, w = e.second;
				if (w > 0) {
					if (w < i) delay[i - w].push_back(v);
					continue;
				}
				if (vis[v]) continue;
				vis[v] = i;
				q.push(v);
			}
		}
	}
}

bool refine() {
	while (true) {
		{
			bool flag = false;
			for (int u = 0; u != N; ++u) {
				for (const auto &e : G[u]) {
					int v = e.first, w = e.second + p[u] - p[v];
					if (w < 0) {
						flag = true;
						break;
					}
				}
				if (flag) break;
			}
			if (!flag) return false;
		}
		if (SCC()) return true;
		memset(len, 0, sizeof(len));
		memset(pref, -1, sizeof(pref));
		memset(sec, -1, sizeof(sec));
		int longest = -1;
		for (int i = 0; i != N2; ++i) {
			int u = topo[i];
			for (const auto &e : G2[u]) {
				int v = e.first, w = e.second;
				if (w > 0) continue;
				if (pref[v] != -1 && len[v] > len[u]) {
					len[u] = len[v];
					pref[u] = pref[v];
				}
			}
			if (mark[u]) {
				int opt = 0;
				for (const auto &e : G2[u]) {
					int v = e.first, w = e.second;
					if (w >= 0) continue;
					if (pref[v] != -1 && len[v] > opt) {
						opt = len[v];
						sec[u] = pref[v];
					}
				}
				if (++opt > len[u]) {
					len[u] = opt;
					pref[u] = u;
				} else sec[u] = -1;
			}
			if (longest == -1 || len[u] > len[longest])
				longest = u;
		}
		int size = 0;
		queue<int> q;
		memset(vis, 0, sizeof(vis));
		for (int i = N2 - 1; i >= 0; --i) {
			int u = topo[i];
			if (vis[u] || !mark[u]) continue;
			++size;
			for (const auto &e : G2[u]) {
				int v = e.first, w = e.second;
				if (!vis[v] && w < 0) {
					vis[v] = true;
					q.push(v);
				}
			}
			while (!q.empty()) {
				int u = q.front(); q.pop();
				for (const auto &e : G2[u]) {
					int v = e.first, w = e.second;
					if (!vis[v] && w <= 0) {
						vis[v] = true;
						q.push(v);
					}
				}
			}
		}
#ifdef ELEGIA_DEBUG
		cerr << "antichain: " << size << ", chain: " << len[longest] << endl;
#endif // ELEGIA_DEBUG
		if (size < len[longest]) chain(longest);
		for (int i = 0; i != N; ++i)
			p[i] -= vis[belongs[i]];
	}
}

bool priceScaling() {
	int L = 0;
	for (int i = 0; i != M; ++i)
		while ((1 << L) <= -w[i]) ++L;
	memset(p, 0, sizeof(p));

	for (int d = L - 1; d >= 0; --d) {
		for (int i = 0; i != N; ++i) G[i].clear();
		for (int i = 0; i != N; ++i) p[i] *= 2;
		for (int i = 0; i != M; ++i)
			G[u[i]].emplace_back(v[i], -((-w[i]) >> d));
		if (refine()) return true;
#ifdef ELEGIA_DEBUG
		cerr << "finish scaling 2^" << d << endl;
#endif // ELEGIA_DEBUG
	}
	return false;
}

int main() {
	ios::sync_with_stdio(false);
	cin.tie(nullptr); cout.tie(nullptr);

	cin >> N >> M;
	for (int i = 0; i != M; ++i) {
		cin >> v[i] >> u[i] >> w[i];
		--u[i]; --v[i];
	}
	if (priceScaling()) {
		cout << "NO\n";
	} else for (int i = 0; i != N; ++i)
		cout << p[i] << " \n"[i == N - 1];

	return 0;
}

posted @ 2022-07-16 19:12  pjykk  阅读(58)  评论(0编辑  收藏  举报