BZOJ

BZOJ 3571

传送门

首先，$N$个配$N$个还带权值，那就显然是用KM了

这道题的思路和最小乘积生成树的思路差不多，可惜BZOJ上一道裸最小乘积生成树被权限了…

考虑把每一种方案抽象成一个点$(\sum A,\sum B)$，然后把所有点扔到一个坐标系里面，然后答案就是使得一个反比例函数$x\cdot y=k$的$k$最小的点

拿KM找出$\sum A$最小的点和$\sum B$最小的点（带权最小匹配就把边权取反算带权最大匹配），最优解一定优于这两个点（或位于两个点之一），也就是在这两个点连成的线段下方

设最优解为$c$，$\sum A$最小的点为$a$,$\sum B$最小的点为$b$

易证明$c$一定在下凸壳上

那么题目就成了求解$a$到$b$的下凸壳。那么问题来了，总点数太多，不可能得到所有的点，怎么求下凸壳？

有一个显然的结论，$\vec{ab}$下方距离$\vec{ab}$最远的点一定在下凸壳上，即使得$S\Delta abc$最大的$c$一定在下凸壳上。

$$S\Delta abc=\frac{\vec{ab}\times\vec{ac}}{2}$$

展开再合起来能得到 （向量箭头不一样高差评）

$$S\Delta abc=\frac{\vec{ab}\times\vec{c}+\vec{a}\times\vec{b}}{2}$$

后面那坨已知，扔掉，即求$\vec{ab}\times\vec{c}$的最大值

找到$c$之后，更新答案，$c$不一定是最优解，我们要将下凸壳遍历，因此继续左右分治下去。

分析时间复杂度，由于要将下凸壳全遍历，设下凸壳上的点数为$S$，每次${\rm KM}$的复杂度是$\Theta(N^3)$，总复杂度是$\Theta(SN^3)$。

在实践中，这个算法是很优秀的，但是还是存在很坏的情况……（跑的还没有暴力快）

上面一大坨实际上是ZYK写的


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
program bzoj_3571;
const inf=maxlongint>>1;
type point=record x,y:longint; end;
var lx,ly,slack,pair:array[1..70]of longint;
    vx,vy:array[1..70]of boolean;
    w,a,b:array[1..70,1..70]of longint;
    t,n,i,j:longint;
    pa,pb:point;

function min(a,b:longint):longint;
begin
  if a<b then exit(a) else exit(b);
end;

function find(x:longint):boolean;
var t,y:longint;
begin
  vx[x]:=true;
  for y:=1 to n do
  begin
    if vy[y] then continue;
    t:=lx[x]+ly[y]-w[x,y];
    if t=0 then
    begin
      vy[y]:=true;
      if (pair[y]=0)or(find(pair[y])) then
      begin
        pair[y]:=x;
        exit(true);
      end;
    end else if t<slack[y] then slack[y]:=t;
  end;
  exit(false);
end;

function KM:point;
var i,j,d:longint;
begin
  fillchar(pair,sizeof(pair),0);
  for i:=1 to n do lx[i]:=-inf;
  fillchar(ly,sizeof(ly),0);
  for i:=1 to n do
  for j:=1 to n do if w[i,j]>lx[i] then lx[i]:=w[i,j];
  for i:=1 to n do
  begin
    fillchar(slack,sizeof(slack),$3f);
    repeat
      fillchar(vx,sizeof(vx),false);
      fillchar(vy,sizeof(vy),false);
      if find(i) then break;
      d:=inf;
      for j:=1 to n do if (not vy[j])and(slack[j]<d) then d:=slack[j];
      for j:=1 to n do if vx[j] then dec(lx[j],d);
      for j:=1 to n do if vy[j] then inc(ly[j],d) else dec(slack[j],d);
    until false;
  end;
  KM.x:=0;
  KM.y:=0;
  for i:=1 to n do inc(KM.x,a[pair[i],i]);
  for i:=1 to n do inc(KM.y,b[pair[i],i]);
end;

function solve(l,r:point):longint;
var t:point;
    i,j:longint;
begin
  for i:=1 to n do
  for j:=1 to n do w[i,j]:=a[i,j]*(r.y-l.y)-b[i,j]*(r.x-l.x);
  t:=KM;
  if ((t.x=l.x)and(t.y=l.y))or((t.x=r.x)and(t.y=r.y)) then exit(min(l.x*l.y,r.x*r.y));
  exit(min(solve(l,t),solve(t,r)));
end;

begin
  readln(t);
  for t:=1 to t do
  begin
    readln(n);
    for i:=1 to n do
    for j:=1 to n do read(a[i,j]);
    for i:=1 to n do
    for j:=1 to n do read(b[i,j]);
    for i:=1 to n do
    for j:=1 to n do w[i,j]:=-a[i,j];
    pa:=KM;
    for i:=1 to n do
    for j:=1 to n do w[i,j]:=-b[i,j];
    pb:=KM;
    writeln(solve(pa,pb));
  end;
end.