bzoj3168[Heoi2013]钙铁锌硒维生素(矩阵求逆,匈牙利算法,二分图字典序最小完备匹配)

Description
银河队选手名单出来了!小林,作为特聘的营养师,将负责银河队选手参加宇宙比赛的饮食。众所周知,前往宇宙
的某个星球,通常要花费好长好长的时间,人体情况在这之间会发生变化,因此,需要根据每天的情况搭配伙食,
来保证营养。小林把人体需要的营养分成了n种,这些营养包括但不限于铁,钙。他准备了2套厨师机器人,一套厨
师机器人有n个,每个厨师机器人只会做一道菜,这道菜一斤能提供第i种营养xi微克。想要吃这道菜的时候,只要
输入一个数,就能吃到对应数量的这道菜了。为防止摄入过量对身体造成的伤害,每个机器人还有防过量摄入药,
只要输入一个数,就能生成一定剂量的药,吃了这些药,就能减少相当于食用对应数目的这道菜提供的营养。小林
之所以准备2套厨师机器人,正是因为旅途漫漫,难以预计,也许某一个厨师机器人在途中坏掉,要是影响了银河
队选手的身体,就不好了。因此,第2套厨师机器人被用来做第1套的备用。小林需要为每一个第1套厨师机器人选
一个第2套厨师机器人作备份,使得当这个机器人坏掉时,用备份顶替,整套厨师机器人仍然能搭配出任何营养需
求,而且,每个第2套厨师机器人只能当一个第1套厨师机器人的备份。
Input
第一行包含一个正整数n。
接下来n行,每行n个整数,表示第1套厨师机器人做的菜每一斤提供的每种营养。
再接下来n行,每行n个整数,表示第2套厨师机器人做的菜每一斤提供的每种营养。
1≤n≤300,所有出现的整数均非负,且不超过10,000。
Output
第一行是一个字符串,如果无法完成任务,输出“NIE”,否则输出“TAK”
并跟着n行,第i行表示第i个第1套机器人的备份是哪一个第2套机器人。
为了避免麻烦,如果有多种可能的答案,请给出字典序最小的那一组。
Sample Input
3
1 0 0
0 1 0
0 0 1
2 3 0
0 7 8
0 0 9
Sample Output
TAK
1
2
3
HINT
Source


比较容易想到我们需要求出每个 a i a_i ai能匹配到的 b j b_j bj

求完之后建个二分图跑字典序最小完备匹配即可

怎么求呢?
首先 { a i } \{a_i\} {ai}是线性无关的
那么 b i b_i bi可以表示为 ∑ i = 1 n c i ∗ a i \sum_{i=1}^nc_i*a_i i=1nciai
c i c_i ci为常数

容易得到,如果 c j ≠ 0 c_j \neq 0 cj̸=0,说明 b i b_i bi不能由其他的 a k a_k ak表示出,那么 b i b_i bi和其他的 a k a_k ak线性无关,也就是 b i b_i bi可以替换 a j a_j aj

故我们只要求出每个 b i b_i bi通过 { a i } \{a_i\} {ai}的线性表示即可

怎么求呢?
假设 x = [ c 1 , c 2 , c 3 . . . . . c n ] x=[c_1,c_2,c_3.....c_n] x=[c1,c2,c3.....cn],这里 x x x是个 n n n阶向量
那么 x ∗ A = b x*A=b xA=b
x = b ∗ A − x=b*A^- x=bA
写个矩阵求逆即可

用浮点数的话需要注意精度
(我随便设的

#include<bits/stdc++.h>
using namespace std;
#define rep(i,j,k) for(int i = j;i <= k;++i)
#define repp(i,j,k) for(int i = j;i >= k;--i)
#define rept(i,x) for(int i = linkk[x];i;i = e[i].n)
#define P pair<int,int>
#define Pil pair<int,ll>
#define Pli pair<ll,int>
#define Pll pair<ll,ll>
#define pb push_back 
#define pc putchar
#define mp make_pair
#define file(k) memset(k,0,sizeof(k))
#define ll long long
#define M matrix
int rd()
{
    int sum = 0;char c = getchar();bool flag = true;
    while(c < '0' || c > '9') {if(c == '-') flag = false;c = getchar();}
    while(c >= '0' && c <= '9') sum = sum * 10 + c - 48,c = getchar();
    if(flag) return sum;
    else return -sum;
}
const double eps = 1e-10;
int n;
double a[310][310],b[310][310];
double f[310][310],g[310][310];
double c[310];
bool d[310][310];
int linkk[610],t,match[610];
bool vis[610],inq[610];
struct node{int n,y;}e[90100];
void mul(int i)
{
	rep(j,1,n)
	{
		c[j] = 0;
		rep(k,1,n) c[j] += b[i][k] * g[k][j];
	}
} 
void get_inv()
{
	rep(i,1,n) g[i][i] = 1.0;
	rep(i,1,n) 
	{
		rep(j,i,n) if(fabs(f[j][i])>eps)
		{
			if(j != i) rep(k,1,n) swap(f[i][k],f[j][k]),swap(g[i][k],g[j][k]);
			break;
		}
		if(fabs(f[i][i]) < eps) {printf("NIE\n");exit(0);}
		double t = f[i][i];
		rep(j,1,n) f[i][j] /= t,g[i][j] /= t;
		rep(j,1,n) if(i != j && fabs(f[j][i]) > eps)
		{
			t = f[j][i];
			rep(k,1,n) f[j][k] -= f[i][k]*t,g[j][k] -= g[i][k]*t;
		}
	}
}
void init()
{
	n = rd();
	rep(i,1,n) rep(j,1,n) a[i][j] = (double) rd();
	rep(i,1,n) rep(j,1,n) b[i][j] = (double) rd();
	rep(i,1,n) rep(j,1,n) f[i][j] = (double) a[i][j];
	get_inv();
}
bool find(int x)
{
	rep(i,1,n) if(d[x][i])
	    if(!vis[i+n])
	    {
	    	int y = i+n;
	    	vis[y] = true;
	    	if(!match[y] || find(match[y]))
	    	{
	    		match[y] = x;
	    		match[x] = y;
	    		return true;
			}
		}
	return false;
}
bool get(int st,int x)
{
	if(x == st) return true;
	rep(y,n+1,2*n) if(match[x] != y && d[x][y-n] && !inq[y] && !vis[y])
	{
		vis[y] = true;
		if(!inq[match[y]] && get(st,match[y]))
		{
			match[x] = y;
			match[y] = x;
			return true;
		}
 	}
	return false;
}
int main() 
{
	init();
	rep(i,1,n)
	{
		mul(i);rep(j,1,n) b[i][j] = c[j];
		rep(j,1,n) if(fabs(b[i][j]) > eps) d[j][i] = true;
	}
	int tot = 0;
	rep(i,1,n)  {file(vis);if(find(i)) tot++;}
	if(tot < n) {printf("NIE\n");return 0;}
	printf("TAK\n");
	rep(i,1,n)
	{
		int tmp = match[i];
		bool fl = false;
		rep(j,n+1,tmp-1) if(d[i][j-n] && !inq[j])
		{
			memset(vis,0,sizeof(vis));vis[j] = true;
			if(!inq[match[j]] && get(i,match[j])) 
			{
			    fl = true;inq[j+n] = inq[i] = false;
			    match[i] = j;match[j] = i;
			    inq[i] = inq[j] = true;
			    break;
			}
		}
		if(!fl) inq[i] = inq[tmp] = true;
	}
	rep(i,1,n) printf("%d\n",match[i]-n);
	return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值