BZOJ 3168 Luogu P4100 [HEOI2013]钙铁锌硒维生素 (矩阵求逆、二分图匹配)

线性代数+图论好题。

题目链接: (bzoj) https://www.lydsy.com/JudgeOnline/problem.php?id=3168

(luogu) https://www.luogu.org/problemnew/show/P4100

题解: 首先有一个结论是,设矩阵\(C\)满足\(CA=B\), 则\(A\)的第\(i\)行可以被\(B\)的第\(j\)行来替代当且仅当\(C_{j,i}\ne 0\).

这个的话就是一定要记住,矩阵初等变换的实质就是给它左乘一个初等矩阵,然后就很好理解了。

那么\(CA=B\)可以推出\(CAA^{-1}=BA^{-1}, C=BA^{-1}\)

数学被各种吊打啊……

然后我们就在\(O(n^3)\)时间内求出了对于每一个\(i,j\), \(A\)中第\(i\)行是否能被\(B\)中第\(j\)行替换

问题转化成了,给一张二分图,保证有完美匹配,求一个完美匹配使得\(A\)中的每个点在\(B\)中的匹配点构成的排列的字典序最小。

这个东西,貌似必须用匈牙利算法。。先跑一边随便求出一个完美匹配,然后从\(1\)号到\(n\)号每个点再匹配尽量小的,如果能找到使当前\(i\)号点更小,且不影响\(i\)号点前面的点的交错路,那么就可以更新答案了。

时间复杂度\(O(n^3)\).

代码
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cassert>
#include<iostream>
#define llong long long
using namespace std;

inline int read()
{
    int x=0; bool f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=0;
    for(; isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^'0');
    if(f) return x;
    return -x;
}

const int N = 300;
const int P = 942030731;
llong quickpow(llong x,llong y)
{
    llong cur = x,ret = 1ll;
    for(int i=0; y; i++)
    {
        if(y&(1ll<<i)) {y-=(1ll<<i); ret = ret*cur%P;}
        cur = cur*cur%P;
    }
    return ret;
}
llong mulinv(llong x) {return quickpow(x,P-2);}
struct Matrix
{
    llong a[N+3][N+3]; int n;
    Matrix() {}
    Matrix(int _n) {n = _n; for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) a[i][j] = 0ll;}
    void read(int _n)
    {
        n = _n;
        for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) scanf("%lld",&a[i][j]);
    }
    void write()
    {
        printf("%d\n",n);
        for(int i=1; i<=n; i++) {for(int j=1; j<=n; j++) printf("%lld ",a[i][j]); puts("");}
    }
    Matrix operator *(const Matrix &arg)
    {
        Matrix ret = Matrix(n);
        for(int i=1; i<=n; i++)
        {
            for(int k=1; k<=n; k++)
            {
                for(int j=1; j<=n; j++)
                {
                    ret.a[i][j] = (ret.a[i][j]+a[i][k]*arg.a[k][j])%P;
                }
            }
        }
        return ret;
    }
    Matrix inv()
    {
        Matrix ret = Matrix(n); for(int i=1; i<=n; i++) ret.a[i][i] = 1ll;
        for(int i=1; i<=n; i++)
        {
            if(a[i][i]==0)
            {
                bool found = false;
                for(int j=i+1; j<=n; j++)
                {
                    if(a[j][i])
                    {
                        for(int k=1; k<=n; k++) {swap(a[i][k],a[j][k]),swap(ret.a[i][k],ret.a[j][k]);}
                        found = true; break;
                    }
                }
                if(found==false) {ret.a[0][0] = P; return ret;}
            }
            for(int j=i+1; j<=n; j++)
            {
                llong coe = (P-a[j][i]*mulinv(a[i][i])%P)%P;
                for(int k=1; k<=n; k++)
                {
                    a[j][k] = (a[j][k]+coe*a[i][k])%P;
                    ret.a[j][k] = (ret.a[j][k]+coe*ret.a[i][k])%P;
                }
            }
        }
        for(int i=1; i<=n; i++)
        {
            llong coe = mulinv(a[i][i]);
            for(int j=1; j<=n; j++) {a[i][j] = a[i][j]*coe%P; ret.a[i][j] = ret.a[i][j]*coe%P;}
        }
//      write();
//      ret.write();
        for(int i=n; i>=1; i--)
        {
            for(int j=n; j>=i+1; j--)
            {
                llong coe = (P-a[i][j]*mulinv(a[j][j])%P)%P;
                a[i][j] = 0ll;
                for(int k=1; k<=n; k++) ret.a[i][k] = (ret.a[i][k]+coe*ret.a[j][k])%P;
            }
        }
        return ret;
    }
} a,b,aux,c;
int g[N+3][N+3];
int vis[N+3];
int match1[N+3],match2[N+3];
int n;

bool dfs1(int u)
{
    for(int i=1; i<=n; i++)
    {
        if(g[i][u]==true && vis[i]==false)
        {
            vis[i] = true;
            if(match2[i]==0 || dfs1(match2[i])==true)
            {
                match2[i] = u; match1[u] = i;
                return true;
            }
        }
    }
    return false;
}

bool dfs2(int u,int u0)
{
    for(int i=1; i<=n; i++)
    {
        if(g[i][u]==true && vis[i]==false)
        {
            vis[i] = true;
            if(match2[i]==u0 || (match2[i]>u0 && dfs2(match2[i],u0)==true))
            {
                match2[i] = u; match1[u] = i;
                return true;
            }
        }
    }
    return false;
}

int main()
{
    scanf("%d",&n);
    a.read(n); b.read(n);
    aux = a.inv();
    if(aux.a[0][0]==P) {printf("NIE"); return 0;}
    c = b*aux;
    for(int i=1; i<=n; i++)
    {
        for(int j=1; j<=n; j++) g[i][j] = c.a[i][j]==0 ? 0 : 1;
    }
//  for(int i=1; i<=n; i++) {for(int j=1; j<=n; j++) printf("%d",g[i][j]); puts("");}
    int mf = 0;
    for(int i=1; i<=n; i++)
    {
        for(int j=1; j<=n; j++) vis[j] = false;
        mf += dfs1(i);
    }
    if(mf<n) {printf("NIE"); return 0;}
    printf("TAK\n");
    for(int i=1; i<=n; i++)
    {
        for(int j=1; j<=n; j++) vis[j] = false;
        dfs2(i,i);
    }
    for(int i=1; i<=n; i++) printf("%d\n",match1[i]);
    return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值