匈牙利算法是为了解决二分图的最大匹配问题。其算法核心，其实就是不停的寻找增广路径，然后进行替换。这里引用Matrix67对其的理解，我觉得很精辟：
说穿了，就是你从二分图中找出一条路径来，让路径的起点和终点都是还没有匹配过的点，并且路径经过的连线是一条没被匹配、一条已经匹配过，再下一条又没匹配这样交替地出现。找到这样的路径后，显然路径里没被匹配的连线比已经匹配了的连线多一条，于是修改匹配图，把路径里所有匹配过的连线去掉匹配关系，把没有匹配的连线变成匹配的，这样匹配数就比原来多1个。不断执行上述操作，直到找不到这样的路径为止。
至于其他的具体概念，大家请自行谷歌。这里我给出匈牙利算法的Matlab和C++代码。网上虽然也有，但是基本都是属于不能用的那种，或者就是出来的结果是错的。
Matlab版：
[sourcecode language=”matlab”]
function assign=assignment(A)
[m,n] = size(A);
M(m,n)=0;
for(i=1:m)
for(j=1:n)
if(A(i,j))
M(i,j)=1;
break;
end
end %求初始匹配 M
if(M(i,j))
break;
end
end %获得仅含一条边的初始匹配 M
while(1)
for(i=1:m)
x(i)=0;
end %将记录X 中点的标号和标记*
for(i=1:n)
y(i)=0;
end %将记录Y 中点的标号和标记*
for(i=1:m)
pd=1; %寻找X 中 M 的所有非饱和点
for(j=1:n)
if(M(i,j))
pd=0;
end;
end
if(pd)
x(i)=-n-1;
end
end %将X 中 M 的所有非饱和点都给以标号0 和标记*, 程序中用 n+1 表示0 标号, 标号为负数时表示标记*
pd=0;
while(1)
xi=0;
for(i=1:m)
if(x(i)<0)
xi=i;
break;
end
end %假如 X 中存在一个既有标号又有标记*的点, 则任 取X 中一个既有标号又有标记*的点xi
if(xi==0)
pd=1;
break;
end %假如X 中所有有标号的点都已去掉了标记*, 算法终止
x(xi)=x(xi)*(-1); %去掉xi 的标记*
k=1;
for(j=1:n )
if(A(xi,j)&y(j)==0)
y(j)=xi;
yy(k)=j;
k=k+1;
end
end %对与 xi 邻接且尚未给标号的 yj 都给以标号i
if(k>1)
k=k-1;
for(j=1:k)
pdd=1;
for(i=1:m)
if(M(i,yy(j)))
x(i)=-yy(j);
pdd=0;
break;
end
end %将yj 在 M 中与之邻接的 点xk (即xkyj ∈M), 给以标号j 和标记*
if(pdd)
break;
end
end
if(pdd)
k=1;
j=yy(j); %yj 不是 M 的饱和点
while(1)
P(k,2)=j;
P(k,1)=y(j);
j=abs(x(y(j))); %任取 M 的一个非饱和点 yj, 逆向返回
if(j==n+1)
break;
end %找到X 中标号为0 的点时结束, 获得 M-增广路
k=k+1;
end
for(i=1:k)
if(M(P(i,1),P(i,2)))
M(P(i,1),P(i,2))=0; %将匹配 M 在增广路 P 中出现的边 去掉
else M(P(i,1),P(i,2))=1;
end
end %将增广路 P 中没有在匹配 M 中出现的边加入 到匹配 M 中
break;
end
end
end
if(pd)
break;
end
end %假如X 中所有有标号的点都已去掉了标记*, 算法终止
assign = M ; %显示最大匹配 M, 程序结束
[/sourcecode]
此版本中矩阵A中不连接的边对应的值请设置成0，下面再给出一个老外写的Matlab版，使用时不连接的边设置成Inf，如下（更多请点击此处）：
[sourcecode language=”matlab”]
function [assignment, cost] = assignmentoptimal(distMatrix)
%ASSIGNMENTOPTIMAL Compute optimal assignment by Munkres algorithm
% ASSIGNMENTOPTIMAL(DISTMATRIX) computes the optimal assignment (minimum
% overall costs) for the given rectangular distance or cost matrix, for
% example the assignment of tracks (in rows) to observations (in
% columns). The result is a column vector containing the assigned column
% number in each row (or 0 if no assignment could be done).
%
% [ASSIGNMENT, COST] = ASSIGNMENTOPTIMAL(DISTMATRIX) returns the
% assignment vector and the overall cost.
%
% The distance matrix may contain infinite values (forbidden
% assignments). Internally, the infinite values are set to a very large
% finite number, so that the Munkres algorithm itself works on
% finite-number matrices. Before returning the assignment, all
% assignments with infinite distance are deleted (i.e. set to zero).
%
% A description of Munkres algorithm (also called Hungarian algorithm)
% can easily be found on the web.
%
% <a href="assignment.html">assignment.html</a> <a href="http://www.mathworks.com/matlabcentral/fileexchange/6543">File Exchange</a> <a href="https://www.paypal.com/cgi-bin/webscr?cmd=_s-xclick&hosted_button_id=EVW2A4G2HBVAU">Donate via PayPal</a>
%
% Markus Buehren
% Last modified 05.07.2011
% save original distMatrix for cost computation
originalDistMatrix = distMatrix;
% check for negative elements
if any(distMatrix(:) < 0)
error(‘All matrix elements have to be non-negative.’);
end
% get matrix dimensions
[nOfRows, nOfColumns] = size(distMatrix);
% check for infinite values
finiteIndex = isfinite(distMatrix);
infiniteIndex = find(~finiteIndex);
if ~isempty(infiniteIndex)
% set infinite values to large finite value
maxFiniteValue = max(max(distMatrix(finiteIndex)));
if maxFiniteValue > 0
infValue = abs(10 * maxFiniteValue * nOfRows * nOfColumns);
else
infValue = 10;
end
if isempty(infValue)
% all elements are infinite
assignment = zeros(nOfRows, 1);
cost = 0;
return
end
distMatrix(infiniteIndex) = infValue;
end
% memory allocation
coveredColumns = zeros(1, nOfColumns);
coveredRows = zeros(nOfRows, 1);
starMatrix = zeros(nOfRows, nOfColumns);
primeMatrix = zeros(nOfRows, nOfColumns);
% preliminary steps
if nOfRows <= nOfColumns
minDim = nOfRows;
% find the smallest element of each row
minVector = min(distMatrix, [], 2);
% subtract the smallest element of each row from the row
distMatrix = distMatrix – repmat(minVector, 1, nOfColumns);
% Steps 1 and 2
for row = 1:nOfRows
for col = find(distMatrix(row,:)==0)
if ~coveredColumns(col)%~any(starMatrix(:,col))
starMatrix(row, col) = 1;
coveredColumns(col) = 1;
break
end
end
end
else % nOfRows > nOfColumns
minDim = nOfColumns;
% find the smallest element of each column
minVector = min(distMatrix);
% subtract the smallest element of each column from the column
distMatrix = distMatrix – repmat(minVector, nOfRows, 1);
% Steps 1 and 2
for col = 1:nOfColumns
for row = find(distMatrix(:,col)==0)’
if ~coveredRows(row)
starMatrix(row, col) = 1;
coveredColumns(col) = 1;
coveredRows(row) = 1;
break
end
end
end
coveredRows(:) = 0; % was used auxiliary above
end
if sum(coveredColumns) == minDim
% algorithm finished
assignment = buildassignmentvector__(starMatrix);
else
% move to step 3
[assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step3__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim); %#ok
end
% compute cost and remove invalid assignments
[assignment, cost] = computeassignmentcost__(assignment, originalDistMatrix, nOfRows);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function assignment = buildassignmentvector__(starMatrix)
[maxValue, assignment] = max(starMatrix, [], 2);
assignment(maxValue == 0) = 0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [assignment, cost] = computeassignmentcost__(assignment, distMatrix, nOfRows)
rowIndex = find(assignment);
costVector = distMatrix(rowIndex + nOfRows * (assignment(rowIndex)-1));
finiteIndex = isfinite(costVector);
cost = sum(costVector(finiteIndex));
assignment(rowIndex(~finiteIndex)) = 0;
% Step 2: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step2__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim)
% cover every column containing a starred zero
maxValue = max(starMatrix);
coveredColumns(maxValue == 1) = 1;
if sum(coveredColumns) == minDim
% algorithm finished
assignment = buildassignmentvector__(starMatrix);
else
% move to step 3
[assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step3__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim);
end
% Step 3: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step3__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim)
zerosFound = 1;
while zerosFound
zerosFound = 0;
for col = find(~coveredColumns)
for row = find(~coveredRows’)
if distMatrix(row,col) == 0
primeMatrix(row, col) = 1;
starCol = find(starMatrix(row,:));
if isempty(starCol)
% move to step 4
[assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step4__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, row, col, minDim);
return
else
coveredRows(row) = 1;
coveredColumns(starCol) = 0;
zerosFound = 1;
break % go on in next column
end
end
end
end
end
% move to step 5
[assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step5__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim);
% Step 4: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step4__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, row, col, minDim)
newStarMatrix = starMatrix;
newStarMatrix(row,col) = 1;
starCol = col;
starRow = find(starMatrix(:, starCol));
while ~isempty(starRow)
% unstar the starred zero
newStarMatrix(starRow, starCol) = 0;
% find primed zero in row
primeRow = starRow;
primeCol = find(primeMatrix(primeRow, :));
% star the primed zero
newStarMatrix(primeRow, primeCol) = 1;
% find starred zero in column
starCol = primeCol;
starRow = find(starMatrix(:, starCol));
end
starMatrix = newStarMatrix;
primeMatrix(:) = 0;
coveredRows(:) = 0;
% move to step 2
[assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step2__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim);
% Step 5: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step5__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim)
% find smallest uncovered element
uncoveredRowsIndex = find(~coveredRows’);
uncoveredColumnsIndex = find(~coveredColumns);
[s, index1] = min(distMatrix(uncoveredRowsIndex,uncoveredColumnsIndex));
[s, index2] = min(s); %#ok
h = distMatrix(uncoveredRowsIndex(index1(index2)), uncoveredColumnsIndex(index2));
% add h to each covered row
index = find(coveredRows);
distMatrix(index, 
= distMatrix(index, + h;
% subtract h from each uncovered column
distMatrix(:, uncoveredColumnsIndex) = distMatrix(:, uncoveredColumnsIndex) – h;
% move to step 3
[assignment, distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows] = step3__(distMatrix, starMatrix, primeMatrix, coveredColumns, coveredRows, minDim);
[/sourcecode]
这两个版本都是我自己运行过绝对可以出结果的。
下面给出C++版（完整源码请看这儿）：
[sourcecode language=”cpp”]
#include "Hungarian.h"
extern int _Verbose;
extern bool _Plot;
void
Hungarian::HungarianAlgo(BipartiteGraph& _bg, vector<VID>& _S, vector<VID>& _T,
vector<VID>& _N, vector<EID>& _EG, vector<EID>& _M){
PlotGraph plot;
if(_Plot){
//remove all previous temp files in /tmp
if(system("rm /tmp/gnuplot-i*")) cout<<endl;
plot.InitPlot();
}
if(_Verbose){
cout<<"The initial labels for vertices:"<<endl;
_bg.DisplayLabels();
}
while(!IsPerfect(_bg)){
VID u = PickFreeAgent(_bg);
if(_Verbose)
cout<<endl<<"Pick an unmatched agent: "<<u<<endl;
InitNewRoot(_bg, u, _S, _T);
while(1){
RefreshBG(_bg, _S, _T, _N, _EG, _M);
if(_Plot){
plot.PlotBipartiteGraph(_bg, _S, _T, _N, _EG, _M);
if(!plot.GetPeriod()){ cout<<"Please press Enter to continue: "<<endl; cin.get(); }
}
//if need relabel, update labels
if(NeedReLabel(_T, _N)){
if(_Verbose)
cout<<"Need relabel. Searching for min_delta…"<<endl;
UpdateLabels(_bg);
if(_Verbose > 1){
cout<<"After updating labels: "<<endl;
_bg.DisplayLabels();
}
}
RefreshBG(_bg, _S, _T, _N, _EG, _M);
if(_Verbose > 1)
DisplayData(1,1,1,1,1);
if(_Plot){
plot.PlotBipartiteGraph(_bg, _S, _T, _N, _EG, _M);
if(!plot.GetPeriod()){ cout<<"Please press Enter to continue: "<<endl; cin.get(); }
}
//if not need relabel
if(!NeedReLabel(_T, _N)){
VID y = PickAvailableTask(_T, _N);
if(_Verbose){
cout<<"No need relabel, pick available task(s) "<<y<<" and put in queue."<<endl;
cout<<"Got task VID "<<y<<endl;
if(!bg.GetTask(y)->GetMatched())
cout<<"Task VID "<<y<<" has NOT been matched, then "<<u<<"–>"<<y
<<" is an augmenting path, searching…"<<endl;
else
cout<<"Task VID "<<y<<" has been matched, and I gonna extend alternating tree."<<endl;
}
//if NOT need relabel and task NOT matched, then found augmenting path
if(!bg.GetTask(y)->GetMatched()){
vector<EID> _path = BFSAugmentingPath(_bg, u, y);
AugmentPath(_bg, _path);
RefreshBG(_bg, _S, _T, _N, _EG, _M);
if(_Verbose > 1){
cout<<"After augmenting path: "<<endl;
DisplayData(1,1,1,1,1);
}
if(_Plot){
plot.PlotAugmentingPath(_bg, _path);
plot.PlotBipartiteGraph(_bg, _S, _T, _N, _EG, _M, y);
if(!plot.GetPeriod()){ cout<<"Please press Enter to continue: "<<endl; cin.get(); }
}
break; //break innter while
}
//if not need relabel and task is *matched*, extend tree
else{
ExtendTree(_bg, y, _S, _T);
if(_Verbose)
cout<<"Extending alternating tree…"<<endl;
if(_Verbose>1){
cout<<"After extending tree: "<<endl;
DisplayData(1,1,1,1,1);
}
}
}//if(!NeedReLabel)
} //while(1)
}//while(perfect)
cout<<endl<<"*****************************************************"<<endl;
cout<<endl<<"Your queried assignment problem is:"<<endl;
DisplayMatrix(*(_bg.GetMatrix()));
cout<<"The matching configuration is:"<<endl;
DisplayConfig(_bg);
cout<<"@_@ Kuhn-Munkres Hungarian Perfect Matching:"<<endl;
DisplayData(_M);
cout<<"Optimization result (weights sum): "<<this->OptimalValue(_bg, _M);
cout<<endl;
cout<<endl<<"*****************************************************"<<endl;
cout<<endl;
if(_Plot){
cout << "Please press Enter to terminate: "<<endl<<endl;
plot.SetPeriod(POS_INF);
cin.get();
plot.ClosePlot();
}
}
void
Hungarian::HungarianAlgo(void){
HungarianAlgo(bg, S, T, N, EG, M);
//DisplayData(0, 0, 0, 1, 1);
}
void
Hungarian::InitNewRoot(BipartiteGraph& _bg, VID root, vector<VID>& _S, vector<VID>& _T){
//init the sets of S and T
_S.clear();
_T.clear();
_S.push_back(root);
//clear history
for(size_t i=0; i<_bg.GetNumAgents(); i++)
_bg.GetAgent(i)->SetColored(false);
for(size_t j=0; j<_bg.GetNumTasks(); j++)
_bg.GetTask(j)->SetColored(false);
//color the root
_bg.GetAgent(root)->SetColored(true);
}
void
Hungarian::RefreshBG(BipartiteGraph& _bg,
const vector<VID>& _S, const vector<VID>& _T,
vector<VID>& _N, vector<EID>& _EG, vector<EID>& _M){
//check feasibility, if not, warning
if(!_bg.CheckFeasibility())
cout<<"Warning: Bipartite Graph is not feasible!"<<endl;
//check admissible edges, update admissible flags for edges, and update set EG
_EG.clear();
for(unsigned int i=0; i<_bg.GetNumAgents(); i++)
for(unsigned int j=0; j<_bg.GetNumTasks(); j++){
if( fabs(_bg.GetMatrix(i,j)->GetWeight() – (_bg.GetAgent(i)->GetLabel() + _bg.GetTask(j)->GetLabel())) <= DOUBLE_EPSILON ){
_bg.GetMatrix(i,j)->SetAdmissibleFlag(true);
_EG.push_back(EID(i,j));
}
else
_bg.GetMatrix(i,j)->SetAdmissibleFlag(false);
} //for j
//check matched edges, update matched flags for vertices, and update set M
_M.clear();
size_t cnt = 0;
for(unsigned int i=0; i<_bg.GetNumAgents(); i++)
_bg.GetAgent(i)->SetMatched(false);
for(unsigned int j=0; j<_bg.GetNumTasks(); j++)
_bg.GetTask(j)->SetMatched(false);
for(unsigned int i=0; i<_bg.GetNumAgents(); i++)
for(unsigned int j=0; j<_bg.GetNumTasks(); j++)
if(_bg.GetMatrix(i,j)->GetMatchedFlag()){
_bg.GetAgent(i)->SetMatched(true);
_bg.GetTask(j)->SetMatched(true);
_M.push_back(EID(i, j));
cnt++;
}//if
//update the variable in BG
_bg.SetNumMatched(cnt);
//update set N
_N.clear();
_N = this->FindNeighbors(_EG, _S);
}
void
Hungarian::ExtendTree(BipartiteGraph& _bg, VID x, vector<VID>& _S, vector<VID>& _T){
//if task x matched some agent
int vid_agent = -1;
if(_bg.GetTask(x)->GetMatched()){
for(unsigned int i=0; i<_bg.GetNumAgents(); i++)
if(_bg.GetMatrix(i, x)->GetMatchedFlag()){
vid_agent = i;
break;
}
if(-1 == vid_agent){
cerr<<"Error: Matched task vertex could not find its matching agent. Stopped."<<endl;
exit(-1);
}
_S.push_back(vid_agent);
_bg.GetAgent(vid_agent)->SetColored(true);
_T.push_back(x);
_bg.GetTask(x)->SetColored(true);
//EG.push_back(EID(vid_agent,x));
}//if
}
//AugmentPath should be followed by RefreshBG(), since after this, config changes.
void
Hungarian::AugmentPath(BipartiteGraph& _bg, vector<EID>& _path){
for(vector<EID>::iterator itr = _path.begin(); itr != _path.end(); itr++){
if(_bg.GetMatrix(itr->first, itr->second)->GetMatchedFlag())
_bg.GetMatrix(itr->first, itr->second)->SetMatchedFlag(false);
else
_bg.GetMatrix(itr->first, itr->second)->SetMatchedFlag(true);
}
//Below has been moved to RefreshBG()
/*
//update matched variable in bg.
int cnt = 0;
M.clear();
for(unsigned int i=0; i<_bg.GetNumAgents(); i++)
for(unsigned int j=0; j<_bg.GetNumTasks(); j++)
if(_bg.bg_matrix(i, j).GetMatchedFlag()){
cnt++;
M.push_back(EID(i,j));
}
_bg.SetNumMatched((size_t)cnt);
*/
}
vector<EID>
Hungarian::BFSAugmentingPath(BipartiteGraph& _bg, VID x, VID y){
size_t tasks_size = _bg.GetNumTasks();
size_t agents_size = _bg.GetNumAgents();
bool found = false;
vector<EID> aug_path;
deque<Vertex> dq;
//before alternating tree, clear all paths
for(unsigned int i=0; i<tasks_size; i++)
_bg.GetTask(i)->path.clear();
for(unsigned int j=0; j<agents_size; j++)
_bg.GetAgent(j)->path.clear();
dq.push_back(*_bg.GetAgent(x));
//cout<<"dq.front: "<<dq.front().GetVID()<<endl;
//loop for searching a path, using queue.
while(dq.size() && !found){
Vertex v = dq.front();
dq.pop_front();
if(v.GetObj() == "TASK"){
if(_Verbose)
cout<<"I am "<<v.GetObj()<<" VID: "<<v.GetVID()<<endl;
size_t agent_index = 0;
for(agent_index=0; agent_index<agents_size; agent_index++)
if(_bg.GetMatrix(agent_index, v.GetVID())->GetMatchedFlag() && !_bg.GetAgent(agent_index)->path.size()){
if(_Verbose)
cout<<"–Found matched AGENT VID: "<<agent_index<<endl;
_bg.GetAgent(agent_index)->path = v.path;
_bg.GetAgent(agent_index)->path.push_back(EID(agent_index, v.GetVID()));
dq.push_back(*_bg.GetAgent(agent_index));
break; //only ONE possible matched edge for ONE vertex.
}
if(agent_index == agents_size){
if(_Verbose)
cout<<"Warning: I found no matched edge during back-tracking from task vertex: "<< v.GetVID()<<endl
<<" It’s possible there exists new other satisfied task vertices. Continue."<<endl;
//exit(-1);
}
}//if "TASK"
else if(v.GetObj() == "AGENT"){
if(_Verbose)
cout<<"I am "<<v.GetObj()<<" VID: "<<v.GetVID()<<endl;
size_t task_index = 0;
for(task_index=0; task_index<tasks_size; task_index++)
if(!_bg.GetMatrix(v.GetVID(), task_index)->GetMatchedFlag()
&& _bg.GetMatrix(v.GetVID(), task_index)->GetAdmissibleFlag()
&& !_bg.GetTask(task_index)->path.size()){
if(_Verbose){
if(task_index == y)
cout<<"–Found target TASK VID: "<<task_index<<", augmenting path found!"<<endl;
else
cout<<"–Found matched TASK VID: "<<task_index<<endl;
}
_bg.GetTask(task_index)->path = v.path;
_bg.GetTask(task_index)->path.push_back(EID(v.GetVID(), task_index));
if(task_index == y){
found = true;
aug_path = _bg.GetTask(y)->path;
if(_Verbose){
cout<<"Display augmenting path: "<<endl;
DisplayData(aug_path);
}
break;
}
else{
dq.push_back(*_bg.GetTask(task_index));
//break;
}
} //if(!_bg.bg_..)
}// if "AGENT"
else{
//wierd situation
cerr<<"Error: Some vertices are found not initialized with Obj…Stopped."<<endl;
exit(-1);
}
}//end while
return aug_path;
}
bool
Hungarian::NeedReLabel(vector<VID>& _T, vector<VID>& _N){
//Below has been done by RefreshBG();
//get the neighbor _N for _S
//_N.clear();
//_N = this->FindNeighbors(_EG, _S);
//then compare the elements in _N and _T
if(_N.size() != _T.size())
return false;
else{
sort(_N.begin(), _N.end());
sort(_T.begin(), _T.end());
if(_N == _T) return true;
else return false;
}
}
vector<VID>
Hungarian::FindNeighbors(const vector<EID>& _EG, const vector<VID>& _S){
vector<VID> _N;
_N.clear();
for(vector<EID>::const_iterator e_itr = _EG.begin(); e_itr != _EG.end(); e_itr++)
for(vector<VID>::const_iterator v1_itr=_S.begin(); v1_itr != _S.end(); v1_itr++)
if(e_itr->first == *v1_itr){
vector<VID>::iterator v2_itr;
for(v2_itr=_N.begin(); v2_itr!= _N.end(); v2_itr++)
if(e_itr->second == *v2_itr) break;
if(v2_itr==_N.end())
_N.push_back(e_itr->second);
} //if
return _N;
}
VID
Hungarian::PickFreeAgent(BipartiteGraph& _bg){
//if still not perfect
int free_agent = -1;
for(size_t i=0; i<_bg.GetNumAgents(); i++)
if(!_bg.GetAgent(i)->GetMatched()){
free_agent = i;
break;
}
if(-1 == free_agent){
cerr<<"Error: No free agent vertex available. Stopped."<<endl;
exit(-1);
}
return (VID)free_agent;
}
VID
Hungarian::PickAvailableTask(vector<VID>& _T, vector<VID>& _N){
int y = -1;
for(vector<VID>::iterator itr1 = _N.begin(); itr1 != _N.end(); itr1++){
vector<VID>::iterator itr2;
for(itr2 = _T.begin(); itr2 != _T.end(); itr2++)
if(*itr2 == *itr1) break;
if(itr2 == _T.end()){
y=*itr1;
break;
}
}
if(-1 == y){
cerr<<"Error: No free vertex can be picked! Stopped."<<endl;
exit(-1);
}
return (VID)y;
}
double
Hungarian::UpdateLabels(BipartiteGraph& _bg){
double delta;
double _min_delta = POS_INF;
//calculate the mininal delta for re-labelling
for(size_t i=0; i<_bg.GetNumAgents(); i++){
if(_bg.GetAgent(i)->GetColored()) //agent is in set of S
for(size_t j=0; j<_bg.GetNumTasks(); j++){
if(!_bg.GetTask(j)->GetColored()){ //task is NOT in set of T
delta = _bg.GetAgent(i)->GetLabel() + _bg.GetTask(j)->GetLabel() – _bg.GetMatrix(i, j)->GetWeight();
if(delta < _min_delta && fabs(delta – _min_delta) > DOUBLE_EPSILON){
_min_delta = delta;
if(_Verbose)
cout<<"("<<i<<","<<j<<")->min_data: "<<delta<<endl;
}
}//if
}//interior for
}//exterior for
//cout<<"min_delta: "<<_min_delta<<endl;
_bg.SetMinDelta(_min_delta);
//agents re-labelling
for(size_t i=0; i<_bg.GetNumAgents(); i++)
//agent is in set of S
if(_bg.GetAgent(i)->GetColored()){
double lb = _bg.GetAgent(i)->GetLabel() – _min_delta;
_bg.GetAgent(i)->SetLabel(lb);
}
//tasks re-labelling
for(size_t j=0; j<_bg.GetNumTasks(); j++)
//task is in set of T
if(_bg.GetTask(j)->GetColored()){
double lb = _bg.GetTask(j)->GetLabel() + _min_delta;
_bg.GetTask(j)->SetLabel(lb);
}
return _min_delta;
}
double
Hungarian::OptimalValue(BipartiteGraph& _bg, vector<EID>& _M){
double sum = 0;
for(vector<EID>::iterator itr=_M.begin(); itr!=_M.end(); itr++)
sum += _bg.GetMatrix(*itr)->GetWeight();
return sum;
}
void
Hungarian::DisplayData(vector<VID>& v){
for(vector<VID>::iterator itr = v.begin(); itr != v.end(); itr++)
cout<<*itr<<" ";
cout<<endl;
}
void
Hungarian::DisplayData(vector<EID>& e){
for(vector<EID>::iterator itr = e.begin(); itr != e.end(); itr++)
cout<<"("<<itr->first<<","<<itr->second<<") ";
cout<<endl;
}
void
Hungarian::DisplayData(bool s, bool t, bool n, bool eg, bool m){
//displaying S
if(s){
cout<<"Set S:"<<endl;
DisplayData(S);
}
//displaying T
if(t){
cout<<"Set T:"<<endl;
DisplayData(T);
}
//displaying N
if(n){
cout<<"Set N:"<<endl;
DisplayData(N);
}
//displaying EG
if(eg){
cout<<"Set EG:"<<endl;
DisplayData(EG);
}
//displaying M
if(m){
cout<<"Set M:"<<endl;
DisplayData(M);
}
}
void
Hungarian::DisplayLabels(vector<Vertex>& vertices){
for(vector<Vertex>::iterator itr = vertices.begin(); itr != vertices.end(); itr++)
cout<<itr->GetLabel()<<" ";
cout<<endl;
}
void
Hungarian::DisplayMatrix(Matrix& _m){
if(_m[0].size() > 30){
cout<<"Matrix is big, no displaying. "<<endl;
return;
}
for(unsigned int i=0; i<_m.size(); i++){
for(unsigned int j=0; j<_m[0].size(); j++)
cout<<" "<<_m[i][j].GetWeight()<<"t";
cout<<endl;
}
cout<<endl;
}
void
Hungarian::DisplayConfig(BipartiteGraph& _bg){
if(_bg.GetNumTasks() > 30){
cout<<"Matrix is big, no displaying. "<<endl;
return;
}
for(unsigned int i=0; i<_bg.GetNumAgents(); i++){
for(unsigned int j=0; j<_bg.GetNumTasks(); j++){
if(_bg.GetMatrix(i,j)->GetMatchedFlag())
cout<<" 1t";
else
cout<<" 0t";
}
cout<<endl;
}
cout<<endl;
}
[/sourcecode]
时间有限，未亲自验证，大家可以自行下载验证。
总的来说匈牙利算法还是很好理解的。