匈牙利算法 Matlab/C++代码

匈牙利算法是为了解决二分图的最大匹配问题。其算法核心,其实就是不停的寻找增广路径,然后进行替换。这里引用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&amp; _bg, vector&lt;VID&gt;&amp; _S, vector&lt;VID&gt;&amp; _T,
vector&lt;VID&gt;&amp; _N, vector&lt;EID&gt;&amp; _EG, vector&lt;EID&gt;&amp; _M){
PlotGraph plot;
if(_Plot){
//remove all previous temp files in /tmp
if(system("rm /tmp/gnuplot-i*")) cout&lt;&lt;endl;
plot.InitPlot();
}

if(_Verbose){
cout&lt;&lt;"The initial labels for vertices:"&lt;&lt;endl;
_bg.DisplayLabels();
}

while(!IsPerfect(_bg)){

VID u = PickFreeAgent(_bg);
if(_Verbose)
cout&lt;&lt;endl&lt;&lt;"Pick an unmatched agent: "&lt;&lt;u&lt;&lt;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&lt;&lt;"Please press Enter to continue: "&lt;&lt;endl; cin.get(); }
}

//if need relabel, update labels
if(NeedReLabel(_T, _N)){
if(_Verbose)
cout&lt;&lt;"Need relabel. Searching for min_delta…"&lt;&lt;endl;
UpdateLabels(_bg);
if(_Verbose &gt; 1){
cout&lt;&lt;"After updating labels: "&lt;&lt;endl;
_bg.DisplayLabels();
}
}
RefreshBG(_bg, _S, _T, _N, _EG, _M);
if(_Verbose &gt; 1)
DisplayData(1,1,1,1,1);
if(_Plot){
plot.PlotBipartiteGraph(_bg, _S, _T, _N, _EG, _M);
if(!plot.GetPeriod()){ cout&lt;&lt;"Please press Enter to continue: "&lt;&lt;endl; cin.get(); }
}

//if not need relabel
if(!NeedReLabel(_T, _N)){
VID y = PickAvailableTask(_T, _N);
if(_Verbose){
cout&lt;&lt;"No need relabel, pick available task(s) "&lt;&lt;y&lt;&lt;" and put in queue."&lt;&lt;endl;
cout&lt;&lt;"Got task VID "&lt;&lt;y&lt;&lt;endl;
if(!bg.GetTask(y)-&gt;GetMatched())
cout&lt;&lt;"Task VID "&lt;&lt;y&lt;&lt;" has NOT been matched, then "&lt;&lt;u&lt;&lt;"–&gt;"&lt;&lt;y
&lt;&lt;" is an augmenting path, searching…"&lt;&lt;endl;
else
cout&lt;&lt;"Task VID "&lt;&lt;y&lt;&lt;" has been matched, and I gonna extend alternating tree."&lt;&lt;endl;
}

//if NOT need relabel and task NOT matched, then found augmenting path
if(!bg.GetTask(y)-&gt;GetMatched()){
vector&lt;EID&gt; _path = BFSAugmentingPath(_bg, u, y);
AugmentPath(_bg, _path);
RefreshBG(_bg, _S, _T, _N, _EG, _M);
if(_Verbose &gt; 1){
cout&lt;&lt;"After augmenting path: "&lt;&lt;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&lt;&lt;"Please press Enter to continue: "&lt;&lt;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&lt;&lt;"Extending alternating tree…"&lt;&lt;endl;
if(_Verbose&gt;1){
cout&lt;&lt;"After extending tree: "&lt;&lt;endl;
DisplayData(1,1,1,1,1);
}
}
}//if(!NeedReLabel)
} //while(1)
}//while(perfect)

cout&lt;&lt;endl&lt;&lt;"*****************************************************"&lt;&lt;endl;
cout&lt;&lt;endl&lt;&lt;"Your queried assignment problem is:"&lt;&lt;endl;
DisplayMatrix(*(_bg.GetMatrix()));
cout&lt;&lt;"The matching configuration is:"&lt;&lt;endl;
DisplayConfig(_bg);

cout&lt;&lt;"@_@ Kuhn-Munkres Hungarian Perfect Matching:"&lt;&lt;endl;
DisplayData(_M);
cout&lt;&lt;"Optimization result (weights sum): "&lt;&lt;this-&gt;OptimalValue(_bg, _M);
cout&lt;&lt;endl;
cout&lt;&lt;endl&lt;&lt;"*****************************************************"&lt;&lt;endl;
cout&lt;&lt;endl;
if(_Plot){
cout &lt;&lt; "Please press Enter to terminate: "&lt;&lt;endl&lt;&lt;endl;
plot.SetPeriod(POS_INF);
cin.get();
plot.ClosePlot();
}

}

&nbsp;

void
Hungarian::HungarianAlgo(void){

HungarianAlgo(bg, S, T, N, EG, M);
//DisplayData(0, 0, 0, 1, 1);

}

&nbsp;

void
Hungarian::InitNewRoot(BipartiteGraph&amp; _bg, VID root, vector&lt;VID&gt;&amp; _S, vector&lt;VID&gt;&amp; _T){

//init the sets of S and T
_S.clear();
_T.clear();
_S.push_back(root);

//clear history
for(size_t i=0; i&lt;_bg.GetNumAgents(); i++)
_bg.GetAgent(i)-&gt;SetColored(false);
for(size_t j=0; j&lt;_bg.GetNumTasks(); j++)
_bg.GetTask(j)-&gt;SetColored(false);

//color the root
_bg.GetAgent(root)-&gt;SetColored(true);

}

&nbsp;

void
Hungarian::RefreshBG(BipartiteGraph&amp; _bg,
const vector&lt;VID&gt;&amp; _S, const vector&lt;VID&gt;&amp; _T,
vector&lt;VID&gt;&amp; _N, vector&lt;EID&gt;&amp; _EG, vector&lt;EID&gt;&amp; _M){

//check feasibility, if not, warning
if(!_bg.CheckFeasibility())
cout&lt;&lt;"Warning: Bipartite Graph is not feasible!"&lt;&lt;endl;

//check admissible edges, update admissible flags for edges, and update set EG
_EG.clear();
for(unsigned int i=0; i&lt;_bg.GetNumAgents(); i++)
for(unsigned int j=0; j&lt;_bg.GetNumTasks(); j++){
if( fabs(_bg.GetMatrix(i,j)-&gt;GetWeight() – (_bg.GetAgent(i)-&gt;GetLabel() + _bg.GetTask(j)-&gt;GetLabel())) &lt;= DOUBLE_EPSILON ){
_bg.GetMatrix(i,j)-&gt;SetAdmissibleFlag(true);
_EG.push_back(EID(i,j));
}
else
_bg.GetMatrix(i,j)-&gt;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&lt;_bg.GetNumAgents(); i++)
_bg.GetAgent(i)-&gt;SetMatched(false);
for(unsigned int j=0; j&lt;_bg.GetNumTasks(); j++)
_bg.GetTask(j)-&gt;SetMatched(false);
for(unsigned int i=0; i&lt;_bg.GetNumAgents(); i++)
for(unsigned int j=0; j&lt;_bg.GetNumTasks(); j++)
if(_bg.GetMatrix(i,j)-&gt;GetMatchedFlag()){
_bg.GetAgent(i)-&gt;SetMatched(true);
_bg.GetTask(j)-&gt;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-&gt;FindNeighbors(_EG, _S);

}
void
Hungarian::ExtendTree(BipartiteGraph&amp; _bg, VID x, vector&lt;VID&gt;&amp; _S, vector&lt;VID&gt;&amp; _T){
//if task x matched some agent
int vid_agent = -1;
if(_bg.GetTask(x)-&gt;GetMatched()){
for(unsigned int i=0; i&lt;_bg.GetNumAgents(); i++)
if(_bg.GetMatrix(i, x)-&gt;GetMatchedFlag()){
vid_agent = i;
break;
}

if(-1 == vid_agent){
cerr&lt;&lt;"Error: Matched task vertex could not find its matching agent. Stopped."&lt;&lt;endl;
exit(-1);
}

_S.push_back(vid_agent);
_bg.GetAgent(vid_agent)-&gt;SetColored(true);
_T.push_back(x);
_bg.GetTask(x)-&gt;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&amp; _bg, vector&lt;EID&gt;&amp; _path){
for(vector&lt;EID&gt;::iterator itr = _path.begin(); itr != _path.end(); itr++){
if(_bg.GetMatrix(itr-&gt;first, itr-&gt;second)-&gt;GetMatchedFlag())
_bg.GetMatrix(itr-&gt;first, itr-&gt;second)-&gt;SetMatchedFlag(false);
else
_bg.GetMatrix(itr-&gt;first, itr-&gt;second)-&gt;SetMatchedFlag(true);
}

//Below has been moved to RefreshBG()
/*
//update matched variable in bg.
int cnt = 0;
M.clear();
for(unsigned int i=0; i&lt;_bg.GetNumAgents(); i++)
for(unsigned int j=0; j&lt;_bg.GetNumTasks(); j++)
if(_bg.bg_matrix(i, j).GetMatchedFlag()){
cnt++;
M.push_back(EID(i,j));
}
_bg.SetNumMatched((size_t)cnt);
*/
}

&nbsp;

vector&lt;EID&gt;
Hungarian::BFSAugmentingPath(BipartiteGraph&amp; _bg, VID x, VID y){
size_t tasks_size = _bg.GetNumTasks();
size_t agents_size = _bg.GetNumAgents();
bool found = false;
vector&lt;EID&gt; aug_path;
deque&lt;Vertex&gt; dq;

//before alternating tree, clear all paths
for(unsigned int i=0; i&lt;tasks_size; i++)
_bg.GetTask(i)-&gt;path.clear();
for(unsigned int j=0; j&lt;agents_size; j++)
_bg.GetAgent(j)-&gt;path.clear();

dq.push_back(*_bg.GetAgent(x));
//cout&lt;&lt;"dq.front: "&lt;&lt;dq.front().GetVID()&lt;&lt;endl;

//loop for searching a path, using queue.
while(dq.size() &amp;&amp; !found){
Vertex v = dq.front();
dq.pop_front();
if(v.GetObj() == "TASK"){
if(_Verbose)
cout&lt;&lt;"I am "&lt;&lt;v.GetObj()&lt;&lt;" VID: "&lt;&lt;v.GetVID()&lt;&lt;endl;
size_t agent_index = 0;
for(agent_index=0; agent_index&lt;agents_size; agent_index++)
if(_bg.GetMatrix(agent_index, v.GetVID())-&gt;GetMatchedFlag() &amp;&amp; !_bg.GetAgent(agent_index)-&gt;path.size()){
if(_Verbose)
cout&lt;&lt;"–Found matched AGENT VID: "&lt;&lt;agent_index&lt;&lt;endl;
_bg.GetAgent(agent_index)-&gt;path = v.path;
_bg.GetAgent(agent_index)-&gt;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&lt;&lt;"Warning: I found no matched edge during back-tracking from task vertex: "&lt;&lt; v.GetVID()&lt;&lt;endl
&lt;&lt;" It’s possible there exists new other satisfied task vertices. Continue."&lt;&lt;endl;
//exit(-1);
}
}//if "TASK"

else if(v.GetObj() == "AGENT"){
if(_Verbose)
cout&lt;&lt;"I am "&lt;&lt;v.GetObj()&lt;&lt;" VID: "&lt;&lt;v.GetVID()&lt;&lt;endl;
size_t task_index = 0;
for(task_index=0; task_index&lt;tasks_size; task_index++)
if(!_bg.GetMatrix(v.GetVID(), task_index)-&gt;GetMatchedFlag()
&amp;&amp; _bg.GetMatrix(v.GetVID(), task_index)-&gt;GetAdmissibleFlag()
&amp;&amp; !_bg.GetTask(task_index)-&gt;path.size()){
if(_Verbose){
if(task_index == y)
cout&lt;&lt;"–Found target TASK VID: "&lt;&lt;task_index&lt;&lt;", augmenting path found!"&lt;&lt;endl;
else
cout&lt;&lt;"–Found matched TASK VID: "&lt;&lt;task_index&lt;&lt;endl;
}
_bg.GetTask(task_index)-&gt;path = v.path;
_bg.GetTask(task_index)-&gt;path.push_back(EID(v.GetVID(), task_index));
if(task_index == y){
found = true;
aug_path = _bg.GetTask(y)-&gt;path;
if(_Verbose){
cout&lt;&lt;"Display augmenting path: "&lt;&lt;endl;
DisplayData(aug_path);
}
break;
}
else{
dq.push_back(*_bg.GetTask(task_index));
//break;
}
} //if(!_bg.bg_..)
}// if "AGENT"

else{
//wierd situation
cerr&lt;&lt;"Error: Some vertices are found not initialized with Obj…Stopped."&lt;&lt;endl;
exit(-1);
}

}//end while

return aug_path;
}
bool
Hungarian::NeedReLabel(vector&lt;VID&gt;&amp; _T, vector&lt;VID&gt;&amp; _N){

//Below has been done by RefreshBG();
//get the neighbor _N for _S
//_N.clear();
//_N = this-&gt;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&lt;VID&gt;
Hungarian::FindNeighbors(const vector&lt;EID&gt;&amp; _EG, const vector&lt;VID&gt;&amp; _S){
vector&lt;VID&gt; _N;
_N.clear();
for(vector&lt;EID&gt;::const_iterator e_itr = _EG.begin(); e_itr != _EG.end(); e_itr++)
for(vector&lt;VID&gt;::const_iterator v1_itr=_S.begin(); v1_itr != _S.end(); v1_itr++)
if(e_itr-&gt;first == *v1_itr){
vector&lt;VID&gt;::iterator v2_itr;
for(v2_itr=_N.begin(); v2_itr!= _N.end(); v2_itr++)
if(e_itr-&gt;second == *v2_itr) break;
if(v2_itr==_N.end())
_N.push_back(e_itr-&gt;second);
} //if
return _N;
}
VID
Hungarian::PickFreeAgent(BipartiteGraph&amp; _bg){

//if still not perfect
int free_agent = -1;

for(size_t i=0; i&lt;_bg.GetNumAgents(); i++)
if(!_bg.GetAgent(i)-&gt;GetMatched()){
free_agent = i;
break;
}

if(-1 == free_agent){
cerr&lt;&lt;"Error: No free agent vertex available. Stopped."&lt;&lt;endl;
exit(-1);
}

return (VID)free_agent;
}
VID
Hungarian::PickAvailableTask(vector&lt;VID&gt;&amp; _T, vector&lt;VID&gt;&amp; _N){
int y = -1;
for(vector&lt;VID&gt;::iterator itr1 = _N.begin(); itr1 != _N.end(); itr1++){
vector&lt;VID&gt;::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&lt;&lt;"Error: No free vertex can be picked! Stopped."&lt;&lt;endl;
exit(-1);
}

return (VID)y;

}
double
Hungarian::UpdateLabels(BipartiteGraph&amp; _bg){
double delta;
double _min_delta = POS_INF;

//calculate the mininal delta for re-labelling
for(size_t i=0; i&lt;_bg.GetNumAgents(); i++){
if(_bg.GetAgent(i)-&gt;GetColored()) //agent is in set of S
for(size_t j=0; j&lt;_bg.GetNumTasks(); j++){
if(!_bg.GetTask(j)-&gt;GetColored()){ //task is NOT in set of T
delta = _bg.GetAgent(i)-&gt;GetLabel() + _bg.GetTask(j)-&gt;GetLabel() – _bg.GetMatrix(i, j)-&gt;GetWeight();
if(delta &lt; _min_delta &amp;&amp; fabs(delta – _min_delta) &gt; DOUBLE_EPSILON){
_min_delta = delta;
if(_Verbose)
cout&lt;&lt;"("&lt;&lt;i&lt;&lt;","&lt;&lt;j&lt;&lt;")-&gt;min_data: "&lt;&lt;delta&lt;&lt;endl;
}
}//if
}//interior for
}//exterior for

//cout&lt;&lt;"min_delta: "&lt;&lt;_min_delta&lt;&lt;endl;
_bg.SetMinDelta(_min_delta);

//agents re-labelling
for(size_t i=0; i&lt;_bg.GetNumAgents(); i++)
//agent is in set of S
if(_bg.GetAgent(i)-&gt;GetColored()){
double lb = _bg.GetAgent(i)-&gt;GetLabel() – _min_delta;
_bg.GetAgent(i)-&gt;SetLabel(lb);
}

//tasks re-labelling
for(size_t j=0; j&lt;_bg.GetNumTasks(); j++)
//task is in set of T
if(_bg.GetTask(j)-&gt;GetColored()){
double lb = _bg.GetTask(j)-&gt;GetLabel() + _min_delta;
_bg.GetTask(j)-&gt;SetLabel(lb);
}

return _min_delta;
}
double
Hungarian::OptimalValue(BipartiteGraph&amp; _bg, vector&lt;EID&gt;&amp; _M){

double sum = 0;
for(vector&lt;EID&gt;::iterator itr=_M.begin(); itr!=_M.end(); itr++)
sum += _bg.GetMatrix(*itr)-&gt;GetWeight();

return sum;

}
void
Hungarian::DisplayData(vector&lt;VID&gt;&amp; v){
for(vector&lt;VID&gt;::iterator itr = v.begin(); itr != v.end(); itr++)
cout&lt;&lt;*itr&lt;&lt;" ";
cout&lt;&lt;endl;
}

void
Hungarian::DisplayData(vector&lt;EID&gt;&amp; e){
for(vector&lt;EID&gt;::iterator itr = e.begin(); itr != e.end(); itr++)
cout&lt;&lt;"("&lt;&lt;itr-&gt;first&lt;&lt;","&lt;&lt;itr-&gt;second&lt;&lt;") ";
cout&lt;&lt;endl;
}

void
Hungarian::DisplayData(bool s, bool t, bool n, bool eg, bool m){

//displaying S
if(s){
cout&lt;&lt;"Set S:"&lt;&lt;endl;
DisplayData(S);
}

//displaying T
if(t){
cout&lt;&lt;"Set T:"&lt;&lt;endl;
DisplayData(T);
}

//displaying N
if(n){
cout&lt;&lt;"Set N:"&lt;&lt;endl;
DisplayData(N);
}

//displaying EG
if(eg){
cout&lt;&lt;"Set EG:"&lt;&lt;endl;
DisplayData(EG);
}

//displaying M
if(m){
cout&lt;&lt;"Set M:"&lt;&lt;endl;
DisplayData(M);
}

}
void
Hungarian::DisplayLabels(vector&lt;Vertex&gt;&amp; vertices){
for(vector&lt;Vertex&gt;::iterator itr = vertices.begin(); itr != vertices.end(); itr++)
cout&lt;&lt;itr-&gt;GetLabel()&lt;&lt;" ";
cout&lt;&lt;endl;
}
void
Hungarian::DisplayMatrix(Matrix&amp; _m){

if(_m[0].size() &gt; 30){
cout&lt;&lt;"Matrix is big, no displaying. "&lt;&lt;endl;
return;
}

for(unsigned int i=0; i&lt;_m.size(); i++){
for(unsigned int j=0; j&lt;_m[0].size(); j++)
cout&lt;&lt;" "&lt;&lt;_m[i][j].GetWeight()&lt;&lt;"t";
cout&lt;&lt;endl;
}
cout&lt;&lt;endl;

}
void
Hungarian::DisplayConfig(BipartiteGraph&amp; _bg){

if(_bg.GetNumTasks() &gt; 30){
cout&lt;&lt;"Matrix is big, no displaying. "&lt;&lt;endl;
return;
}

for(unsigned int i=0; i&lt;_bg.GetNumAgents(); i++){
for(unsigned int j=0; j&lt;_bg.GetNumTasks(); j++){
if(_bg.GetMatrix(i,j)-&gt;GetMatchedFlag())
cout&lt;&lt;" 1t";
else
cout&lt;&lt;" 0t";
}
cout&lt;&lt;endl;
}
cout&lt;&lt;endl;

}

[/sourcecode]

时间有限,未亲自验证,大家可以自行下载验证。

总的来说匈牙利算法还是很好理解的。

发表评论

电子邮件地址不会被公开。 必填项已用*标注