LightOJ 1356 Prime Independence 二分图最大独立集,HK算法

时间:2023-03-09 02:11:48
LightOJ 1356 Prime Independence 二分图最大独立集,HK算法

这个题唯一需要说的就是普通的匈牙利算法是O(nm)的,过不了

然后HK算法可以O(n^0.5m),这个算法可以每次找很多同样长度的最短增广路

分析见:http://www.hardbird.net/lightoj-1356-prime-independence%E6%9C%80%E5%A4%A7%E7%8B%AC%E7%AB%8B%E9%9B%86-hopcroft-carp%E7%AE%97%E6%B3%95/

#include <cstdio>
#include <iostream>
#include <ctime>
#include <vector>
#include <cmath>
#include <map>
#include <queue>
#include <algorithm>
#include <cstring>
using namespace std;
typedef long long LL;
const int N=4e4+;
const int INF=0x3f3f3f3f;
int prime[],cnt,T,n;
int a[N],p[][];
vector<int>fac[N],g[N];
void getprime(){
bool v[];
memset(v,,sizeof(v));
for(int i=;i*i<=;++i)
if(!v[i])
for(int j=i*i;j<=;j+=i)
v[j]=;
for(int i=;i<=;++i)
if(!v[i])prime[++cnt]=i;
}
int Mx[N],My[N],dx[N],dy[N],dis;
bool used[N];
bool SearchP(){
queue<int>q;
dis=INF;
memset(dx,-,sizeof(dx));
memset(dy,-,sizeof(dy));
for(int i=;i<=n;++i)
if(Mx[i]==-){
q.push(i);
dx[i]=;
}
while(!q.empty()){
int u=q.front();
q.pop();
if(dx[u]>dis)break;
for(int i=;i<g[u].size();++i){
int v=g[u][i];
if(dy[v]!=-)continue;
dy[v]=dx[u]+;
if(My[v]==-)dis=dy[v];
else{
dx[My[v]]=dy[v]+;
q.push(My[v]);
}
}
}
return dis!=INF;
}
bool dfs(int u){
for(int i=;i<g[u].size();++i){
int v=g[u][i];
if(!used[v]&&dy[v]==dx[u]+){
used[v]=true;
if(My[v]!=-&&dy[v]==dis)continue;
if(My[v]==-||dfs(My[v])){
My[v]=u;
Mx[u]=v;
return true;
}
}
}
return false;
}
int MaxMatch(){
int res=;
memset(Mx,-,sizeof(Mx));
memset(My,-,sizeof(My));
while(SearchP()){
memset(used,,sizeof(used));
for(int i=;i<=n;++i)
if(Mx[i]==-&&dfs(i))++res;
}
return res;
}
int main()
{
getprime();
int cas=;
scanf("%d",&T);
while(T--){
scanf("%d",&n);
memset(p,-,sizeof(p));
for(int i=;i<=n;++i){
scanf("%d",&a[i]);
g[i].clear();
fac[i].clear();
int tot=,t=a[i];
for(int j=;prime[j]*prime[j]<=t;++j){
if(t%prime[j]==){
fac[i].push_back(prime[j]);
while(t%prime[j]==)t/=prime[j],++tot;
}
}
if(t>)fac[i].push_back(t),++tot;
p[tot&][a[i]]=i;
}
for(int i=;i<=n;++i){
if(p[][a[i]]!=-){
for(int j=;j<fac[i].size();++j){
int tmp=a[i]/fac[i][j];
if(p[][tmp]==-)continue;
g[i].push_back(p[][tmp]);
}
}
else{
for(int j=;j<fac[i].size();++j){
int tmp=a[i]/fac[i][j];
if(p[][tmp]==-)continue;
g[p[][tmp]].push_back(i);
}
}
}
printf("Case %d: %d\n",++cas,n-MaxMatch());
}
return ;
}