首页 Web开发 Windows程序 编程语言 数据库 移动开发 系统相关 微信 其他好文

POJ 3621--Sightseeing Cows(0-1规划求最大密度)

时间：2015-01-30 09:06:51   收藏：0   阅读：338

题意：给定一个<V,E>规模分别为L和P的有向图，每个点有点权，每个有向边有边权，请找出一个环路，使此环的点权和与边权和的除积最大。

题解：此类求两个集合和的最大除积，类似于0-1规划问题可以先转化为线性式子，然后利用二分或者Dinkelbach算法求得最大或者最小值。

  1. 当y < y_max，对于最大除积的环路，Z(y_max) < 0，即我们利用新生成的边权w_new = W[u][v]*y-U[v]，可以推出此环为一负环，对于其他环路可能存在负环，也可能不存在。
  2. 当y >= y_max，根据函数式可以推出所有环路都是不存在负环的。
注意：负环可以用spfa和bellman_ford算法，可以首先初始化所有点距离为0，即把所有点都当成初始点在全图中寻找负环，加入简化的并查集优化可以加快速度。
  • Dinkelbach算法：
  1. 假设图中有k个环，其中的除积按升序排列为：y1 < y2 < ... < yk。
  2. 则我们初始化一个可知的下界y0 < y1，显然对于所有环路都是能得到负环的，而在寻找负环的过程中我们是可以利用并查集加上点权和数组找出这个负环的U和W的，也就可以找到这个环路所对应的yi = U/W。
  3. 注意到我们每次找到负环求得的yi是逐渐递增直到yk，但是每次出现哪个负环是不可知，这样就能在平均意义上加快搜索速度。
  • 负环搜索次数简单分析：假设剩下n个负环，接下来我们搜索到其中每个负环的概率都是相等的。令x[i]为遇到yi后还需要搜索的期望次数，i = 0,1,2,...,k，其中x[0]就是我们要求的平均次数。
  1. 易知x[k] = 0，x[k-1] = 1。
  2. x[i] = {x[i+1]+1+x[i+2]+1+...+x[k]+1}/(k-i) = 1+sum{x[j],j = i+1,i+2,...,k}/(k-i)。
  3. 得到递推公式：x[i]-x[i+1] = 1/(k-i)。
  4. x[0] = 1+1/2+1/3+...+1/k。前k个倒数和近似于k的对数复杂度。
注意：因为精度问题，若直接用yi作为下一个求负环的输入，则仍然有可能得到此负环的除积yi，而导致错误的结果，可以加上一个小于输出精度的数，使上一个负环不成立并且不影响最终结果。

bellman_ford&二分：

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;
#define maxN 1005
#define maxP 5005
#define eps (1e-4)

struct node
{
    short x,y;
    short val;
};

class solve
{
private:
    short funVal[maxN];
    node edges[maxP];
    int L,P;
    double maxFunUnit;
public:
    solve(int l,int p):L(l),P(p)
    {
        processIn();
        maxFunUnit = two_partion();
        if(maxFunUnit != -1)
        {
            printf("%.2lf\n",maxFunUnit);
        }
        else
        {
            printf("0.00\n");
        }
    }
    int processIn();
    double two_partion();
    int bellman_ford(double funUnit);
};

int solve::bellman_ford(double funUnit)
{
    short fa[maxN];
    double dis[maxN];
    double tmpDis;
    int i,j;
    char IsRelax;
    short x,y;
    for(i = 1;i <= L;i++)
    {
        dis[i] = 1e9;
    }
    for(i = 1;i <= L;i++)
    {
        fa[i] = i;
    }
    for(i = 0;i < L-1;i++)
    {
        IsRelax = false;
        for(j = 0;j < P;j++)
        {
            x = edges[j].x;
            y = edges[j].y;
            tmpDis = dis[x]+((double)edges[j].val)*funUnit-(double)funVal[y];
            if(dis[y] > tmpDis)
            {
                if(fa[x] == y)
                    return 1;
                dis[y] = tmpDis;
                fa[y] = fa[x];
                IsRelax = true;
            }
        }
        if(!IsRelax)
            return 0;
    }
    for(j = 0;j < P;j++)
    {
        x = edges[j].x;
        y = edges[j].y;
        tmpDis = dis[x]+((double)edges[j].val)*funUnit-(double)funVal[y];
        if(dis[y] > tmpDis)
        {
            return 1;
        }
    }
    return 0;
}

double solve::two_partion()
{
    if(!bellman_ford(-1.0))
        return -1;
    double left,right,mid;
    left = 1e-3;
    right = 1000;
    while(right-left > eps)
    {
        mid = (left+right)/2;
        if(bellman_ford(mid))      //有负环
        {
            left = mid;
        }
        else
        {
            right = mid;
        }
    }
    return right;
}

int solve::processIn()
{
    int i;
    for(i = 1;i <= L;i++)
    {
        scanf("%d",funVal+i);
    }
    for(i = 0;i < P;i++)
    {
        scanf("%d%d%d",&(edges[i].x),&(edges[i].y),&(edges[i].val));
    }
    return 0;
}

int main()
{
    int l,p;
    while(~scanf("%d%d",&l,&p))
    {
        solve sightseeing(l,p);
    }
    return 0;
}

spfa&二分：

#include<iostream>
#include<cstdio>
#include<cstring>
#include<vector>
#include<queue>
using namespace std;
#define maxN 1005
#define eps (1e-4)

struct node
{
    int adjNode;
    int val;
};

class solve
{
private:
    short funVal[maxN];
    int L,P;
    double maxFunUnit;
    vector<vector<node> > adj;
public:
    solve(int l,int p):L(l),P(p)
    {
        processIn();
        maxFunUnit = two_partion();
        if(maxFunUnit != -1)
        {
            printf("%.2lf\n",maxFunUnit);
        }
        else
        {
            printf("0.00\n");
        }
    }
    int processIn();
    double two_partion();
    int spfa(double funUnit);
};

int solve::spfa(double funUnit)
{
    char vis[maxN];
    short fa[maxN];
    short depth[maxN];
    double dis[maxN];
    double tmpDis;
    int adjSize;
    int i;
    short now,next;
    queue<short> q;
    for(i = 1;i <= L;i++)
    {
        dis[i] = 1e9;
        q.push(i);
        fa[i] = i;
    }
    memset(vis,1,sizeof(vis));
    memset(depth,0,sizeof(depth));
    while(!q.empty())
    {
        now = q.front();
        q.pop();
        vis[now] = 0;
        adjSize = adj[now].size();
        for(i = 0;i < adjSize;i++)
        {
            next = adj[now][i].adjNode;
            tmpDis = dis[now]+((double)(adj[now][i].val))*funUnit-(double)funVal[next];
            if(dis[next] > tmpDis)
            {
                if(fa[now] == next)
                    return 1;
                dis[next] = tmpDis;
                fa[next] = fa[now];
                if(!vis[next])
                {
                    vis[next] = 1;
                    q.push(next);
                    depth[next]++;
                    if(depth[next] > L)
                    {
                        return 1;
                    }
                }
            }
        }
    }
    return 0;
}

double solve::two_partion()
{
    if(!spfa(-1.0))
        return -1;
    double left,right,mid;
    left = 1e-3;
    right = 1000;
    while(right-left > eps)
    {
        mid = (left+right)/2;
        if(spfa(mid))      //有负环
        {
            left = mid;
        }
        else
        {
            right = mid;
        }
    }
    return right;
}

int solve::processIn()
{
    int i;
    int a;
    node tmpNode;
    for(i = 1;i <= L;i++)
    {
        scanf("%d",funVal+i);
    }
    adj.resize(L+1);
    for(i = 0;i < P;i++)
    {
        scanf("%d%d%d",&a,&(tmpNode.adjNode),&(tmpNode.val));
        adj[a].push_back(tmpNode);
    }
    return 0;
}

int main()
{
    int l,p;
    while(~scanf("%d%d",&l,&p))
    {
        solve sightseeing(l,p);
    }
    return 0;
}

spfa&Dinkelbach：

#include<iostream>
#include<cstdio>
#include<cstring>
#include<vector>
#include<queue>
#include<cmath>
using namespace std;
#define maxN 1005
#define eps (1e-3)

struct node
{
    int adjNode;
    int val;
};

class solve
{
private:
    short funVal[maxN];
    double sumVal[maxN];
    int L,P;
    double maxFunUnit;
    double nextFunUnit;
    vector<vector<node> > adj;
public:
    solve(int l,int p):L(l),P(p)
    {
        processIn();
        maxFunUnit = two_partion();
        printf("%.2lf\n",maxFunUnit);
    }
    int processIn();
    double two_partion();
    int spfa(double funUnit);
};

int solve::spfa(double funUnit)
{
    char vis[maxN];
    short fa[maxN];
    short depth[maxN];
    double dis[maxN];
    double tmpDis;
    int adjSize;
    int i;
    short now,next;
    queue<short> q;
    for(i = 1;i <= L;i++)
    {
        q.push(i);
        fa[i] = i;
    }
    memset(dis,0,sizeof(dis));
    memset(sumVal,0,sizeof(sumVal));
    memset(vis,1,sizeof(vis));
    memset(depth,0,sizeof(depth));
    while(!q.empty())
    {
        now = q.front();
        q.pop();
        vis[now] = 0;
        adjSize = adj[now].size();
        for(i = 0;i < adjSize;i++)
        {
            next = adj[now][i].adjNode;
            tmpDis = dis[now]+((double)(adj[now][i].val))*funUnit-(double)funVal[next];
            if(dis[next] > tmpDis)
            {
                if(fa[now] == next)
                {
                    nextFunUnit = funUnit/((tmpDis/(sumVal[now]+funVal[next]))+1.0);
                    //计算此负环的除积
                    return 1;
                }
                dis[next] = tmpDis;
                fa[next] = fa[now];
                sumVal[next] = sumVal[now]+funVal[next];
                //记录从祖先结点到此节点的所有新生成的边权和
                if(!vis[next])
                {
                    vis[next] = 1;
                    q.push(next);
                    depth[next]++;
                    if(depth[next] > L)
                    {
                        nextFunUnit = funUnit/((tmpDis/(sumVal[now]+funVal[next]))+1.0);
                        return 1;
                    }
                }
            }
        }
    }
    return 0;
}

double solve::two_partion()
{
    double left,pre_Left;
    int tmpRe;
    left = 1e-3;
    pre_Left = 100;
    while(fabs(left-pre_Left) > eps)
    {
        tmpRe = spfa(left);
        pre_Left = left;
        if(tmpRe)
        {
            left = nextFunUnit+0.00001;		
            //加上精度防止再找到此负环
        }
    }
    return left;
}

int solve::processIn()
{
    int i;
    int a;
    node tmpNode;
    for(i = 1;i <= L;i++)
    {
        scanf("%d",funVal+i);
    }
    adj.resize(L+1);
    for(i = 0;i < P;i++)
    {
        scanf("%d%d%d",&a,&(tmpNode.adjNode),&(tmpNode.val));
        adj[a].push_back(tmpNode);
    }
    return 0;
}

int main()
{
    int l,p;
    while(~scanf("%d%d",&l,&p))
    {
        solve sightseeing(l,p);
    }
    return 0;
}


评论（0）
© 2014 mamicode.com 版权所有 京ICP备13008772号-2  联系我们:gaon5@hotmail.com
迷上了代码！