bzoj3571

同样的最小乘积XXX模型，这里显然是二分图带权匹配

我不会写KM……于是写了个费用流，由于是稠密图，会退化到n^4

然后本地跑了56s，交上去过了………………一定是我电脑太慢……

改天写个KM吧

  1 const inf=14000*14000;
  2 type node=record
  3        po,next,flow,cost:longint;
  4      end;
  5      point=record
  6        x,y:longint;
  7      end;
  8
  9 var e:array[0..2000010] of node;
 10     ma,pre,cur,p,d:array[0..145] of longint;
 11     v:array[0..145] of boolean;
 12     a,b,c:array[0..75,0..75] of longint;
 13     q:array[0..2000010] of longint;
 14     j,len,i,tt,t,n,ans:longint;
 15     p1,p2:point;
 16
 17 procedure add(x,y,f,c:longint);
 18   begin
 19     inc(len);
 20     e[len].po:=y;
 21     e[len].flow:=f;
 22     e[len].cost:=c;
 23     e[len].next:=p[x];
 24     p[x]:=len;
 25   end;
 26
 27 procedure build(x,y,f,c:longint);
 28   begin
 29     add(x,y,f,c);
 30     add(y,x,0,-c);
 31   end;
 32
 33 function spfa:boolean;
 34   var f,r,x,y,i,j:longint;
 35   begin
 36     fillchar(v,sizeof(v),false);
 37     for i:=1 to t do
 38       d[i]:=inf;
 39     d[0]:=0;
 40     f:=1;
 41     r:=1;
 42     q[1]:=0;
 43     while f<=r do
 44     begin
 45       x:=q[f];
 46       v[x]:=false;
 47       i:=p[x];
 48       while i<>-1 do
 49       begin
 50         y:=e[i].po;
 51         if e[i].flow>0 then
 52           if d[x]+e[i].cost<d[y] then
 53           begin
 54             d[y]:=d[x]+e[i].cost;
 55             cur[y]:=i;
 56             pre[y]:=x;
 57             if not v[y] then
 58             begin
 59               inc(r);
 60               q[r]:=y;
 61               v[y]:=true;
 62             end;
 63           end;
 64         i:=e[i].next;
 65       end;
 66       inc(f);
 67     end;
 68     if d[t]=inf then exit(false) else exit(true);
 69   end;
 70
 71 procedure change(j:longint);
 72   begin
 73     dec(e[j].flow);
 74     inc(e[j xor 1].flow);
 75   end;
 76
 77 function get:point;
 78   var i,j,x,y:longint;
 79   begin
 80     len:=-1;
 81     fillchar(p,sizeof(p),255);
 82     for i:=1 to n do
 83     begin
 84       build(0,i,1,0);
 85       build(i+n,t,1,0);
 86     end;
 87     x:=1; y:=1;
 88     for i:=1 to n do
 89       for j:=1 to n do
 90       begin
 91         build(i,j+n,1,c[i,j]);
 92         if c[i,j]<c[x,y] then
 93         begin
 94           x:=i;
 95           y:=j;
 96         end;
 97       end;
 98     ma[y+n]:=x;
 99     i:=p[x];
100     while i<>-1 do
101     begin
102       if (e[i].po=y+n) then change(i);
103       if e[i].po=0 then change(i xor 1);
104       i:=e[i].next;
105     end;
106     i:=p[y+n];
107     while i<>-1 do
108     begin
109       if e[i].po=t then
110       begin
111         change(i);
112         break;
113       end;
114       i:=e[i].next;
115     end;
116     while spfa do
117     begin
118       i:=t;
119       while i<>0 do
120       begin
121         j:=cur[i];
122         change(j);
123         if (i>n) and (i<>t) then
124           ma[i]:=pre[i];
125         i:=pre[i];
126       end;
127     end;
128     get.x:=0; get.y:=0;
129     for i:=1 to n do
130     begin
131       j:=ma[i+n];
132      // write(j,‘ ‘);
133       inc(get.x,a[j,i]);
134       inc(get.y,b[j,i]);
135     end;
136    // writeln;
137     if get.x*get.y<ans then ans:=get.x*get.y;
138   end;
139
140 function cross(a,b,c:point):longint;
141   begin
142     exit((a.x-c.x)*(b.y-c.y)-(a.y-c.y)*(b.x-c.x));
143   end;
144
145 procedure work(p1,p2:point);
146   var p:point;
147       i,j:longint;
148   begin
149     for i:=1 to n do
150       for j:=1 to n do
151         c[i,j]:=a[i,j]*(p1.y-p2.y)+b[i,j]*(p2.x-p1.x);
152     p:=get;
153     if cross(p2,p,p1)>=0 then exit;
154     work(p1,p);
155     work(p,p2);
156   end;
157
158 begin
159   readln(tt);
160   while tt>0 do
161   begin
162     dec(tt);
163     readln(n);
164     t:=2*n+1;
165     for i:=1 to n do
166       for j:=1 to n do
167         read(a[i,j]);
168     for i:=1 to n do
169       for j:=1 to n do
170         read(b[i,j]);
171     for i:=1 to n do
172       for j:=1 to n do
173         c[i,j]:=a[i,j];
174     ans:=14000*140000;
175     p1:=get;
176     for i:=1 to n do
177       for j:=1 to n do
178         c[i,j]:=b[i,j];
179     p2:=get;
180     work(p1,p2);
181     writeln(ans);
182   end;
183 end.