匈牙利算法求解指派问题（C++代码）

前言

匈牙利算法能精确求解指派问题，获取最优分配方案。匈牙利算法求解指派问题基于以下原理：在一个成本矩阵中，对某一行或者某一列加上或减去一个数，最优的分配方案不变。基于此原理，我们可以对成本矩阵进行变换，直到使用试指派能够找到最优解（对一个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文件运行结果如下所示：