poj 2516--Minimum Cost

32 篇文章 0 订阅
30 篇文章 0 订阅

最小费用最大流问题(看了一周《算法导论》的单纯性算法已经长跪不起了!):在最大流的基础上,使各个路径的花费和最小。

第一种方法:参考最大流

简要分析一下两种问题的异同,最大流问题建图时只有有向边的容量capacity,而最小费用最大流还有有向边的单位容量的花费cost,所以找增广路径时略有不同。

  1. 简单地说,最大流的增广路径的可行边是有流量的边就可以走,能到终点就是增广路径。
  2. 而最小费用最大流的增广路径必须是:有容量且单位容量花费和最小的路径
为什么是单位容量最小花费和,而不是增广路径最小总花费和(路径最小容量x路径单位容量花费和),其实易知第二个值是不容易求的,因为它不是最短路而是最小乘法加权值之和。具体原理待证明。
算法具体实现:对每件商品建图求最小费用最大流
  1. 建图:新建超级源节点到所有供应商的容量为供应商的存储量,同时供应商到商店的容量初值也设置为供应商存储量,新建超级汇点,商店到汇点容量为商店需求量;供应商到商店花费为单位商品运费,同时取其相反数赋给商店到供应商花费(相当于退货)。
  2. 找图中增广路径,若失败则说明已最优;反之根据路径更改容量(同最大流方法)。
求最小单位容量花费路径可以使用SPFA,Bellman_Ford。

第二种方法: KM算法
KM算法是一种求最大加权二分匹配的方法,运用于此题,首先注意到二分匹配图中权值可以等效于单位容量花费,但是其中没有边的容量概念。但是我们可以把商店和供应商分解成若干(需求量或者供应量)相当于容量为1的小商店以及小供应商,就可以等价于求二分最大加权匹配了。
算法具体实现:
  1. 根据需求量把商店分解到二分图中X集合中的点,同理根据供应量把供应商分解为二分图中Y集合中的点,各边权值仍然为原商店与原供应商的单位花费的相反数。可知X中的点一定是有匹配的,即需求是必须满足的,但Y中的点不一定有匹配(供大于求)。
  2. 根据X中点序,直到能找到一条增广路径再继续X中下一个点,若找不到则修改左右标签值。

SPFA:
#include<iostream>
#include<cstring>
#include<cstdio>
#include<queue>
using namespace std;
#define INF 0X7F7F
#define maxN 55
#define maxK 55

int N,M,K;
short need[maxK][maxN];
short storage[maxK][maxN];
short capacity[(maxN)<<1][(maxN)<<1];
short cost[(maxN)<<1][(maxN)<<1];
short minCost;
short dest;

short spfa(short* fa)
{
    short i,now;
    short minDis[(maxN)<<1];
    char vis[(maxN)<<1];
    memset(vis,0,sizeof(vis));
    memset(minDis,127,sizeof(minDis));
    minDis[0] = 0;
    queue<short> q;
    q.push(0);
    vis[0] = true;
    while(!q.empty())
    {
        now = q.front();
        q.pop();
        vis[now] = false;
        for(i = 1;i <= dest;i++)
        {
            if(capacity[now][i] > 0&&minDis[now]+cost[now][i] < minDis[i])     //找出有容量且单位商品的最小花费路径
            {
                fa[i] = now;
                minDis[i] = minDis[now]+cost[now][i];
                if(!vis[i])
                {
                    q.push(i);
                    vis[i] = true;
                }
            }
        }
    }
    return minDis[dest];
}

int minCost_maxFlow()
{
    short x;
    short minCapacity = INF;
    short fa[(maxN)<<1];
    memset(fa,0,sizeof(fa));
    while(spfa(fa) != INF)
    {
        x = dest;
        while(x)       //找出增广路径的最小容量
        {
            minCapacity = min(minCapacity,capacity[fa[x]][x]);
            x = fa[x];
        }
        x = dest;
        while(x)    //修改增广路径中容量
        {
            minCost += minCapacity*cost[fa[x]][x];
            capacity[fa[x]][x] -= minCapacity;
            capacity[x][fa[x]] += minCapacity;
            x = fa[x];
        }
        memset(fa,0,sizeof(fa));
    }
    return 0;
}

int main()
{
    int i,j,k;
    short sumNeed[maxK];
    short sumStorage[maxK];
    char needEnough;
    while(~scanf("%d%d%d",&N,&M,&K)&&N)
    {
        minCost = 0;
        dest = N+M+1;
        needEnough = true;
        memset(cost,0,sizeof(cost));
        memset(sumNeed,0,sizeof(sumNeed));
        memset(sumStorage,0,sizeof(sumStorage));
        for(i = 1;i <= N;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&need[k][i]);   //第k个商家需要多少第j种商品
                sumNeed[k] += need[k][i];
            }
        }
        for(i = 1;i <= M;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&storage[k][i]);  //第k个供货商有多少j的库存
                sumStorage[k] += storage[k][i];
            }
        }
        for(k = 1;k <= K;k++)
        {
            if(sumNeed[k] > sumStorage[k])     //第k中商品的总供应量少于需求
            {
                needEnough = false;
                break;
            }
        }
        for(k = 1;k <= K;k++)
        {
            memset(capacity,0,sizeof(capacity));
            for(i = 1;i <= N;i++)
            {
                capacity[M+i]dest] = need[k][i];      //商店到汇点的容量为商店的需求
                for(j = 1;j <= M;j++)
                {
                    capacity[j][M+i] = capacity[0][j] = storage[k][j];  //源点到供应商以及供应商到商店的容量为供应商存储量
                    scanf("%d",&cost[j][M+i]);
                    cost[M+i][j] = -1*cost[j][M+i];     //商店到供应商的花费为供应商到商店花费的相反数(相当于退回商品)
                }
            }
            if(needEnough)
            minCost_maxFlow();
        }
        if(!needEnough)
        {
            printf("-1\n");
        }
        else
        printf("%d\n",minCost);
    }
    return 0;
}

Bellman_Ford:
#include<iostream>
#include<cstring>
#include<cstdio>
#include<queue>
using namespace std;
#define INF 0X7F7F
#define maxN 55
#define maxK 55
#define maxE 5100

int N,M,K;
short need[maxK][maxN];
short storage[maxK][maxN];
short capacity[(maxN)<<1][(maxN)<<1];
short cost[(maxN)<<1][(maxN)<<1];
short edges[maxE][2];
short edgeNum;
short minCost;
short dest;

int generateEdge()     //每次容量改变后都重新生成边集合
{
    edgeNum = 0;
    short i,j;
    for(i = 0;i < dest;i++)
    {
        for(j = 1;j <= dest;j++)
        {
            if(capacity[i][j] > 0)
            {
                edges[edgeNum][0] = i;
                edges[edgeNum++][1] = j;
            }
        }
    }
    return 1;
}

short bellman_ford(short* fa)
{
    short i,j;
    short x,y;
    short minDis[(maxN)<<1];
    char IsRelax;
    memset(minDis,0x7F,sizeof(minDis));
    minDis[0] = 0;
    for(i = 0;i < dest;i++)
    {
        IsRelax = false;
        for(j = 0;j < edgeNum;j++)
        {
            x = edges[j][0];
            y = edges[j][1];
            if(minDis[y] > minDis[x]+cost[x][y])
            {
                fa[y] = x;
                minDis[y] = minDis[x]+cost[x][y];
                IsRelax = true;
            }
        }
        if(!IsRelax)
            break;
    }
    return minDis[dest];
}

int minCost_maxFlow()
{
    short x;
    short minCapacity;
    short fa[(maxN)<<1];
    while(generateEdge()&&bellman_ford(fa) != INF)
    {
        minCapacity = INF;
        x = dest;
        while(x)
        {
            minCapacity = min(minCapacity,capacity[fa[x]][x]);
            x = fa[x];
        }
        x = dest;
        while(x)
        {
            minCost += minCapacity*cost[fa[x]][x];
            capacity[fa[x]][x] -= minCapacity;
            capacity[x][fa[x]] += minCapacity;
            x = fa[x];
        }
        memset(fa,0,sizeof(fa));
    }
    return 0;
}

int main()
{
    int i,j,k;
    short sumNeed[maxK];
    short sumStorage[maxK];
    char needEnough;
    while(~scanf("%d%d%d",&N,&M,&K)&&N)
    {
        minCost = 0;
        dest = N+M+1;
        needEnough = true;
        memset(cost,0,sizeof(cost));
        memset(sumNeed,0,sizeof(sumNeed));
        memset(sumStorage,0,sizeof(sumStorage));
        for(i = 1;i <= N;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&need[k][i]);   //第k个商家需要多少第j种商品
                sumNeed[k] += need[k][i];
            }
        }
        for(i = 1;i <= M;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&storage[k][i]);  //第k个供货商有多少j的库存
                sumStorage[k] += storage[k][i];
            }
        }
        for(k = 1;k <= K;k++)
        {
            if(sumNeed[k] > sumStorage[k])
            {
                needEnough = false;
                break;
            }
        }
        for(k = 1;k <= K;k++)
        {
            memset(capacity,0,sizeof(capacity));
            for(i = 1;i <= N;i++)
            {
                capacity[M+i][M+N+1] = need[k][i];
                for(j = 1;j <= M;j++)
                {
                    capacity[j][M+i] = capacity[0][j] = storage[k][j];
                    scanf("%d",&cost[j][M+i]);
                    cost[M+i][j] = -1*cost[j][M+i];
                }
            }
            if(needEnough)
            minCost_maxFlow();
        }
        if(!needEnough)
        {
            printf("-1\n");
        }
        else
        printf("%d\n",minCost);
    }
    return 0;
}



双向循环链表&Bellman_Ford(不用每次都生成边的集合):

#include<iostream>
#include<cstring>
#include<cstdio>
#include<queue>
using namespace std;
#define INF 0X7F7F
#define maxN 55
#define maxK 55

struct node
{
    short nodeNo;
    node* prev;
    node* next;
};
node* adj[(maxN)<<1];
int N,M,K;
short need[maxK][maxN];
short storage[maxK][maxN];
short capacity[(maxN)<<1][(maxN)<<1];
short cost[(maxN)<<1][(maxN)<<1];
short minCost;
short dest;

int adjInsert(short a,short b)
{
    node* newNode = new node;
    newNode->nodeNo = b;
    adj[a]->next->prev = newNode;
    newNode->next = adj[a]->next;
    newNode->prev = adj[a];
    adj[a]->next = newNode;
    return 0;
}

int adjInit()
{
    short i,j;
    for(i = 0;i < dest;i++)     //初始化双向循环链表
    {
        adj[i] = new node;
        adj[i]->prev = adj[i]->next = adj[i];
    }
    for(i = 0;i < dest;i++)
    {
        for(j = 1;j <= dest;j++)
        {
            if(capacity[i][j] > 0)
            {
                adjInsert(i,j);
            }
        }
    }
    return 0;
}

int adjDelete(short a,short b)
{
    node* tmpNode = adj[a]->next;
    while(tmpNode->nodeNo != b)
    {
        tmpNode = tmpNode->next;
    }
    tmpNode->prev->next = tmpNode->next;
    tmpNode->next->prev = tmpNode->prev;
    return 0;
}

short bellman_ford(short* fa)
{
    short i,j;
    short x;
    node* tmpNode;
    short minDis[(maxN)<<1];
    char IsRelax;
    memset(minDis,0x7F,sizeof(minDis));
    minDis[0] = 0;
    for(i = 0;i < dest;i++)
    {
        IsRelax = false;
        for(j = 0;j < dest;j++)
        {
            for(tmpNode = adj[j]->next;tmpNode != adj[j];tmpNode = tmpNode->next)
            {
                x = tmpNode->nodeNo;
                if(minDis[x] > minDis[j]+cost[j][x])
                {
                    fa[x] = j;
                    minDis[x] = minDis[j]+cost[j][x];
                    IsRelax = true;
                }
            }
        }
        if(!IsRelax)
            break;
    }
    return minDis[dest];
}

int minCost_maxFlow()
{
    short x;
    short minCapacity;
    short fa[(maxN)<<1];
    while(bellman_ford(fa) != INF)
    {
        minCapacity = INF;
        x = dest;
        while(x)
        {
            minCapacity = min(minCapacity,capacity[fa[x]][x]);
            x = fa[x];
        }
        x = dest;
        while(x)
        {
            minCost += minCapacity*cost[fa[x]][x];
            if((capacity[fa[x]][x] -= minCapacity) == 0)
            {
                adjDelete(fa[x],x);
            }
            if(x != dest)
            {
                adjInsert(x,fa[x]);
            }
            capacity[x][fa[x]] += minCapacity;
            x = fa[x];
        }
    }
    return 0;
}

int main()
{
    int i,j,k;
    short sumNeed[maxK];
    short sumStorage[maxK];
    char needEnough;
    while(~scanf("%d%d%d",&N,&M,&K)&&N)
    {
        minCost = 0;
        dest = N+M+1;
        needEnough = true;
        memset(cost,0,sizeof(cost));
        memset(sumNeed,0,sizeof(sumNeed));
        memset(sumStorage,0,sizeof(sumStorage));
        for(i = 1;i <= N;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&need[k][i]);   //第k个商家需要多少第j种商品
                sumNeed[k] += need[k][i];
            }
        }
        for(i = 1;i <= M;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&storage[k][i]);  //第k个供货商有多少j的库存
                sumStorage[k] += storage[k][i];
            }
        }
        for(k = 1;k <= K;k++)
        {
            if(sumNeed[k] > sumStorage[k])
            {
                needEnough = false;
                break;
            }
        }
        for(k = 1;k <= K;k++)
        {
            memset(capacity,0,sizeof(capacity));
            for(i = 1;i <= N;i++)
            {
                capacity[M+i][M+N+1] = need[k][i];
                for(j = 1;j <= M;j++)
                {
                    capacity[j][M+i] = capacity[0][j] = storage[k][j];
                    scanf("%d",&cost[j][M+i]);
                    cost[M+i][j] = -1*cost[j][M+i];
                }
            }
            if(needEnough)
            {
                adjInit();
                minCost_maxFlow();
            }
        }
        if(!needEnough)
        {
            printf("-1\n");
        }
        else
        printf("%d\n",minCost);
    }
    return 0;
}

KM算法:

#include<iostream>
#include<cstring>
#include<cstdio>
using namespace std;
#define INF 0X7F7F
#define maxN 55
#define maxK 55
#define maxM 155

short cost[maxM][maxM];         //二分图左右权值
short lx[maxM],ly[maxM];        //二分图左右标签
char visx[maxM],visy[maxM];     //增广路中访问情况
short match[maxM];              //二分图匹配情况
short slack[maxM];              //松弛便签
short nx,ny;
short minCost;

bool DFS(short x)
{
    short i,diff;
    visx[x] = true;
    for(i = 1;i <= ny;i++)
    {
        if(!visy[i])
        {
            diff = lx[x]+ly[i]-cost[x][i];
            if(diff == 0)          //增广路中可行边
            {
                visy[i] = true;
                if(match[i] == 0||DFS(match[i]))      //直到最后遇到Y中一点没有匹配才返回为1
                {
                    match[i] = x;
                    return 1;
                }
            }
            else if(slack[i] > diff)     //i不属于增广路中Y的集合,更改其便签松弛量
            {
                slack[i] = diff;
            }
        }
    }
    return 0;
}

int KM()
{
    short i,j;
    short d;
    memset(match,0,sizeof(match));
    memset(lx,0X80,sizeof(lx));
    memset(ly,0,sizeof(ly));
    for(i = 1;i <= nx;i++)          //初始化标签满足权值要求
    {
        for(j = 1;j <= ny;j++)
        {
            if(lx[i] < cost[i][j])
            {
                lx[i] = cost[i][j];
            }
        }
    }
    for(i = 1;i <= nx;i++)
    {
        memset(slack,0X7F,sizeof(slack));      //初始化松弛量为极大值
        while(1)
        {
            memset(visx,0,sizeof(visx));
            memset(visy,0,sizeof(visy));
            if(DFS(i))          //找到增广路径就跳出循环,没找到改标号
            {
                break;
            }
            d = INF;
            for(j = 1;j <= ny;j++)     //找出使Y集合不在增广路中的点能满足要求的最大可松弛量d
            {
                if(!visy[j]&&slack[j] < d)
                {
                    d = slack[j];
                }
            }
            for(j = 1;j <= nx;j++)     //X集合中在增广路中的点减去d
            {
                if(visx[j])
                {
                    lx[j] -= d;
                }
            }
            for(j = 1;j <= ny;j++)
            {
                if(visy[j])     //Y中集合在增广路中点加上d
                {
                    ly[j] += d;
                }
                else         //X中集合在增广路的点减去d,则Y中不在增广路中点松弛量还可减去d
                {
                    slack[j] -= d;
                }
            }
        }
    }
    for(i = 1;i <= ny;i++)
    {
        if(match[i])        //可能供大于求,Y中的点不一定有匹配
        {
            minCost += cost[match[i]][i];
        }
    }
    return 0;
}

int main()
{
    int N,M,K;
    int i,j,k;
    int tmpCost[maxN][maxN];
    int need[maxK][maxN];
    int storage[maxK][maxN];
    short node2shop[maxM],node2supply[maxM];
    short sumNeed[maxK];
    short sumStorage[maxK];
    char needEnough;
    while(~scanf("%d%d%d",&N,&M,&K)&&N)
    {
        minCost = 0;
        needEnough = true;
        memset(sumNeed,0,sizeof(sumNeed));
        memset(sumStorage,0,sizeof(sumStorage));
        for(i = 1;i <= N;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&need[k][i]);   //第k个商家需要多少第j种商品
                sumNeed[k] += need[k][i];
            }
        }
        for(i = 1;i <= M;i++)
        {
            for(k = 1;k <= K;k++)
            {
                scanf("%d",&storage[k][i]);  //第k个供货商有多少j的库存
                sumStorage[k] += storage[k][i];
            }
        }
        for(k = 1;k <= K;k++)
        {
            if(sumNeed[k] > sumStorage[k])
            {
                needEnough = false;
                break;
            }
        }
        for(k = 1;k <= K;k++)
        {
            for(i = 1;i <= N;i++)
            {
                for(j = 1;j <= M;j++)
                {
                    scanf("%d",&tmpCost[i][j]);
                }
            }
            if(needEnough)
            {
                nx = ny = 1;
                for(i = 1;i <= N;i++)
                {
                    for(j = 1;j <= need[k][i];j++)      //按照需求量将商店分解成若干需求为1的商店
                    {
                        node2shop[nx++] = i;        //记录新生成的商店原来属于哪个商店
                    }
                }
                for(i = 1;i <= M;i++)
                {
                    for(j = 1;j <= storage[k][i];j++)   //按照供应量将供应商分解成若干供应量为1的供应商
                    {
                        node2supply[ny++] = i;   //记录新生成的供应商原来属于哪个供应商
                    }
                }
                nx--;
                ny--;
                for(i = 1;i <= nx;i++)
                {
                    for(j = 1;j <= ny;j++)
                    {
                        cost[i][j] = -1*tmpCost[node2shop[i]][node2supply[j]];  //转换成相反数求最大二分加权匹配
                    }
                }
                KM();
            }
        }
        if(!needEnough)
        {
            printf("-1\n");
        }
        else
        {
            printf("%d\n",-1*minCost);
        }
    }
    return 0;
}


  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值