最优比例生成树(0/1分数规划)

首先这是要解决什么问题:一个带权完全图,每条边都有自己的花费值cost[i]和收益值benifit[i],如果用x[i]来代表一条边取或不取,那么求一个生成树。要求:r=(∑cost[i]*x[i] ) / (∑benifit[i]*x[i] )最小。

经典题目:POJ2728 - Desert King

如何来求解:这里用到了0-1分数规划思想,对于上式可以变形为 z(r)=∑cost[i]*x[i] -r*∑benifit[i]*x[i]。而z(r)=0为我们所求。

这里有个非常重要的结论:z(r)为单调递减函数,因此是线性的。于是"我们可以兴高采烈地把z(r)看做以 cost[i]-r*benifit[i] 为边权的最小生成树的总权值"。显然,我们所求即为z(max(r))=0。这里max(r)是由z=0确定的。

于是有了两种算法:

1、Dinkelbanch算法,即迭代法。取r的一个初值,一般为0,得到一个生成树之后,令r=(∑benifit[i]*x[i] ) / (∑cost[i]*x[i] ),哐哐迭代之后……r是趋于最大值的,而z趋于0。

2、二分法,z(r)=0是我们所求,而它又是单调函数,那么二分r即可。据说二分比迭代慢了很多。

 

 

模板Dinkelbach POJ 2728

#include <iostream>
#include <cstdio>
#include <cmath>
#include <vector>
#include <cstring>
#include <algorithm>
#include <string>
#include <set>
#include <ctime>
#include <queue>
#include <map>
#include <sstream>

#define CL(arr, val)    memset(arr, val, sizeof(arr))
#define REP(i, n)       for((i) = 0; (i) < (n); ++(i))
#define FOR(i, l, h)    for((i) = (l); (i) <= (h); ++(i))
#define FORD(i, h, l)   for((i) = (h); (i) >= (l); --(i))
#define L(x)    (x) << 1
#define R(x)    (x) << 1 | 1
#define MID(l, r)   (l + r) >> 1
#define Min(x, y)   x < y ? x : y
#define Max(x, y)   x < y ? y : x
#define E(x)    (1 << (x))

const double eps = 1e-6;
typedef long long LL;
using namespace std;

const int N = 1<<10;
const int inf = ~0u>>2;

struct node {
    int x, y, z;
} point[N];

double cost[N][N];
double benifit[N][N];
double mi[N];
int pre[N];
bool vis[N];
double c, d;
int n;

double prim(double ans) {
    int i, j, f;
    double mx;
    c = d = 0;
    for(i = 1; i < n; ++i) {
        mi[i] = cost[0][i] - ans*benifit[0][i];
        pre[i] = 0;
        vis[i] = false;
    }
    mi[0] = 0;
    vis[0] = true;
    for(i = 0; i < n - 1; ++i) {
        mx = inf; f = 0;
        for(j = 0; j < n; ++j) {
            if(!vis[j] && mx > mi[j]) {
                mx = mi[j];
                f = j;
            }
        }
        c += cost[pre[f]][f];
        d += benifit[pre[f]][f];
        vis[f] = true;

        for(j = 1; j < n; ++j) {
            double val = cost[f][j] - ans*benifit[f][j];
            if(!vis[j] && mi[j] > val) {
                pre[j] = f;
                mi[j] = val;
            }
        }
    }
    return c/d;
}

int iabs(int x) {
    return x < 0 ? -x : x;
}

void init() {
    int i, j;
    for(i = 0; i < n; ++i) {
        for(j = 0; j < n; ++j) {
            cost[i][j] = iabs(point[i].z - point[j].z);
            benifit[i][j] = sqrt(1.*(point[i].x - point[j].x)*(point[i].x - point[j].x) +
                                 1.*(point[i].y - point[j].y)*(point[i].y - point[j].y));
        }
    }
}

void solve() {
    double ans = 0, tmp;
    while(1) {
        tmp = prim(ans);
        if(fabs(ans - tmp) < eps)   break;
        else    ans = tmp;
    }
    printf("%.3f\n", ans);
}

int main() {
    freopen("data.in", "r", stdin);

    int i;
    while(scanf("%d", &n), n) {
        for(i = 0; i < n; ++i) {
            scanf("%d%d%d", &point[i].x, &point[i].y, &point[i].z);
        }
        init();
        solve();
    }
    return 0;
}
posted @ 2012-07-17 21:44  AC_Von  阅读(702)  评论(0编辑  收藏  举报