POJ-2976 Dropping tests(01分数规划入门)

01分数规划问题：

所谓的01分数规划问题就是指这样的一类问题，给定两个数组，a[i]表示选取i的收益，b[i]表示选取i的代价。如果选取i，定义x[i]=1否则x[i]=0。每一个物品只有选或者不选两种方案，求一个选择方案使得R=sigma(a[i]*x[i])/sigma(b[i]*x[i])取得最值，即所有选择物品的总收益/总代价的值最大或是最小。

分析(详)：

大体概括一下，就是定义一个函数F(L) = sigma(a[i]*x[i])-L*sigma(b[i]*x[i])，其实就是将R换为自变量L后上述式子的变形，分离参数，得到F(L)=sigma((a[i]-L*b[i])*x[i])。如果我们固定L在某个方案下能导致F(L)大于0，这代表存在更优的L，所以我们不断变换L，寻找一个临界L使得不存在一种方案，能够使F(L)>0。另d[i]=(a[i]-L*b[i])，显然d数组是随着L的增大而单调减的，所以根据d数组去判断边界，如果d数组的相应方案的和>=0则表明可能存在更优解，否则不存在。

所以根据性质，很自然想到用二分去求解。二分的下界一般都是取为0(a[i],b[i]一般都是大于0的)，但是上界怎么求呢？假设两个方案a1,b1和a2,b2。它们分别的比率是a1/b1,a2/b2，组合的比率是(a1+a2)/(b1+b2)

a1/b1 - (a1+a2)/(b1+b2)化简后是(a1b2-a2b1)/(b1(b1+b2))

a2/b2 - (a1+a2)/(b1+b2)化简后是(a2b1-a1b2)/(b2(b1+b2))

由于正负取决于分子，所以可以看出两者组合的比率必定小于其中一个单独的比率。所以上界就是寻找一个最大的单独比率。

二分是一个非常通用的办法，但是二分的时候我们只是用到了F(L)>0这个条件，而对于使得F(L)>0的这组解所求到的R值没有使用。因为F(L)>0,我们已经知道了R是一个更优的解，与其漫无目的的二分，为什么不将解移动到R上去呢？求01分数规划的另一个方法就是Dinkelbach算法，它就是基于这样的一个思想，他并不会去二分答案，而是先随便给定一个答案，然后根据更优的解不断移动答案，逼近最优解。由于对每次判定使用的更加充分，所以它比二分会快上很多。但是，它的弊端就是需要保存这个解，而我们知道，有时候验证一个解和求得一个解的复杂度是不同的。二分和Dinkelbach算法写法都非常简单，各有长处，要根据题目谨慎使用。

POJ-2976题意：

给定A数组B数组，从中选择N-K个使得R最大，输出Round(100*R);

代码1(二分方法 79ms)：

#include <algorithm>#include <cstdio>#define LL long longusing namespace std;const int maxn = 1005;const double eps = 1e-6;LL a[maxn], b[maxn];double d[maxn];int n, k;LL _round(double x){if(x-(LL)x < 0.5) return (LL)x;else return (LL)x+1;}bool jg(double v){double tmp = 0;for(int i = 1; i <= n; ++i)d[i] = a[i]-b[i]*v;sort(d+1, d+n+1);for(int i = 0; i < k; ++i)tmp += d[n-i];return tmp >= 0;}void work(){double l = 0, r = 0, mid;for(int i = 1; i <= n; ++i)if(1.0*a[i]/b[i] > r) r = 1.0*a[i]/b[i];while(r-l >= eps){mid = (l+r)/2.0;if(jg(mid)) l = mid;else r = mid;}printf("%lld\n", _round(l*100));}int main(){while(scanf("%d %d", &n, &k), (n||k)){for(int i = 1; i <= n; ++i)scanf("%lld", &a[i]);for(int i = 1; i <= n; ++i)scanf("%lld", &b[i]);k = n-k; work();}return 0;}

代码2(Dinkelbach算法 32ms)

#include <algorithm>#include <cstdio>#include <cmath>#define LL long longusing namespace std;const int maxn = 1005;const double eps = 1e-6;struct node{int a, b;double d;} p[maxn];int cmp(node x, node y) {return x.d < y.d;}int n, k;LL _round(double x){if(x-(LL)x < 0.5) return (LL)x;else return (LL)x+1;}void work(){double ans = 1.0, tmp, x, y;do{tmp = ans;for(int i = 1; i <= n; ++i)p[i].d = p[i].a-p[i].b*tmp;sort(p+1, p+n+1, cmp);x = y = 0.0;for(int i = 0; i < k; ++i)x += p[n-i].a, y += p[n-i].b;ans = x/y;}while(fabs(tmp-ans) >= eps);printf("%lld\n", _round(ans*100));}int main(){while(scanf("%d %d", &n, &k), (n||k)){for(int i = 1; i <= n; ++i)scanf("%d", &p[i].a);for(int i = 1; i <= n; ++i)scanf("%d", &p[i].b);k = n-k; work();}return 0;}

继续加油~