POJ 2728 Desert King

主要参考：http://blog.csdn.net/sdj222555/article/details/7490797

http://www.cnblogs.com/lotus3x/archive/2009/03/21/1418480.html

概念

有带权图G, 对于图中每条边e[i], 都有benifit[i](收入)和cost[i](花费), 我们要求的是一棵生成树T, 它使得 ∑(benifit[i]) / ∑(cost[i]), i∈T 最大(或最小).

这显然是一个具有现实意义的问题.

 

解法之一 0-1分数规划

设x[i]等于1或0, 表示边e[i]是否属于生成树.

则我们所求的比率 r = ∑(benifit[i] * x[i]) / ∑(cost[i] * x[i]), 0≤i<m .

为了使 r 最大, 设计一个子问题---> 让 z = ∑(benifit[i] * x[i]) - l * ∑(cost[i] * x[i]) = ∑(d[i] * x[i]) 最大 (d[i] = benifit[i] - l * cost[i]) , 并记为z(l). 我们可以兴高采烈地把z(l)看做以d为边权的最大生成树的总权值.

然后明确两个性质:

　1.  z单调递减

　　证明: 因为cost为正数, 所以z随l的减小而增大.

　2.  z( max(r) ) = 0

　　证明: 若z( max(r) ) < 0, ∑(benifit[i] * x[i]) - max(r) * ∑(cost[i] * x[i]) < 0, 可化为 max(r) < max(r). 矛盾;

　　        若z( max(r) ) >= 0, 根据性质1, 当z = 0 时r最大. 

其实原题就是求 MIN（ ∑CiXi / ∑DiXi ） Xi∈{0，1} ，对每个生成树，设其比率r=∑CiXi / ∑DiXi ，可得∑CiXi - ∑DiXi * r=0（条件1）

那么对于所有的生成树，显然∑CiXi - ∑DiXi * min(r) >= 0，当 ∑CiXi / ∑DiXi = min(r)时，等号成立。 而我们现在不知道min(r)是多少，只好进行枚举，对每个枚举的r ,构建新的权值（Ci-Di*r），然后求最小生成树，  为什么求最小呢？ 我的理解就是这是为了寻找使得生成树的总权值为0的可能性，因为只有当其等于0 的时候，才满足了条件1 这个条件， 说明这个r是可行的，并且如果r枚举到值为min(r)时，其最小生成树的的总权值必然恰好等于0，但是如果不能等于0， 比如大于0， 显然是对该r值，所有的生成树上无论如何也满足不了条件1，说明r值就是偏小了。同理如果小于0，r值是偏大的,说明可能存在某些生成树使得满足条件1，而我们的目标是在满足条件1的情况下使得r最小。

根据这个我们可以发现，实际上r的值是可以进行二分查找的。 而也有人给出了更为高效的迭代方法。

第一种是二分的方法

#include <iostream>
#include <cstdio>
#include <cstring>
#include <vector>
#include <stack>
#include <cmath>
#define INF 1000000000
#define eps 1e-7
using namespace std;
const int N = 1005;
int n;
double e[N][N],low[N];
int  nearn[N];
struct point{
	double x,y,z;
}p[N];
double getdis(point a, point b){
	return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
 bool prim(double l){
    double sum = 0;
    for(int i = 1; i <= n; i++){
        nearn[i] = 1;
        low[i] = abs(p[1].z - p[i].z) - e[1][i] * l;
    }
    nearn[1] = -1;
    for(int i = 1; i < n; i++){
        double mi = INF;
        int v = -1;
        for(int j = 1; j <= n; j++)
            if(nearn[j] != -1 && low[j] < mi){
                v = j;
                mi = low[j];
            }
        if(v == -1)continue;
		nearn[v] = -1;
		sum += low[v];
		for(int j = 1; j <= n; j++){
			double tmp = abs(p[v].z - p[j].z) - e[v][j] * l;
			if(nearn[j] != -1 && tmp < low[j])
			{
				low[j] = tmp;
				nearn[j] = v;
			}
		}
    }
    return sum >= 0;
}
int main() {
    while(~scanf("%d", &n) && n) {
        for(int i = 1; i <= n; i++)
            scanf("%lf%lf%lf", &p[i].x, &p[i].y, &p[i].z);
        for(int i = 1; i <= n; i++)
            for(int j = 1; j <= n; j++)
                e[i][j] = getdis(p[i], p[j]);
        double l = 0.0, r = 50.0, mid;
        while(r - l > eps) {
            mid = (l + r) / 2;
            if(prim(mid)) l = mid;
            else r = mid;
        }
        printf("%.3f\n", r);
    }
    return 0;
}

第二种是迭代

#include <iostream>
#include <cstdio>
#include <cstring>
#include <vector>
#include <stack>
#include <cmath>
#define INF 1000000000
#define eps 1e-7
using namespace std;
const int N = 1005;
int n;
double e[N][N],low[N];
int  nearn[N];
struct point{
	double x,y,z;
}p[N];
double getdis(point a, point b){
	return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
double prim(double l){
    double cost = 0, len = 0;
    double sum = 0;
    for(int i = 1; i <= n; i++){
        nearn[i] = 1;
        low[i] = abs(p[1].z - p[i].z) - e[1][i] * l;
    }
    nearn[1] = -1;
    for(int i = 1; i < n; i++){
        double mi = INF;
        int v = -1;
        for(int j = 1; j <= n; j++)
            if(nearn[j] != -1 && low[j] < mi){
                v = j;
                mi = low[j];
            }
        if(v == -1)continue;
		cost += abs(p[nearn[v]].z - p[v].z);
		len += e[nearn[v]][v];
		nearn[v] = -1;
		sum += low[v];
		for(int j = 1; j <= n; j++){
			double tmp = abs(p[v].z - p[j].z) - e[v][j] * l;
			if(nearn[j] != -1 && tmp < low[j])
			{
				low[j] = tmp;
				nearn[j] = v;
			}
		}
    }
    return cost/len;
}
int main() {
    while(~scanf("%d", &n) && n) {
        for(int i = 1; i <= n; i++)
            scanf("%lf%lf%lf", &p[i].x, &p[i].y, &p[i].z);
        for(int i = 1; i <= n; i++)
            for(int j = 1; j <= n; j++)
                e[i][j] = getdis(p[i], p[j]);
        double a = 0, b;
		while(true){
            b = prim(a);
            if(fabs(a-b) < eps) break;
            a = b;
        }
        printf("%.3f\n", a);
    }
    return 0;
}