矩阵快速幂加速最短路

矩阵快速幂加速最短路

通常用来优化Floyd的实现

[NOI Online #1 入门组] 魔法

题目描述

C 国由 \(n\) 座城市与 \(m\) 条有向道路组成，城市与道路都从 \(1\) 开始编号，经过 \(i\) 号道路需要 \(t_i\) 的费用。

现在你要从 \(1\) 号城市出发去 \(n\) 号城市，你可以施展最多 \(k\) 次魔法，使得通过下一条道路时，需要的费用变为原来的相反数，即费用从 \(t_i\) 变为 \(-t_i\)。请你算一算，你至少要花费多少费用才能完成这次旅程。

注意：使用魔法只是改变一次的花费，而不改变一条道路自身的 \(t_i\)；最终的费用可以为负，并且一个城市可以经过多次（包括 \(n\) 号城市）。

输入格式

输入的第一行有三个整数，分别代表城市数 \(n\)，道路数 \(m\) 和魔法次数限制 \(k\)。

第 \(2\) 到第 \((m + 1)\) 行，每行三个整数。第 \((i + 1)\) 行的整数 \(u_i, v_i, t_i\) 表示存在一条从 \(u_i\) 到 \(v_i\) 的单向道路，经过它需要花费 \(t_i\)。

输出格式

输出一行一个整数表示最小的花费。

样例 #1

样例输入 #1

4 3 2
1 2 5
2 3 4
3 4 1

样例输出 #1

-8

样例 #2

样例输入 #2

2 2 2
1 2 10
2 1 1

样例输出 #2

-19

提示

输入输出样例 1 解释

依次经过 \(1\) 号道路、\(2\) 号道路、\(3\) 号道路，并在经过 \(1,2\) 号道路前使用魔法。

输入输出样例 2 解释

依次经过 \(1\) 号道路、\(2\) 号道路、\(1\) 号道路，并在两次经过 \(1\) 号道路前都使用魔法。

数据规模与约定

本题共 \(20\) 个测试点，各测试点信息见下表。

测试点编号	\(n \leq\)	\(m \leq\)	\(k \leq\)	无环
\(1 \sim 2\)	\(5\)	\(20\)	\(0\)	不保证
\(3 \sim 4\)	\(10\)	\(20\)	\(50\)	不保证
\(5 \sim 6\)	\(10\)	\(20\)	\(0\)	不保证
\(7 \sim 8\)	\(20\)	\(200\)	\(50\)	是
\(9 \sim 10\)	\(20\)	\(200\)	\(0\)	不保证
\(11 \sim 12\)	\(100\)	\(200\)	\(50\)	是
\(13 \sim 14\)	\(100\)	\(200\)	\(50\)	不保证
\(15 \sim 18\)	\(100\)	\(2500\)	\(1000\)	不保证
\(19 \sim 20\)	\(100\)	\(2500\)	\(10^6\)	不保证

对于【无环】一栏为“是”的测试点，保证给出的图是一张有向无环图，否则不对图的形态做任何保证。

对于全部的测试点，保证：

• \(1 \leq n \leq 100\)，\(1 \leq m \leq 2500\)，\(0 \leq k \leq 10^6\)。
• \(1 \leq u_i, v_i \leq n\)，\(1 \leq t_i \leq 10^9\)。
• 给出的图无重边和自环，且至少存在一条能从 \(1\) 到达 \(n\) 的路径。

首先有70分的做法，即使用动态规划。

将边存入数组中，\(g[i].u\)和\(g[i].v\)表示边的两端，g[i].w表示边权

设\(dp[i][j][k]\)表示从i到j使用k次魔法时的最短路。

当k=0时:\(dp[i][j][0]=min(dp[i][j][0],dp[i][k][0]+dp[k][j][0])\)

当k=1时:\(dp[i][j][1]=min(dp[i][j][1],dp[i][g[k].u][0]+dp[g[k].v][j][0]-g[k].w)\)

使用k次魔法:\(dp[i][j][k]=min(dp[i][j][k],dp[i][l][k-1]+dp[l][j][1])\)

\(k-1+1\)刚好为k，表示魔法的使用次数。

最后输出\(dp[1][n][k]\)即为答案。

但是这个做法的时间复杂度为\(O(n^3k)\),空间复杂度为\(O(n^2k)\),无法AC。

这里就需要用一个工具——矩阵加速来降低复杂度了。

先对处理样例1时dp数组进行打表:

k=0时:

0 5 9 10
\ 0 4 5 
\ \ 0 1
\ \ \ 0

k=1时:

\ -5 -1  0
\  0 -4 -3
\  \  0 -1
\  \  \  0

k=2时:

0 -5 -9 -8
\  0 -4 -5
\  \  0 -1
\  \  \  0

按照题目要求，可以把矩阵乘法进行修改:

\(\begin{bmatrix}a_{1,1}&a_{1,2}&a_{1,3}&a_{1,4}\\a_{2,1}&a_{2,2}&a_{2,3}&a_{2,4}\\a_{3,1}&a_{3,2}&a_{3,3}&a_{3,4}\\a_{4,1}&a_{4,2}&a_{4,3}&a_{4,4}\\\end{bmatrix} *\begin{bmatrix}b_{1,1}&b_{1,2}&b_{1,3}&b_{1,4}\\b_{2,1}&b_{2,2}&b_{2,3}&b_{2,4}\\b_{3,1}&b_{3,2}&b_{3,3}&b_{3,4}\\b_{4,1}&b_{4,2}&b_{4,3}&b_{4,4}\\\end{bmatrix}=\begin{bmatrix}min(a_{1,k},b_{k,1})&min(a_{1,k},b_{k,2})&min(a_{1,k},b_{k,3})&min(a_{1,k},b_{k,4})\\min(a_{2,k},b_{k,1})&min(a_{2,k},b_{k,2})&min(a_{2,k},b_{k,3})&min(a_{2,k},b_{k,4})\\min(a_{3,k},b_{k,1})&min(a_{3,k},b_{k,2})&min(a_{3,k},b_{k,3})&min(a_{3,k},b_{k,4})\\min(a_{4,k},b_{k,1})&min(a_{4,k},b_{k,2})&min(a_{4,k},b_{k,3})&min(a_{4,k},b_{k,4})\\\end{bmatrix}\)

验证后可得到

另初始矩阵为\(E_{1}\)，第二个矩阵为\(E_{2}\),可得到\(E_{k}=E_{1}*E_{2}^k\)。

\(E_{1}\)和\(E_{2}\)可用动态规划(Floyd) 直接求得，然后在用ksm即可快速得到答案。

时间复杂度为\(O(n^2m+n^3log(k))\),空间复杂度为\(O(n^2)\),完全可通过此题。

#include<bits/stdc++.h>
#define inf 1000000000000000
using namespace std;
const int N=110;
struct Matrix{		//矩阵
    long long a[N][N];
    int n,m;
    Matrix(int n=0,int m=0):n(n),m(m){}
    void init(){
        for(int i=1;i<=n;i++){
            for(int j=1;j<=n;j++){
                if(i==j) a[i][j]=0;
                else a[i][j]=inf;
            }
        }
    }
    Matrix operator*(const Matrix&A) const{
        Matrix ret(n,A.n);
        for(int i=1;i<=n;i++){
            for(int j=1;j<=A.m;j++){
                long long tmp=inf;
                for(int k=1;k<=m;k++) tmp=min(tmp,a[i][k]+A.a[k][j]);
                ret.a[i][j]=tmp;
            }
        }
        return ret;
    }
};
Matrix ksm(Matrix A,int b){		//快速幂
    Matrix ret(A.n,A.n);
    ret.init();
    while(b){
        if(b&1) ret=ret*A;
        b>>=1;
        A=A*A;
    }
    return ret;
}
int main(){
    int n,m,k;
    scanf("%d%d%d",&n,&m,&k);
    Matrix D(n,n),K(n,n);
    D.init();
    K.init();
    for(int i=1;i<=m;i++){
        int u,v,w;
        scanf("%d%d%d",&u,&v,&w);
        D.a[u][v]=w;
        K.a[u][v]=-w;
    }
    Matrix G=ksm(D,n-1);
    G=G*ksm(K*G,k);
    printf("%lld",G.a[1][n]);
    return 0;
}