魔术球问题 2011-12-29

算法实现题8-4 魔术球问题（习题 8-14）
´问题描述：
假设有n根柱子，现要按下述规则在这n根柱子中依次放入编号为 1，2，3，¼的球。
（1）每次只能在某根柱子的最上面放球。
（2）在同一根柱子中，任何2个相邻球的编号之和为完全平方数。
试设计一个算法，计算出在n根柱子上最多能放多少个球。例如，在4 根柱子上最多可
放11个球。
´编程任务：
对于给定的n，计算在 n根柱子上最多能放多少个球。
´数据输入：
由文件input.txt提供输入数据。文件第1 行有 1个正整数n，表示柱子数。
´结果输出:
程序运行结束时，将 n 根柱子上最多能放的球数以及相应的放置方案输出到文件
output.txt中。文件的第一行是球数。接下来的n行，每行是一根柱子上的球的编号。

输入文件示例
input.txt
4

输出文件示例
output.txt
11
1 8 
2 7 9 
3 6 10 
4 5 11

————————————————————

将每个点 i拆成 xi和yi两个点

对于i,j 有i+j为完全平方数，且j>i 则建边 ｛xi,yj}

｛s,xi},{yi,t} 权值均为1。

依次枚举答案，通过网络最大流判断是否可行。

————————————————————

  1 Program Stone;
  2 var n,le,flow,m,s,t:longint;
  3     head,vh,dis,cur,pre,last:array[-1..5000]of longint;
  4     next,date,point:array[-50000..50000]of longint;
  5  procedure add(x,y:longint);
  6   begin
  7      inc(le);
  8      date[le]:=1;
  9      point[le]:=y;
 10      next[le]:=head[x];
 11      head[x]:=le;
 12      point[-le]:=x;
 13      next[-le]:=head[y];
 14      head[y]:=-le;
 15   end;
 16  function min(a,b:longint):longint;
 17   begin
 18     if a<b then min:=a else min:=b;
 19   end;
 20  function aug(x,nf:longint):longint;
 21  var i,j,l,d,minh,ins:longint;
 22   begin
 23     if x=t then exit(nf);
 24     l:=nf;
 25     i:=cur[x];
 26     while i<>0 do
 27      begin
 28        if (date[i]>0)and(dis[point[i]]+1=dis[x]) then
 29           begin
 30               cur[x]:=i;
 31               d:=aug(point[i],min(l,date[i]));
 32               if (d>0)and(x mod 2<>0) then begin
 33                                                       last[x]:=point[i];
 34                                                       pre[point[i]]:=x;
 35                                                    end;
 36               dec(date[i],d);
 37               inc(date[-i],d);
 38               dec(l,d);
 39               if (dis[s]=m*2+2)or(l=0) then exit(nf-l);
 40           end;
 41        i:=next[i];
 42      end;
 43     if l=nf then
 44      begin
 45         minh:=m*2+2;
 46         i:=head[x];
 47         while i<>0 do
 48          begin
 49            if (date[i]>0)and(dis[point[i]]<minh) then begin minh:=dis[point[i]];ins:=i;end;
 50            i:=next[i];
 51          end;
 52         cur[x]:=ins;
 53         dec(vh[dis[x]]);
 54         if vh[dis[x]]=0 then dis[s]:=m*2+2;
 55         dis[x]:=minh+1;
 56         inc(vh[dis[x]]);
 57      end;
 58     aug:=nf-l;
 59   end;
 60  procedure print;
 61  var i,j,k:longint;
 62   begin
 63     writeln(m-1);
 64     for i:=1 to m-1 do
 65      if pre[i*2]=0 then
 66       begin
 67         j:=i*2-1;
 68         while true do
 69          begin
 70            write((j+1)div 2,‘ ‘);
 71            if last[j]=0 then break;
 72            j:=last[j]-1;
 73          end;
 74         writeln;
 75       end;
 76     close(input);close(output);
 77     halt;
 78   end;
 79  procedure main;
 80  var i,j,k:longint;
 81   begin
 82     s:=0;t:=-1;m:=0;
 83     while true do
 84      begin
 85        inc(m);
 86        add(s,m*2-1);add(m*2,t);
 87        for i:=1 to m-1 do
 88         if sqr(trunc(sqrt(i+m)))=i+m then add(2*i-1,2*m);
 89        fillchar(vh,sizeof(vh),0);
 90        fillchar(dis,sizeof(dis),0);
 91        vh[0]:=m*2+2;
 92        for i:=-1 to m*2 do cur[i]:=head[i];
 93        while dis[s]<m*2+2 do inc(flow,aug(s,maxint));
 94        if m-flow>n then print;
 95      end;
 96
 97   end;
 98
 99 Begin
100  assign(input,‘prog84.in‘);assign(output,‘prog84.out‘);
101  reset(input);rewrite(output);
102    readln(n);
103    main;
104 end.
105
106