Desert King POJ - 2728 (最优比率生成树) 二分

题意：
在这么一个图中求一棵生成树,这棵树的单位长度的花费最小是多少？
思路
题目可以写成如下公式：其中E指的是 生成树中的边

R=sigma(a[i]*x[i])/sigma(b[i]*x[i])

我们先定义一个函数F(L):=sigma(a[i]*x[i])-L*sigma(b[i]*x[i])，显然这只是对目标式的一个简单的变形。分离参数，得到F(L):=sigma((a[i]-L*b[i])*x[i])。这时我们就会发现，如果L已知的话，a[i]-L*b[i]就是已知的，当然x[i]是未知的。记d[i]=a[i]-L*b[i]，那么F(L):=sigma(d[i]*x[i])

L就是目标式中的R，最大化R也就是最大化L。

F(L)=sigma(a[i]*x[i])-L*sigma(b[i]*x[i])>0即sigma(a[i]*x[i])/sigma(b[i]*x[i])>L也就是说，如果一个方案使得F(L)>0说明了这组方案可以得到一个比现在的L更优的一个L，既然有一个更优的解，那么为什么不用呢？

显然，d数组是随着L的增大而单调减的。也就是说，存在一个临界的L使得不存在一种方案，能够使F(L)>0. 当F(L)=0使，对应方案的R值恰好等于此时的L值。

可以用二分 ：

求01分数规划的另一个方法就是Dinkelbach算法
它本质上是观察直线
R=cx’-Ldx’
于函数F(L)在L=L’处相切，这里x’是子问题Q(L’)的最优解。因此，-dx’是F(L)在L’处的斜率。而且很容易看出上面的直线与L轴相交与L=cx’/dx’.

Dinkelbach算法：
步骤1：设L=L1,使 L*<=L1<=nC
步骤2：解决子问题Q(L)并得到最优解x
步骤3：如果z(L)=0，那么输出x并终止。否则，设L=cx/dx跳到步骤2
为了初始化L1,将用到nC，因此充分挖掘拓展问题的结构将能做出更好的选择。

二分代码

/*
 * Author       :  Echo
 * Email        :  1666424499@qq.com  
 * Description  :   
 * Created Time :  2017/10/10 19:25:26
 * Last Modify  :  2017/10/12 20:11:14
 * File Name    :  \yanga11ang\write.cpp
 */
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <string>
#include <queue>
using namespace std;
const int maxn=1100;
const double INF=1e9;
struct node{
    double x,y,h;
}an[maxn];
double cx[maxn][maxn];
double dx[maxn][maxn];
double mp[maxn][maxn];
bool vis[maxn];
double dis[maxn];
int n;
double _abs(double a){
    return a>=0 ? a:-a;
}
double prim(){//普利姆算法 最短路
    for(int i=1;i<=n;i++) 
        dis[i]=mp[1][i],vis[i]=0;
    vis[1]=1;
    double res=0;
    for(int i=1;i<n;i++){
        double minnum=INF;
        int u;
        for(int j=2;j<=n;j++){
            if(vis[j])continue;
            if(minnum>dis[j]){
                minnum=dis[j];
                u=j;
            }
        }
        res+=dis[u];
        vis[u]=1;
        for(int j=2;j<=n;j++){
            if(vis[j])continue;
            if(mp[u][j]<dis[j]){
                dis[j]=mp[u][j];
            }
        }
    }
    return res;
}
double dist(node &a,node &b){
    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)); 
}
int main(){
    while(scanf("%d",&n),n){
        for(int i=1;i<=n;i++){
            scanf("%lf%lf%lf",&an[i].x,&an[i].y,&an[i].h);
        }
        double l=1e-9,r=1e9,mid;
        for(int i=1;i<n;i++)
            for(int j=i+1;j<=n;j++){
                cx[i][j]=cx[j][i]= _abs(an[i].h-an[j].h);
                dx[i][j]=dx[j][i]=dist(an[i],an[j]);
            }
        while(r-l>1e-4){
            mid=(l+r)/2;
            for(int i=1;i<n;i++)
                for(int j=i+1;j<=n;j++)
                    mp[i][j]=mp[j][i]=cx[i][j]-mid*dx[i][j];
            if(prim()<0) r=mid;
            else l=mid;
        }
        printf("%.3f\n",l);
    }
    return 0;
}
/*

4
0 0 0
0 1 1
1 1 2
1 0 3


1.000
 */