「Bzoj 1706」「usaco2007 Nov」relays 奶牛接力跑 (矩阵快速幂+图上DP+离散化/倍增Floyd)

BZOJ 1706
Luogu P2886
无向连通图，求$S$到$E$经过$k$条边的最短路

矩阵经典问题，看这里

图中DP：设$dp(i,x)$为走了$i$步在$x$点位置，则转移$dp(i,x)=\sum dp(i-1, y)+w(y, x)$，其中$(x,y) \in E$
同上面不好算。
我们定义$dp(k,i,j)$为$i$到$j$经过$k$个点(走$k$步)的最短路。则
$$dp(k,i,j)=min(dp(k-1,i,k)+dp(k-1,k,j))$$
压掉第一维，则$dp(i,j)=min(dp(i,k)+dp(k,j))$
我们发现和上面式子很像，那么我们可以试图修改一下矩阵乘法定义如下

$$C[i,j]=min(A[i,k]+B[k,i])$$

易推导结合律。可以用矩阵快速幂加速。

其实这也体现出 Floyd 的倍增思想。相当于对一个图做几次 Floyd 就是走了几步。矩阵快速幂就是倍增的过程。
这题在图上DP思想实际上就是分层图思想，每一步对应一张新图。与 NOIP2017Day1T3 有点像

知识点：
1、图上的 DP
2、矩阵快速幂(乘法)运用(修改)
3、矩阵快速幂的矩阵大小不要太大，尽可能离散化
4、倍增 Floyd 思想。
5、矩阵乘法修改后单位矩阵记得改

#include<cstdio> 
#include<cstring>
#include<algorithm>
#include<iostream>
#include<vector>
#include<set>
#include<cmath>
#define ms(i, j) memset(i, j, sizeof i)
#define LL long long
#define db double
#define fir first
#define sec second
#define mp make_pair
using namespace std;

namespace flyinthesky {

    const int MAXN = 100 + 10;

    struct matrix {
        int a[MAXN][MAXN];
        matrix() {ms(a, 60);} 
    }whw;

    int n, k, m, s, e, G[MAXN][MAXN], num[1005], sz;

    int gt(int x) {return num[x] ? num[x] : num[x] = ++sz;}

    matrix mul(matrix a, matrix b) {
        matrix ret;
        for (int i = 1; i <= sz; ++i) { // A i 行
            for (int j = 1; j <= sz; ++j) { // B j 列 
                for (int l = 1; l <= sz; ++l) { // foreach 
                    ret.a[i][j] = min(ret.a[i][j], a.a[i][l] + b.a[l][j]);
                }
            }
        }
        return ret;
    }
    matrix ksm(matrix a, int b) {
        matrix ans, bs = a;
        int fl = 0;
        while (b) {
            if (b & 1) {
                if (fl) ans = mul(ans, bs); else fl = 1, ans = bs;
            }
            bs = mul(bs, bs);
            b >>= 1;
        }
        return ans;
    }

    void clean() {
        ms(num, 0), sz = 0, ms(G, 0);
    }
    int solve() {
        scanf("%d%d%d%d", &k, &m, &s, &e);
        n = 1000;
        clean();
        for (int w, x, y, i = 1; i <= m; ++i) {
            scanf("%d%d%d", &w, &x, &y);
            x = gt(x), y = gt(y);
            whw.a[x][y] = G[x][y] = w;
            whw.a[y][x] = G[y][x] = w;
        }
        whw = ksm(whw, k);
        printf("%d\n", whw.a[gt(s)][gt(e)]);
        return 0; 
    }
}
int main() {
    flyinthesky::solve();
    return 0;
}