网络流24题——魔术球问题

发表于2022-06-14更新于2022-06-14阅读次数

题意：

有 n 个柱子，现往上依次放入编号为 1, 2, 3, 4, … 的球，要求：

  1. 每次只能在某跟柱子的最上面放球
  2. 同一根柱子中，相邻两个球的和为完全平方数

问你最多能放多少个球，并输出方案。

这道题我第一次做网络流 24 的时候（大概是 17 年）就看到了，当时应该是看了题解之后仍然一知半解，之后又尝试过几次，都没做出来。结果 5 年后随便看了一眼，居然一下子就做出来了，属实有点以外。是我思考能力真的增强了呢，还是只是恰好灵稽一动呢，我不好说 233 。

这道题要求最大的球数，乍一看好像要 run 最大流，最大的流量就表示能放的球数。但是从这个思路出发的话会发现很难设计约束。

如果转念一想，尝试能不能用二分等手段把问题转化为可行性问题，再用网络流判定，问题就简单多了。

假设现在已经确定球数，要判断 m 个球是否可行，那么只需抓住关键性质：完全平方数的限制只和相邻的球有关。如果我们能对每一个球都找到一个邻居，并且每个球最多有两个邻居，那就必然能组成若干条链。

继续往下思考，就可以发现一个柱子本质就是一个链表，每个球都指向自己下面的球。柱子本身可以抽象成一个节点，任何人都可以都往上连，表示新开一个节点。除了柱子就只能连向值比自己小，而且和为完全平方数的点。

这样左边是 m 个点，右边是 m 个点 + n 个柱子点，左边往右边可以放的点连边，权值为 1 。s 到左边，右边到 t 都连权值为 1 的边，表示一个点只能用一次。如果跑满流就表示存在合法解，这其实也是个完美匹配问题。

接下来可以二分答案，然后判断答案合法性，但是其实这道题点数还是很多的（答案最大有 1500），边数也很恐怖，所以感觉二分会很慢。我用的方法是每次在残余网络的基础上额外左右加两个点和若干边，然后再执行 dinic ，如果流量增量为 1 就表示此答案合法，可以继续往上加。dinic 算法的原理应该是保证即使是在跑了一半的残余网络上，也可以修正出最大流的。但是之前的流量已经被统计过了，所以 dinic 的结果实际上是流量增量。

输出方案可以检查中间的边，然后就可以建出若干个链表，表示每一个柱子。

据说这题有贪心做法，证明在 loj 讨论区上，有兴趣的大神可以观摩一下。

代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
using namespace std;
const int oo=1000000007;
struct edg{int y,nxt;  int v;}e[5100000];  int lk[5100],ltp=1;
void ist(int x,int y,int z){
	e[++ltp]=(edg){y,lk[x],z};  lk[x]=ltp;
	e[++ltp]=(edg){x,lk[y],0};  lk[y]=ltp;
}
int n,m;  int s,t;
int lvl[5100];
int q[5100],hd=0;
int crt[5100];
int aq[5100];
bool gtlvl(){
	memset(lvl,0,sizeof(lvl));
	lvl[q[hd=1]=s]=1;
	for(int k=1;k<=hd;++k){
		crt[q[k]]=lk[q[k]];
		for(int i=lk[q[k]];i;i=e[i].nxt)if(e[i].v && !lvl[e[i].y]){
			q[++hd]=e[i].y;
			lvl[e[i].y]=lvl[q[k]]+1;
		}
	}
	return lvl[t];
}
int mxflw(int x,int y){
	if(x==t)  return y;
	int bwl=0,flw=0;
	for(int i=crt[x];i && bwl<y;i=e[i].nxt)
		if(lvl[e[i].y]==lvl[x]+1 && e[i].v)
			if((flw=mxflw(e[i].y,min(y-bwl,e[i].v)))){
				e[i].v-=flw,e[i^1].v+=flw;
				bwl+=flw;
				if(!e[i].v)  crt[x]=e[i].nxt;
			}
	if(!bwl)  lvl[x]=0;
	return bwl;
}
int dnc(){
	int bwl=0,flw=0;
	while(gtlvl())while((flw=mxflw(s,oo)))  bwl+=flw;
	return bwl;
}
bool chc(int x,int y){
	int q=sqrt(x+y);
	for(int i=-2;i<=2;++i)
		if((q+i)*(q+i)==x+y)  return true;
	return false;
}
void svans(int x){
	memset(aq,0,sizeof(aq));
	for(int k=1;k<=x;++k){
		int y=n+k*2;
		for(int i=lk[y];i;i=e[i].nxt)
			if(e[i].y!=s && e[i].v==0)
				aq[e[i].y]=k;
	}
}
void otans(){
	for(int i=1;i<=n;++i){
		int x=aq[i];
		while(x!=0){
			printf("%d",x);
			x=aq[n+x*2+1];
			if(x!=0)  printf(" ");
		}
		printf("\n");
	}
}
int main(){
	scanf("%d",&n);
	s=5000,t=5001;
	for(int i=1;i<=n;++i)  ist(i,t,1);
	int ans=0;
	for(int i=1;;++i){
		ist(s,n+i*2,1);
		ist(n+i*2+1,t,1);
		for(int j=1;j<=n;++j)  ist(n+i*2,j,1);
		for(int j=1;j<i;++j)if(chc(i,j))
			ist(n+i*2,n+j*2+1,1);
		int tmp=dnc();
		if(tmp==0){
			ans=i-1;
			break;
		}
		else{
			svans(i);
		}
	}
	printf("%d\n",ans);
	otans();
	return 0;
}