BZOJ4386[POI2015] Wycieczki

BZOJ4386[POI2015] Wycieczki

Description

给定一张n个点m条边的带权有向图，每条边的边权只可能是1，2，3中的一种。

将所有可能的路径按路径长度排序，请输出第k小的路径的长度，注意路径不一定是简单路径，即可以重复走同一个点。

Input

第一行包含三个整数n,m,k(1<=n<=40，1<=m<=1000，1<=k<=10^18)。

接下来m行，每行三个整数u,v,c(1<=u,v<=n，u不等于v，1<=c<=3)，表示从u出发有一条到v的单向边，边长为c。

可能有重边。

Output

包含一行一个正整数，即第k短的路径的长度，如果不存在，输出-1。

Sample Input

6 6 11

1 2 1

2 3 2

3 4 2

4 5 1

5 3 1

4 6 3

Sample Output

4

HINT

长度为1的路径有1->2，5->3，4->5。

长度为2的路径有2->3，3->4，4->5->3。

长度为3的路径有4->6，1->2->3，3->4->5，5->3->4。

长度为4的路径有5->3->4->5。

Solution：

这是一道矩阵乘法的题。

首先来看边权是一的情况：

在$n*n$矩阵中$A$中，$A_{ij}$表示从$i$到$j$的路径条数。设$A_{k}[i][j]$表示长度为$k$的从$i$到$j$的路径条数。显然根据矩阵乘法$A_{m+n}=A_m*A_n$。$A_1$即是原图。（其实对于图上来说这与倍增Floyd的想法差不多）

要求任意两点之间的路径条数，我们可以添加一个虚拟源点，将它向其他点连边，并给它自己连一条自环（用于矩阵乘法转移）。然后$A_k[0][i]$就能表示以$i$结尾的长度为$k$的路径条数。

然后来看这一道题：因为矩阵乘法只能适用于边权为1的情况，而它的边权又只是1或2或3，所以我们把一个点拆成三个点，再进行连边，如图：

然后我们在这个图上做矩阵乘法（倍增Floyd），可以求出长度为$2^k$的路径个数。因为我们拆了点，有些点是虚点，所以把每个点作为最后的倒数第二个点处理，它的贡献就是$A[0][i]*cnt[i]$，即以 i <script type="math/tex" id="MathJax-Element-119">i</script>结尾的路径个数乘以它的出度（即将原来的点按出边长度分开处理，可以避免虚点造成的干扰）。

然后就使用倍增的试探方式求出这个答案就可以了。

本题有个神坑：

矩阵乘法（统计答案）的时候有可能爆long long，乘之前应要先判断有没有超过K，否则说明已经超过了K，就标记一下，不更新答案，直接return就可以了。

#include<stdio.h>
#include<string.h>
#include<iostream>
#define ll long long
#define M 125 
#define Exit {res.m[0][0]=-1;return res;}//标记返回 
using namespace std;
void Rd(int &res){
    char c;res=0;
    while(c=getchar(),!isdigit(c));
    do{
        res=(res<<1)+(res<<3)+(c^48);
    }while(c=getchar(),isdigit(c));
}
int n,m,tot,id[M][3],cnt[M];
ll K;
struct Matrix{
    ll m[M][M];
    Matrix(){
        memset(m,0,sizeof(m));
    }
    Matrix operator *(const Matrix &A)const{
        Matrix res;
        for(int i=0;i<=tot;i++)
            for(int j=0;j<=tot;j++)
                if(m[i][j]<0||A.m[i][j]<0)Exit;
        for(int i=0;i<=tot;i++)
            for(int j=0;j<=tot;j++)
                if(m[i][j]){
                    for(int k=0;k<=tot;k++)
                        if(A.m[j][k]){
                            if(m[i][j]>K/A.m[j][k])Exit;
                            res.m[i][k]+=m[i][j]*A.m[j][k];
                            if(res.m[i][j]>K)Exit;
                        }
                }
        return res;
    }
}a[65],b,c;
bool check(){
    if(b.m[0][0]<0)return false;
    ll k=0;
    for(int i=1;i<=tot;i++){
        if(b.m[0][i]&&cnt[i]){
            if(b.m[0][i]<0)return false;
            if(b.m[0][i]>K/cnt[i])return false;
            k+=b.m[0][i]*cnt[i];
            if(k>=K)return false;
        }
    }
    return k<K;
}
int main(){
    Rd(n);Rd(m);cin>>K;
    for(int i=1;i<=n;i++)
        for(int j=0;j<3;j++)
            id[i][j]=++tot;
    for(int i=1;i<=n;i++){
        for(int j=0;j<2;j++)
            a[0].m[id[i][j]][id[i][j+1]]++;
        a[0].m[0][id[i][0]]++;
    }
    a[0].m[0][0]++;
    while(m--){
        int x,y,z;
        Rd(x);Rd(y);Rd(z);
        z--;
        cnt[id[x][z]]++;
        a[0].m[id[x][z]][id[y][0]]++;
    }
    ll Max=K*3,ans=0;
    int len=0;
    while((1LL<<len)<=Max)len++;
    for(int i=1;i<len;i++)
        a[i]=a[i-1]*a[i-1];
    c.m[0][0]=1;
    for(int i=len-1;i>=0;i--){
        b=c*a[i];
        if(check()){
            ans|=(1LL<<i);
            memcpy(c.m[0],b.m[0],sizeof(b.m[0]));
        }
    }
    ans++;
    if(ans>=Max)ans=-1;
    cout<<ans<<endl;
    return 0;
}