Dijkstra算法与SPFAlyh's blog

Dijkstra算法

        Dijkstra算法是一种基于贪心思想的求单源最短路算法，时间复杂度为$O(n^2)$。Dijkstra算法要求边权为正。
步骤:
        ① 一个数组dis[]保存最短路径长的结果，源点值置为0，其余为inf。
        ② 一个集合S保存所有已求出最短路径长的节点，初始化S为空。
        ③ 从不在S中的节点集合中找到dis最小的节点，用这个节点更新与之相邻的的节点dis值(这个过程叫松弛)，并将这个节点加入S。
        重复③，直到所有节点都加入S中。

        我们证明，不在S中的节点集合里面dis最小的节点，设为min，它的dis最小值必为当前值。这是因为若存在另外一条路径使得点dis值更小，必定会经过其余的某一个点，假设该点dis值为k，则应该有k+p < min，然而p ≥ 0且k ≥ min，可知这是不可能的。从这个原理上就可以明白Dijkstra要求边权非负的原因。
        在实际操作中可以用小根堆来维护最小值，这样时间复杂度降至稳定的$O((m+n)logn)$。这里值得注意的是，由于堆中经常进行调整，对于一个dis值已经被修改然而已在堆中的节点，我们应该对其进行调整。所以构造一个从节点序号到堆中具体位置的映射，来更新节点在堆中的位置。具体细节见示例代码。
        本示例代码采用链式前向星+小根堆优化的方法。

#include<iostream>
#include<cstring>

using namespace std;
int n, m, from;
const int MAX = 100000;
struct edge {
    int to, value, next;
} op[2 * MAX + 1];//链式前向星存图
int head[MAX + 1], dis[MAX + 1];
int heap[MAX + 1], size = 0;//小根堆使用
int ID[MAX + 1] = {0};//构造映射

inline void solve(int x) {
    int x1 = x << 1, x2 = x << 1 | 1, minn = x;
    if (x1 <= size && dis[heap[x1]] < dis[heap[minn]])minn = x1;
    if (x2 <= size && dis[heap[x2]] < dis[heap[minn]])minn = x2;
    if (minn != x) {
        swap(ID[heap[minn]], ID[heap[x]]);
        swap(heap[minn], heap[x]);
        solve(minn);
    }
}

inline int top() {
    int r = heap[1];
    heap[1] = heap[size--];
    ID[heap[1]] = 1;
    solve(1);
    return r;
}

inline void up(int x) {
    if (x == 1)return;
    if (dis[heap[x / 2]] > dis[heap[x]])
        ID[heap[x / 2]] = x, ID[heap[x]] = x / 2, swap(heap[x / 2], heap[x]), up(x / 2);
}

inline void add(int x) {
    heap[++size] = x;
    ID[x] = size;
    up(size);
}

int main() {
    ios::sync_with_stdio(false);
    cin >> n >> m >> from;
    memset(head, -1, sizeof(head));
    for (register int i = 1; i <= n; i++)dis[i] = 2147483647;
    dis[from] = 0;
    for (register int i = 1; i <= m; i++) {
        int x, y, z;
        cin >> x >> y >> z;
        op[i].next = head[x], op[i].to = y, op[i].value = z;
        head[x] = i;
    }
    add(from);//将源点加入堆
    while (size > 0) {
        int t = top();//取堆中dis值最小节点,并将其移出堆,此后这个节点不可能再被加入堆中
        for (register int i = head[t]; i != -1; i = op[i].next) {//找到该节点邻接节点
            if (dis[t] + op[i].value < dis[op[i].to]) {
                dis[op[i].to] = dis[t] + op[i].value;//修改该邻接点dis值
                if (ID[op[i].to] != 0) {//已在堆中,更新其在堆中的位置
                    up(ID[op[i].to]);
                } else {
                    add(op[i].to);//不在堆中,将其加入堆
                }
            }
        }
    }
    for (register int i = 1; i <= n; i++)cout << dis[i] << " ";
    return 0;
}

SPFA

        接下来介绍SPFA，SPFA本质上是基于队列的Bellman-Ford算法，支持负边权并且可以判断负环。该算法时间复杂度达$O(km)$。该算法不稳定，在某些情况下时间复杂度会达到朴素Bellman-Ford级别。
算法思想:
        ① 将源点加入队列
        ② 弹出队首元素，用该元素松弛其邻接点，若松弛后的节点不在队列中，让该节点入队。
        重复②，直到队列为空为止。
        当图存在负环时，最短路径是不存在的。可以证明存在负环的充要条件是某一个节点入队次数大于n次。利用这个结论可以方便地判定负环，这是BFS版的SPFA判负环方法。SPFA还有DFS版本，即用栈代替队列，只要需要松弛的点在栈中就可以判断有负环。

#include<iostream>
#include<cstring>
#include<queue>

using namespace std;
int n, m, from;
const int MAX = 100000;
struct edge {
    int to, value, next;
} op[2 * MAX + 1];//链式前向星存图
int head[MAX + 1];
int dis[MAX + 1];
bool vis[MAX + 1] = {false};
int num[MAX + 1] = {0};

int main() {
    ios::sync_with_stdio(false);
    queue<int> que;
    cin >> n >> m >> from;
    memset(head, -1, sizeof(head));
    for (register int i = 1; i <= n; i++)dis[i] = 2147483647;
    dis[from] = 0;
    for (register int i = 1; i <= m; i++) {
        int x, y, z;
        cin >> x >> y >> z;
        op[i].next = head[x], op[i].to = y, op[i].value = z;
        head[x] = i;
    }
    que.push(from), vis[from] = true, num[from]++;
    while (!que.empty()) {
        int r = que.front();
        que.pop();
        vis[r] = false;
        for (int i = head[r]; i != -1; i = op[i].next) {
            if (dis[r] + op[i].value < dis[op[i].to]) {
                dis[op[i].to] = dis[r] + op[i].value;
                if (!vis[op[i].to])que.push(op[i].to), vis[op[i].to] = true, num[op[i].to]++;
                if (num[op[i].to] > n) {
                    cout << "ERR" << endl;
                    return 0;
                }
            }
        }
    }
    for (register int i = 1; i <= n; i++)cout << dis[i] << " ";
    return 0;
}

        经过洛谷上的测试，Dijkstra+heap优化可以通过P3371和强化版的P4779。SPFA可以通过前者，但由于在后者被卡，无法通过。
        另外，Dijkstra+heap与SPFA的速度比较存在争议，对上面算法的正确性和算法选择请读者自行理解。推荐在没有负边时使用堆优化Dijkstra算法，仅在含负边或者需要判负环时使用SPFA。
        关于SPFA的两个优化：

SLF（Small Label First）优化：设当前队首元素为f（要取出队首），待入队元素为x，若dict[x]<dict[f]，则将x插入队首，否则插入队尾。
LLL（Large Label First）优化：若队首为f，如果dict[f]大于当前队列中dict的平均值，则取出f放到队尾，跳过该元素判断下一个，直到队首元素的dict不大于平均值，将其拿出进行松弛。实际操作中我们需要不停记录队列中的元素个数和dict和来维护平均值。

对最短路算法的进一步理解

        Dijkstra算法中，每一个结点仅入堆一次，这样可以保证时间复杂度在较低的水平，并且能够保证算法的稳定性。但是如果我们每次取的点并不一定是最小dist的点又会如何？
        倘若我们不加堆，每次取一个随机的结点，这样最终也是可以完成松弛操作，同样可以求出最短路。但是那样结点会多次入堆，造成效率的大大降低。
        现在来看SPFA的算法过程，可以发现SPFA就如同随机从点中取一个进行松弛，只不过加了一个队列而已。这样SPFA的效率应该比Dijkstra堆优化版本更低。但是在某些随机的情况下，SPFA可以做到每一个结点仅入队一次或很少次，这样做到了和Dijkstra相似的算法过程，还省略了堆优化的log级别的复杂度，总体效率会优于Dijkstra。但是在某些构造数据下，SPFA算法过程中结点会多次入队，直接被卡到甚至$O(nm)$级别。但相同情况下，Dijkstra通过加堆优化，强行先松弛dist最小的结点的出点，保证结点入堆次数仅有一次，实现了算法稳定性的提升。
        到这里就可以更深层次地理解Dijkstra要求边权非负的原因：利用边权非负来保证结点入堆次数仅有一次以提高效率。
        对于有负权的图如何应用Dijkstra算法？一种方法是用Johnson重标记法来解决这个问题，以后再介绍。
        对于上面的问题和讨论，可以发现Dijkstra其实也可以处理不含负环的负边权问题，只需允许结点多次入堆即可，代码如下：

//该代码未经测试!!!!!
#include<iostream>
#include<cstring>

using namespace std;
int n, m, from;
const int MAX = 100000;
struct edge {
    int to, value, next;
} op[2 * MAX + 1];
int head[MAX + 1], dis[MAX + 1];
int heap[MAX + 1], size = 0;
int ID[MAX + 1] = {0};

inline void solve(int x) {
    int x1 = x << 1, x2 = x << 1 | 1, minn = x;
    if (x1 <= size && dis[heap[x1]] < dis[heap[minn]])minn = x1;
    if (x2 <= size && dis[heap[x2]] < dis[heap[minn]])minn = x2;
    if (minn != x) {
        swap(ID[heap[minn]], ID[heap[x]]);
        swap(heap[minn], heap[x]);
        solve(minn);
    }
}

inline int top() {
    int r = heap[1];
    ID[r] = 0;//映射置为0表示可以重新加入
    if (size > 1)heap[1] = heap[size--], ID[heap[1]] = 1, solve(1);//size>1更新信息
    else size = 0;//否则size=0
    return r;
}

inline void up(int x) {
    if (x == 1)return;
    if (dis[heap[x / 2]] > dis[heap[x]])
        ID[heap[x / 2]] = x, ID[heap[x]] = x / 2, swap(heap[x / 2], heap[x]), up(x / 2);
}

inline void add(int x) {
    heap[++size] = x;
    ID[x] = size;
    up(size);
}

int main() {
    ios::sync_with_stdio(false);
    cin >> n >> m >> from;
    memset(head, -1, sizeof(head));
    for (register int i = 1; i <= n; i++)dis[i] = 2147483647;
    dis[from] = 0;
    for (register int i = 1; i <= m; i++) {
        int x, y, z;
        cin >> x >> y >> z;
        op[i].next = head[x], op[i].to = y, op[i].value = z;
        head[x] = i;
    }
    add(from);
    while (size > 0) {
        int t = top();
        for (register int i = head[t]; i != -1; i = op[i].next) {
            if (dis[t] + op[i].value < dis[op[i].to]) {
                dis[op[i].to] = dis[t] + op[i].value;
                if (ID[op[i].to] != 0) {
                    up(ID[op[i].to]);
                } else {
                    add(op[i].to);
                }
            }
        }
    }
    for (register int i = 1; i <= n; i++)cout << dis[i] << " ";
    return 0;
}

分层图最短路

考虑一个问题，现在给定一个图，允许将其中不超过k条边改为0边权，问某点的单源最短路最小为多少？这就是分层图最短路问题，洛谷有模板题。
方法是改进后的Dijkstra算法。规定dict[i][j]为源点到i且将j条边权值置为0时的最短路，根据此进行松弛。这个算法在Dijkstra算法基础上给dict添加了第二个指标元素来描述以更改边权边的条数。类似的方法也可以处理将边权改成其它值的情况。

#include<bits/stdc++.h>

#define inf 1l<<50
using namespace std;

inline int read() {
    char e = getchar();
    int s = 0;
    while (e < '-')e = getchar();
    while (e > '-')s = (s << 1) + (s << 3) + (e & 15), e = getchar();
    return s;
}

struct Edge {
    int next, to, v;
} edge[50005 << 1];

struct Node {
    int r, k;

    Node(int x = 0, int y = 0) : r(x), k(y) {}
} heap[10005 * 24];

int n, m, head[10005], cnt = 1, k, size = 0, rk[10005][24];
long long dp[10005][24];

inline void add(int x, int y, int z) {
    edge[cnt].to = y, edge[cnt].next = head[x], edge[cnt].v = z, head[x] = cnt++;
}

void solve(int x) {
    int x1 = x << 1, x2 = x << 1 | 1, minn = x;
    if (x1 <= size && dp[heap[x1].r][heap[x1].k] < dp[heap[minn].r][heap[minn].k])minn = x1;
    if (x2 <= size && dp[heap[x2].r][heap[x2].k] < dp[heap[minn].r][heap[minn].k])minn = x2;
    if (x != minn)
        swap(heap[minn], heap[x]), swap(rk[heap[minn].r][heap[minn].k], rk[heap[x].r][heap[x].k]), solve(minn);
}

void up(int x) {
    if (x == 1)return;
    if (dp[heap[x >> 1].r][heap[x >> 1].k] > dp[heap[x].r][heap[x].k]) {
        swap(heap[x >> 1], heap[x]), swap(rk[heap[x >> 1].r][heap[x >> 1].k], rk[heap[x].r][heap[x].k]), up(x >> 1);
    }
}

Node top() {
    Node t = heap[1];
    rk[t.r][t.k] = 0, rk[heap[size].r][heap[size].k] = 1, heap[1] = heap[size--], solve(1);
    return t;
}

void add(Node z) {
    rk[z.r][z.k] = ++size, heap[size] = z, up(size);
}

int main() {
    n = read(), m = read(), k = read();
    for (int i = 1; i <= m; i++) {
        int x = read(), y = read(), z = read();
        add(x, y, z), add(y, x, z);
    }
    for (int i = 1; i <= n; i++)
        for (int j = 0; j <= k; j++)dp[i][j] = inf;
    for (int i = 0; i <= k; i++)dp[1][i] = 0, add(Node(1, i));
    while (size > 0) {
        Node t = top();
        for (int i = head[t.r]; i; i = edge[i].next) {
            if (dp[edge[i].to][t.k] > dp[t.r][t.k] + edge[i].v) {//第一种松弛
                dp[edge[i].to][t.k] = dp[t.r][t.k] + edge[i].v;
                if (rk[edge[i].to][t.k])up(rk[edge[i].to][t.k]);
                else add(Node(edge[i].to, t.k));
            }
            if (t.k < k && dp[edge[i].to][t.k + 1] > dp[t.r][t.k]) {//第二种松弛
                dp[edge[i].to][t.k + 1] = dp[t.r][t.k];
                if (rk[edge[i].to][t.k + 1])up(rk[edge[i].to][t.k + 1]);
                else add(Node(edge[i].to, t.k + 1));
            }
        }
    }
    long long ans = inf;
    for (int i = 0; i <= k; i++)ans = min(ans, dp[n][i]);
    cout << ans;
    return 0;

}

关于分层图的理解：将原图赋值若干份，平行排布，下一层与上一层的点之间有边（与原图边对应）相连，边权为0。这样分层图就有了所有的决策情况，求一个最短路即可。分层图模型在很多问题上都有应用。