BZOJ 3168: [Heoi2013]钙铁锌硒维生素 [线性基 Hungary 矩阵求逆]

3168: [Heoi2013]钙铁锌硒维生素

题意：给一个线性无关组A，再给一个B，要为A中每个向量在B中选一个可以代替的向量，替换后仍然线性无关。判断可行和求字典序最小的解

---

PoPoQQQ orz

显然是一个二分图匹配的模型

A是一个线性基，用它把B中每个向量表示出来，那么\(B_i\)可以替换\(A_j\)当且仅当表示\(B_i\)用到了\(A_j\)

可是A并不是每一位独立，怎么求表示啊？

A和B可以看成两个矩阵(横向量组)
\(C*A=B \rightarrow C=B*A^{-1}\)
\(C_{i,j}=1\)说明表示\(B_i\)用到了\(A_j\)，那么\(C^T\)就是这个二分图的邻接矩阵了

求矩阵的逆

这里说一种方法，对A进行高斯约当消元，右面的常数列换成单位矩阵。校园后，左面变成了单位矩阵，右面就是\(A^{-1}\)

二分图匹配字典序最小的解

求任意一个完美匹配，然后从1到n贪心选择字典序最小的解，方法和hungary类似，但是要比较匹配点和当前点的字典序

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
using namespace std;
typedef unsigned long long ll;
const int N=305, P=1e9+7;
inline int read() {
	char c=getchar(); int x=0, f=1;
	while(c<'0' || c>'9') {if(c=='-')f=-1; c=getchar();}
	while(c>='0' && c<='9') {x=x*10+c-'0'; c=getchar();}
	return x*f;
}

inline ll Pow(ll a, int b) {
	ll ans=1;
	for(; b; b>>=1, a=a*a%P)
		if(b&1) ans=ans*a%P;
	return ans;
}
inline void mod(int &x) {if(x<0) x+=P; else if(x>=P) x-=P;}

int n, g[N][N]; char s[N];
struct Matrix {
	int a[N][N];
	Matrix(){memset(a, 0, sizeof(a));}
	int* operator [](int x) {return a[x];}
	inline void im() {for(int i=1; i<=n; i++) a[i][i]=1;}
	void print() {for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) printf("%d%c",a[i][j],j==n?'\n':' ');}
}a, b, c;

Matrix inverse(Matrix a) {
	Matrix c; c.im();
	for(int i=1; i<=n; i++) {
		int r;
		for(r=i; r<=n; r++) if(a[r][i]) break;
		// r != n+1
		if(r!=i) for(int j=1; j<=n; j++) 
			swap(a[i][j], a[r][j]), swap(c[i][j], c[r][j]);
		ll inv = Pow(a[i][i], P-2); 
		for(int j=1; j<=n; j++) 
			a[i][j] = a[i][j]*inv%P, c[i][j] = c[i][j]*inv%P;
		for(int k=1; k<=n; k++) if(k!=i) {
			ll t = a[k][i]%P;
			for(int j=1; j<=n; j++) 
				mod(a[k][j] -= a[i][j]*t%P), mod(c[k][j] -= c[i][j]*t%P);
		}
	}
	return c;
}
Matrix operator *(Matrix a, Matrix b) {
	Matrix c;
	for(int i=1; i<=n; i++)
		for(int k=1; k<=n; k++) if(a[i][k])
			for(int j=1; j<=n; j++) 
				mod(c[i][j] += (ll)a[i][k]*b[k][j]%P);
	return c;
}

int vis[N], le[N];
bool dfs(int u) {
	for(int v=1; v<=n; v++) 
		if(!vis[v] && g[u][v]) {
			vis[v]=1;
			if(!le[v] || dfs(le[v])) {
				le[v]=u;
				return true;
			}
		}
	return false;
}
bool dfs(int u, int now) {
	for(int v=1; v<=n; v++) 
		if(!vis[v] && g[u][v]) {
			vis[v]=1;
			if(le[v]==now || (le[v]>now && dfs(le[v], now))) {
				le[v]=u;
				return true;
			}
		}
	return false;
}
int main() {
	freopen("in","r",stdin);
	//freopen("ferrous.in","r",stdin);
	//freopen("ferrous.out","w",stdout);
	n=read();
	for(int i=1; i<=n; i++) 
		for(int j=1; j<=n; j++) a[i][j] = read();
	for(int i=1; i<=n; i++) 
		for(int j=1; j<=n; j++) b[i][j] = read();
	c = b * inverse(a); //puts("c");c.print();
	//Matrix t = a * inverse(a); puts("t"); t.print();
	for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) if(c[i][j]) g[j][i]=1;
	for(int i=1; i<=n; i++) {
		memset(vis, 0, sizeof(vis));
		if(!dfs(i)) {puts("NIE"); return 0;}
	}
	puts("TAK");
	for(int i=1; i<=n; i++) {
		memset(vis, 0, sizeof(vis));
		dfs(i, i);
	}
	for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) if(le[j]==i) printf("%d\n",j);
}