前言
匈牙利算法能精确求解指派问题,获取最优分配方案。匈牙利算法求解指派问题基于以下原理:在一个成本矩阵中,对某一行或者某一列加上或减去一个数,最优的分配方案不变。基于此原理,我们可以对成本矩阵进行变换,直到使用试指派能够找到最优解(对一个n*n的成本矩阵而言,找到n个独立0元素)。
一、指派问题
实际中,我们会经常碰到此类问题:有n项任务需要均分给n个工人完成,工人i完成任务j的成本为cij,我们要找到一种分配方案,使得总成本最小。
如下图是一个4员工4任务的指派问题:
二、匈牙利算法
对于一个规模为n的指派问题而言,使用匈牙利算法求解的核心是找出n个独立0元素。匈牙利算法的大致流程如下图所示:
以前面介绍指派问题给出的数据为例,对匈牙利算法的流程做一个简单介绍,首先给出初始成本矩阵:
对矩阵进行行归约,每行减去该行的最小值,保证每行都出现0元素:
对矩阵进行列归约,每列减去该列的最小值,保证每列都出现0元素:
接下来是试指派,检验是否能找到n个独立0元素,否则需要对矩阵进行进一步变换以增加矩阵中的0元素个数。试指派从0元素最少的行或列开始,当确定了cij=0为所选元素,划去第i行和第j列,即不在考虑该区域内的元素。从0元素最少的行开始指派,因此指派的顺序是(行):2->3->1->4,通过试指派找到了最优分配,即员工1执行任务4、员工2执行任务2、员工3执行任务1、员工4执行任务3。
并不是所有的情况都能通过试指派找到最优分配方案(找到n个独立0元素),此时需要对矩阵进行变换,目的是增加0矩阵中0元素。假设最后得到的矩阵如下所示:
矩阵变换的目的是使得成本矩阵中出现更多的0元素,我们执行以下操作:
一:打勾操作
(1)对没有独立0元素的行打勾
(2)对打勾行含0元素的列打勾
(3)对打勾列含独立0元素的行打勾
(4)重复(2)(3)直到没有新的勾出现
此步值得注意的是要分清0元素和独立0元素,标红的为本次选择的独立0元素。
二:画线操作
(1)对没有打勾的行和打勾的列进行画线操作
三:矩阵变换
(1)找到未被线覆盖的最小元素
(2)没有画线的行减去该元素
(3)画线的列加上该元素
变换后的矩阵如下所示,对其进行试指派,可以找到n(此处为4)个独立0元素。
因此最后得到的分配方案是:员工1执行任务1、员工2执行任务2、员工3执行任务4、员工4执行任务3。
三 、C++代码
为了方便读者参考和使用,本文将匈牙利代码写在一个.cpp文件中,并将员工数量n和成本矩阵c的值存储在.dat文件中,对于不同数据,读者根据具体问题更改.dat文档中的数据即可。
数据文件
数据文件为一个后缀名为.dat的文档,格式如下图所示:
其中4为n,其后数据为一个4*4的成本矩阵。
代码
#include<iostream>
#include<algorithm>
#include<fstream>
using namespace std;
int n;//元素个数
int *assign;//分配结果
int **mat;//代价矩阵
int **matRcd;//代价矩阵
int totalCost;//总成本
//-------------------------------------数据读取
bool read(const char* ad)
{
//数据读取
ifstream iff;
iff.open(ad);
if(!iff)return false;
iff>>n;
assign=new int[n];for(int i=0;i<n;i++)assign[i]=-1;
mat=new int*[n];for(int i=0;i<n;i++)mat[i]=new int[n];
matRcd=new int*[n];for(int i=0;i<n;i++)matRcd[i]=new int[n];
for(int i=0;i<n;i++){
for(int j=0;j<n;j++){
iff>>mat[i][j];
matRcd[i][j]=mat[i][j];
}
}
iff.close();
totalCost=0;
return true;
}
//-------------------------------------行归约
void rowSub()
{
int *minEmt=new int[n];for(int i=0;i<n;i++)minEmt[i]=(int)1e8;
for(int i=0;i<n;i++)for(int j=0;j<n;j++)if(mat[i][j]<minEmt[i])minEmt[i]=mat[i][j];
for(int i=0;i<n;i++)for(int j=0;j<n;j++)mat[i][j]-=minEmt[i];
delete []minEmt;
}
//-------------------------------------列归约
void columnSub()
{
int *minEmt=new int[n];for(int j=0;j<n;j++)minEmt[j]=(int)1e8;
for(int j=0;j<n;j++)for(int i=0;i<n;i++)if(mat[i][j]<minEmt[j])minEmt[j]=mat[i][j];
for(int j=0;j<n;j++)for(int i=0;i<n;i++)mat[i][j]-=minEmt[j];
delete []minEmt;
}
//-------------------------------------检验最优
bool isOptimal(int *assign)
{
int *tAssign=new int[n];for(int i=0;i<n;i++)tAssign[i]=-1;
int *nZero=new int[n];
bool *rowIsUsed=new bool[n];
bool *columnIsUsed=new bool[n];
for(int i=0;i<n;i++)rowIsUsed[i]=columnIsUsed[i]=0;
int nLine=0;
while(nLine<n){
for(int i=0;i<n;i++)nZero[i]=0;
for(int i=0;i<n;i++){
if(rowIsUsed[i]==1)continue;
for(int j=0;j<n;j++){
if(columnIsUsed[j]!=1&&mat[i][j]==0)nZero[i]++;
}
}
int minZeros=n;
int rowId=-1;
for(int i=0;i<n;i++){
if(rowIsUsed[i]==0&&nZero[i]<minZeros&&nZero[i]>0){
minZeros=nZero[i];
rowId=i;
}
}
if(rowId==-1)break;
for(int j=0;j<n;j++){
if(mat[rowId][j]==0&&columnIsUsed[j]==0){
rowIsUsed[rowId]=1;
columnIsUsed[j]=1;
tAssign[rowId]=j;
break;
}
}
nLine++;
}
for(int i=0;i<n;i++)assign[i]=tAssign[i];
delete []tAssign;
delete []nZero;
delete []rowIsUsed;
delete []columnIsUsed;
for(int i=0;i<n;i++)if(assign[i]==-1)return false;
return true;
}
//-------------------------------------矩阵变换
void matTrans()
{
bool *rowTip=new bool[n];
bool *columnTip=new bool[n];
bool *rowLine=new bool[n];
bool *columnLine=new bool[n];
for(int i=0;i<n;i++)rowTip[i]=columnTip[i]=rowLine[i]=columnLine[i]=0;
//打勾
for(int i=0;i<n;i++)if(assign[i]==-1)rowTip[i]=1;
while(1){
int preTip=0;
for(int i=0;i<n;i++)preTip=preTip+rowTip[i]+columnTip[i];
for(int i=0;i<n;i++){
if(rowTip[i]==1){
for(int j=0;j<n;j++){
if(mat[i][j]==0)columnTip[j]=1;
}
}
}
for(int j=0;j<n;j++){
if(columnTip[j]==1){
for(int i=0;i<n;i++){
if(assign[i]==j)rowTip[i]=1;
}
}
}
int curTip=0;
for(int i=0;i<n;i++)curTip=curTip+rowTip[i]+columnTip[i];
if(preTip==curTip)break;
}
//画线
for(int i=0;i<n;i++){
if(rowTip[i]==0)rowLine[i]=1;
if(columnTip[i]==1)columnLine[i]=1;
}
//找最小元素
int minElmt=(int)1e8;
for(int i=0;i<n;i++)for(int j=0;j<n;j++)if(rowLine[i]==0&&columnLine[j]==0&&mat[i][j]<minElmt)minElmt=mat[i][j];
//变换
for(int i=0;i<n;i++)if(rowTip[i]==1)for(int j=0;j<n;j++)mat[i][j]-=minElmt;
for(int j=0;j<n;j++)if(columnTip[j]==1)for(int i=0;i<n;i++)mat[i][j]+=minElmt;
delete []rowTip;
delete []columnTip;
delete []rowLine;
delete []columnLine;
}
//-------------------------------------匈牙利算法
void hungary()
{
read("cost.dat");//读取数据
rowSub();//行归约
columnSub();//列归约
//如果不能找到n个独立的0元素,则对矩阵进行变换
while(!isOptimal(assign)){
matTrans();
}
for(int i=0;i<n;i++)totalCost+=matRcd[i][assign[i]];
for(int i=0;i<n;i++)delete []mat[i];delete []mat;
for(int i=0;i<n;i++)delete []matRcd[i];delete []matRcd;
}
int main()
{
//调用匈牙利算法
hungary();
//输出结果
cout<<"总成本为"<<totalCost<<endl;
for(int i=0;i<n;i++)cout<<"员工"<<i+1<<"-->"<<"任务"<<assign[i]+1<<endl;
cin.get();
}
代码中mat为成本矩阵,n为员工或任务数,assign数组为最后的分配方案,totalCost为此分配下的总成本。在主函数中调用匈牙算法函数hungary(),最后输出总成本以及分配方案。接口和结果都清晰明了,且代码按照上文叙述编写,结合注释,应该可以完全看懂。代码编译环境为VS2019,上述.dat文件运行结果如下所示: