【SOJ#824】【SR#12】传播者

Author Avatar
mrsrz 2月 24, 2021

网络流，LGV 引理，线段树维护矩阵，高斯消元

题目大意：

给定一个具有 $n$ 层的分层图，每一层有 $k$ 个点，第 $i$ 层和第 $i+1$ 层之间连了若干条有向边。

设 $S,T$ 为两个点集，满足 $S$ 中所有点在第 $l$ 层，$T$ 中所有点在第 $r$ 层，且 $l<r$。

记 $cut(S,T)$ 为 $S$ 和 $T$ 的最小点割 —— 即至少要去掉多少个点 (及其关联的边)，才能保证不存在一条从 $s$ 到 $t$ 的有向路径，其中 $s\in S$，$t\in T$，且 $s,t$ 均未被去掉。

题解：

本题在《2021 年信息学奥林匹克 中国国家集训队论文》的《浅谈线性代数与图论的关系》中出现。

考虑如何处理单组询问，只需要拆点后跑最大流即可。

考察这里的最大流有什么性质。可以发现，若最大流为 $t$，则它实际上找到了 $t$ 条起点在 $S$ 中，终点在 $T$ 中，且不经过重复点（包括起点和终点）的路径。

先考虑一个比较简单的问题，我们假设 $S$ 和 $T$ 都是全集，并且只需要判断答案是否为 $k$。

使用 Lindström–Gessel–Viennot lemma（LGV 引理）即可解决这个问题。

该引理可用来解决有向无环图上不相交路径计数问题，其内容如下。

定义 $\omega(P)$ 表示路径 $P$ 上的所有边权的乘积，$e(u,v)$ 表示 $u$ 到 $v$ 的所有路径 $P$ 的 $\omega(P)$ 之和，即

给定大小为 $n$ 的起点集合 $S_1,S_2,\ldots,S_n$ 以及终点集合 $T_1,T_2,\ldots,T_n$。

一组 $S$ 到 $T$ 的不相交路径 $P$ 满足对于 $i\in[1,n]$，都有一条 $S_i$ 到 $T_{\sigma(P)_i}$ 的路径 $P_i$，且对于 $i\neq j$，$S_i$ 和 $S_j$ 出发的路径没有公共点。

记 $N(p)$ 表示排列 $p$ 的逆序对数。

记

则有

也就是说，我们只需要知道 $S$ 中每个点到 $T$ 中每个点的路径条数，就可以使用高斯消元在 $O(k^3)$ 的时间复杂度里求出 $S$ 到 $T$ 的不相交匹配路径的带权和。

这里我们只需要知道是否存在这样的路径，所以我们可以选择一个质数 $P$ 作为模数进行行列式求值，若结果非 $0$，则一定存在这样的路径。

由于这里求的是带权和，所以可能出现行列式值为 $0$，实际上存在方案的情况。解决方案是对每条边赋一个随机权值，然后进行行列式求值。可以认为这样做的出错概率为 $1-\frac 1 P$，取 $10^9$ 左右的模数可以保证非常高的正确率。

至于 $e(S_i,T_i)$ 的求法，可以使用矩阵乘法来解决。修改和询问相当于对单个矩阵的修改，以及询问区间的矩阵乘积，使用线段树维护即可。

上述矩阵的乘积恰好就是 LGV 引理中的行列式 $M$。

现在考虑原问题，我们只保留与点集 $S,T$ 中的点有关的值，可以得到一个 $|S|\times|T|$ 的矩阵。

我们需要找一个 $t$ 满足：存在一个 $t\times t$ 的子矩阵，行列式不为 $0$，且任意一个 $(t+1)\times(t+1)$ 的子矩阵，行列式都为 $0$。

因此矩阵的秩就是答案，使用高斯消元可以 $O(k^3)$ 求矩阵的秩。

算法的总时间复杂度为 $O((n+q\log n)k^3)$，空间复杂度 $O(nk^2)$。

Code：

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
typedef unsigned long long u64;
const int md=1004535809,N=8300,W=24;
mt19937 mt(time(0)+(size_t)new int);
int n,k,q,A[W+2][W+2];
char mp[W+2],gl[W+2],gr[W+2];
inline int pow(int a,int b){
    int res=1;
    for(;b;b>>=1,a=(u64)a*a%md)if(b&1)res=(u64)res*a%md;
    return res;
}
struct Matrix{
    int A[W][W];
    inline Matrix(){memset(A,0,sizeof A);}
    inline Matrix operator*(const Matrix&b)const{
        Matrix c;
        for(int i=0;i<24;++i)
        for(int j=0;j<24;++j){
            u64 tmp1=(u64)A[i][ 0]*b.A[ 0][j]+(u64)A[i][ 1]*b.A[ 1][j]+(u64)A[i][ 2]*b.A[ 2][j]+(u64)A[i][ 3]*b.A[ 3][j]+
            (u64)A[i][ 4]*b.A[ 4][j]+(u64)A[i][ 5]*b.A[ 5][j]+(u64)A[i][ 6]*b.A[ 6][j]+(u64)A[i][ 7]*b.A[ 7][j];
            u64 tmp2=(u64)A[i][ 8]*b.A[ 8][j]+(u64)A[i][ 9]*b.A[ 9][j]+(u64)A[i][10]*b.A[10][j]+(u64)A[i][11]*b.A[11][j]+
            (u64)A[i][12]*b.A[12][j]+(u64)A[i][13]*b.A[13][j]+(u64)A[i][14]*b.A[14][j]+(u64)A[i][15]*b.A[15][j];
            u64 tmp3=(u64)A[i][16]*b.A[16][j]+(u64)A[i][17]*b.A[17][j]+(u64)A[i][18]*b.A[18][j]+(u64)A[i][19]*b.A[19][j]+
            (u64)A[i][20]*b.A[20][j]+(u64)A[i][21]*b.A[21][j]+(u64)A[i][22]*b.A[22][j]+(u64)A[i][23]*b.A[23][j];
            c.A[i][j]=(tmp1%md+tmp2%md+tmp3)%md;
        }
        return c;
    }
}qwq[N],D[N*2],K;
void build(int l,int r,int o){
    if(l==r)D[o]=qwq[l];else{
        const int mid=(l+r)/2;
        build(l,mid,o<<1),build(mid+1,r,o<<1|1);
        D[o]=D[o<<1]*D[o<<1|1];
    }
}
void modify(int l,int r,int o,const int&pos,const int&u,const int&v){
    if(l==r){
        if(int&g=D[o].A[u][v])g=0;
        else g=mt()%(md-1)+1;
        return;
    }
    const int mid=(l+r)/2;
    if(pos<=mid)modify(l,mid,o<<1,pos,u,v);
    else modify(mid+1,r,o<<1|1,pos,u,v);
    D[o]=D[o<<1]*D[o<<1|1];
}
void query(int l,int r,int o,const int&L,const int&R){
    if(L<=l&&r<=R)K=K*D[o];else{
        const int mid=(l+r)/2;
        if(L<=mid)query(l,mid,o<<1,L,R);
        if(mid<R)query(mid+1,r,o<<1|1,L,R);
    }
}
int gauss(int n,int m){
    int rk=0;
    for(int i=0;i<m&&rk<n;++i){
        for(int j=rk;j<n;++j)if(A[j][i]){
            if(rk!=j)swap(A[rk],A[j]);
            break;
        }
        if(!A[rk][i]){
            continue;
        }
        int inv=pow(A[rk][i],md-2);
        for(int j=rk+1;j<n;++j){
            int d=(u64)inv*A[j][i]%md;
            for(int k=i;k<m;++k)
            A[j][k]=(A[j][k]-(LL)A[rk][k]*d%md+md)%md;
        }
        ++rk;
    }
    return rk;
}
int solve(){
    static int a[W+2],b[W+2],ta,tb;
    ta=tb=0;
    for(int i=0;i<k;++i){
        if(gl[i]=='1')a[ta++]=i;
        if(gr[i]=='1')b[tb++]=i;
    }
    for(int i=0;i<ta;++i)
    for(int j=0;j<tb;++j)A[i][j]=K.A[a[i]][b[j]];
    return gauss(ta,tb);
}
int main(){
    ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    cin>>n>>k>>q;
    for(int i=1;i<n;++i){
        for(int s=0;s<k;++s){
            cin>>mp;
            for(int t=0;t<k;++t)
            if(mp[t]=='1')qwq[i].A[s][t]=mt()%(md-1)+1;
        }
    }
    build(1,n-1,1);
    while(q--){
        char op;
        cin>>op;
        if(op=='T'){
            int x,u,v;
            cin>>x>>u>>v,--u,--v;
            modify(1,n-1,1,x,u,v);
        }else{
            int l,r;
            cin>>l>>r>>gl>>gr;
            K=Matrix();
            for(int i=0;i<k;++i)K.A[i][i]=1;
            query(1,n-1,1,l,r-1);
            cout<<solve()<<'\n';
        }
    }
    return 0;
}
Powered by Hexo
Theme - Material