分配問題與Hungarian算法

分配問題與Hungarian算法

分配問題

指派問題

匈牙利算法

匈牙利方法是一種可以在多項式時間內解決分配問題(assignment problem)的組合優化算法。它由Harold Kuhn 與1955年發展並提出，因爲該算法很大程度上依賴於先前兩位匈牙利數學家:Denes Konig 和 Jeno Egervary,因此被命名爲「匈牙利方法」。
1957年James Munkres從新審視了這個方法，證實發現該方法是嚴格polynomial的，因此以後該方法也被稱爲Kuhn-Munkres 算法或者Munkres分配算法。原始的匈牙利算法的時間複雜度是,然而以後Edmonds和Karp，以及Tomizawa獨立發現通過必定的修改，該算法能改達到的時間複雜度。 Ford和Fulkerson將該方法擴展到通常運輸問題的求解上。 2006年，研究發現Carl Custav Jacobi在19實際就解決了assignment問題，而且在其逝世後的1890年求解過程被以拉丁語形式發表。。。

指派問題

匈牙利法解決的指派問題應該具備兩個約束條件

• workes 和tasks的數目應該相同，即o2o問題。

• 求解的是最小化問題，如工做時間的最小化、費用的最小化等等

指派問題示例：
有三個workers: Jim, Steve和Alan，如今有3個工做：clean the bathroom, sweep the floors和wash the windows須要交給這三我的，每一個人只能完成一個任務，對應的cost matrix以下

---	Clean bathroom	Sweep floors	Wash windows
Jim	$2	$3	$3
Steve	$3	$2	$3
Alan	$3	$3	$2

那麼如何分配任務是開銷最小就是一個指派問題

匈牙利法步驟

問題： 假定某單位有甲、乙、丙、丁、戊五個員工，現須要完成A、B、C、D、E五項任務，每一個員工完成某項任務的時間以下圖所示，應該如何分配任務，才能保證完成任務所須要的時間開銷最小？

1476015762594.jpg

解：

1. 寫出係數矩陣

2. 更新系數矩陣，使係數矩陣的每一行每一列都減去該行該列的最小值，保證每一行每一列都有0元素出現,參見定理2.

3. 選擇只有一個0元素的行或列將該0元素標註爲獨立0元素，並將該0元素所在的列或行中0元素劃掉，直至找不到知足條件的行或列，須要注意的是在循環時，劃掉的0元素再也不視爲0元素。好比咱們找下面這個係數矩陣的標註0元素和劃掉的0元素

那麼咱們用1表示0元素，0表示非0元素，用2標註獨立0元素，-2表示劃掉0元素，依次獲得的中間結果爲

能夠發現咱們在標註第一行第2個0元素時，並無把已經劃掉的第一個0看成0元素看待。

4. 劃蓋0線。
   a. 首先找到不含有獨立0元素的行，將行標註 (4)
   b. 找到標註行中劃掉的0元素所在的列，將列標註 (2,3)
   c. 將標註列中獨立0元素所在的行標註 (4, 1, 2)
   d. 重複b,c步驟知道全部的標註行中再也不存在劃掉的0元素
   (行：4 ,1, 2, 3 ; 列: 2, 3, 1)
   e. 在全部標註的列上作蓋0線，在全部未被標註的行上作蓋0線，咱們在矩陣中用3表示被蓋0線覆蓋，則

   能夠發現全部的0元素都被覆蓋掉，因此稱爲蓋0線。
   根據定理1，若是此時蓋0線的個數等於矩陣的維數，至關於找到n個獨立0元素，則跳到步驟7，不然步驟5更新系數矩陣。

5. 更新系數矩陣。
   a. 找到未被蓋0線覆蓋的元素中的最小值
   b. 全部未被蓋0線覆蓋的元素減去最小值
   c. 全部蓋0線交叉處的元素值加上最小值

6. 重複步驟4，5

7. 計算最優解。
   相似步驟3中找獨立0元素的方法在C中找到n個不一樣行不一樣列的0元素所對應的位置。

定理1. 係數矩陣C中獨立零元素的最多個數等於能覆蓋全部零元素的最少線數。 —— D. Konig
定理2. 若將分配問題的係數矩陣每一行及每一列分別減去各行及各列的最小元素，則新的分配問題與原分配問題有相同的最優解，只是最優值差一個常數。

固然，注意找蓋〇線的方法並非惟一的，好比下述方法

1. 同上一方法

2. 同上一方法

3. 在C中尋找未被蓋〇線覆蓋的存在0元素且數組最少的行(列)，標註一個0元素做爲獨立0元素，該0元素所在的列(行)作蓋0線。重複此步驟，直至找不到符合條件的行或列，已經被找過的行(列)就再也不找了！

4. 如果蓋0線數組等於矩陣維度，則跳到7

5. 同上一方法

6. 重複3，4，5

7. 同上一方法

固然還有別的方法，好比從包含最多0元素的行或列開始作蓋0線直到將全部的0元素覆蓋掉等

這裏代碼實現了上面兩個方法，註釋掉的是第二個方法的核心

1. function [cost,CMatrix]=Assignment(C,ismin) 

   ---

2. % Assignment problem solved by hungarian method. 

   ---

3. % 

   ---

4. % input: 

   ---

5. % C - 係數矩陣，能夠適應workers和tasks數目不一樣的情形 

   ---

6. % ismin - 1表示最小化問題，0表示最大化問題 

   ---

7. % ouput: 

   ---

8. % cost - 最終花費代價 

   ---

9. % CMatrix - 對應的匹配矩陣，元素1所在位置c_{ij}表示j task分配給 i worker。 

   ---

10. % 

    ---

11.  

    ---

12. [m,n]=size(C); 

    ---

13. if ismin==0 

    ---

14. C=max(C(:))-C; 

    ---

15. end 

    ---

16.  

    ---

17. %workes 和tasks數目不相同 

    ---

18. if m<n 

    ---

19. C=[C;zeros(n-m,n)]; 

    ---

20. elseif m>n 

    ---

21. C=[C zeros(m,m-n)]; 

    ---

22. end 

    ---

23. copyC=C; 

    ---

24. d=max(m,n);% 最終係數矩陣的維度 

    ---

25. C=C-repmat(min(C,[],2),1,d); 

    ---

26. C=C-repmat(min(C,[],1),d,1); 

    ---

27.  

    ---

28. %% 方法一 

    ---

29. while 1 

    ---

30. A=int8((C==0)); 

    ---

31. nIZeros=0; % 獨立0元素的個數 

    ---

32. while 1 

    ---

33. r=sum(A==1,2); % 每一行0元素的個數 

    ---

34. [~,idr]=find(r'==1);%找到只有一個0元素的行 

    ---

35. if ~isempty(idr) % 若是找到這樣的行 

    ---

36. tr=A(idr(1),:); 

    ---

37. [~,idc]=find(tr==1);%找到0元素所在列 

    ---

38. A(idr(1),idc)=2;%標註獨立元素 

    ---

39. tc=A(:,idc); 

    ---

40. tc(idr(1))=2; 

    ---

41. [~,idr]=find(tc'==1);%找到獨立0元素所在列的其餘0元素 

    ---

42. A(idr,idc)=-2;%劃掉獨立0元素所在列的其他0元素 

    ---

43. nIZeros=nIZeros+1; 

    ---

44. else 

    ---

45. c=sum(A==1,1); % 每一列0元素的個數 

    ---

46. [~,idc]=find(c==1);%找到只含有一個0元素的列 

    ---

47. if ~isempty(idc)% 找到這樣的列 

    ---

48. tc=A(:,idc(1)); 

    ---

49. [~,idr]=find(tc'==1);%0元素所在的行 

    ---

50. A(idr,idc(1))=2;%標註獨立0元素 

    ---

51. tr=A(idr,:); 

    ---

52. tr(idc(1))=2; 

    ---

53. [~,idc]=find(tr==1);%獨立0元素所在行的其餘0元素 

    ---

54. A(idr,idc)=-2;%劃掉獨立0元素所在行的其他0元素 

    ---

55. nIZeros=nIZeros+1; 

    ---

56. else 

    ---

57. break; 

    ---

58. end 

    ---

59. end 

    ---

60. end 

    ---

61.  

    ---

62. if nIZeros==d 

    ---

63. %計算最優解 

    ---

64. CMatrix=(A==2); 

    ---

65.  

    ---

66. if ismin==1 

    ---

67. cost=sum(copyC(:).*CMatrix(:)); 

    ---

68. else 

    ---

69. cost = sum((max(copyC(:))-copyC(:)).*CMatrix(:)); 

    ---

70. end 

    ---

71. CMatrix=CMatrix(1:m,1:n); 

    ---

72. break;%找到d個獨立0元素，則跳出循環 

    ---

73. else% 獨立0元素個數不足，就要找蓋0線了 

    ---

74. r=sum(A==2,2); 

    ---

75. [~,idr]=find(r'==0);%不含有獨立0元素的行 

    ---

76. idrr=idr; 

    ---

77. idcc=[]; 

    ---

78. while 1 

    ---

79. tr=A(idrr,:); 

    ---

80. [~,idc]=find(tr==-2);%不含獨立0元素的行中劃掉的0元素所在列 

    ---

81. if isempty(idc),break;end 

    ---

82. tc=A(:,unique(idc)); 

    ---

83. [idrr,~]=find(tc==2);%這些列中標註的0元素所在行 

    ---

84. idr=[idr,idrr']; 

    ---

85. idcc=[idcc,idc]; 

    ---

86. end 

    ---

87. idry=1:d; 

    ---

88. idry(idr)=[];%蓋0線所在的行索引 

    ---

89. TempC=C;%存儲非覆蓋元素 

    ---

90. TempC(idry,:)=[]; 

    ---

91. TempC(:,idcc)=[]; 

    ---

92. minUnOverlap=min(TempC(:)); 

    ---

93. %更新系數矩陣 

    ---

94. C=C-minUnOverlap; 

    ---

95. C(idry,:)=C(idry,:)+minUnOverlap; 

    ---

96. C(:,idcc)=C(:,idcc)+minUnOverlap; 

    ---

97. end 

    ---

98. end 

    ---

99. %% 方法二 

    ---

100. % while 1 

     ---

101. % CMatrix=zeros(d); 

     ---

102. % nLines=0; 

     ---

103. % A=(C==0); 

     ---

104. % idx=[]; 

     ---

105. % idy=[]; 

     ---

106. % sr=[]; 

     ---

107. % sc=[]; 

     ---

108. % while 1 

     ---

109. % r=sum(A,2); 

     ---

110. % c=sum(A,1); 

     ---

111. % r(sr)=0; 

     ---

112. % c(sc)=0; 

     ---

113. % trc=[r(:);c(:)]; 

     ---

114. % [trc,idtrc]=sort(trc,1,'ascend'); 

     ---

115. % [~,idn0]=find(trc'>0); 

     ---

116. % if ~isempty(idn0) 

     ---

117. % id=idtrc(idn0(1)); 

     ---

118. % if id>d 

     ---

119. % tc=A(:,id-d); 

     ---

120. % [~,idr]=find(tc'==1); 

     ---

121. % A(idr(1),:)=0;  

     ---

122. % nLines=nLines+1; 

     ---

123. % idy=[idy,idr(1)]; 

     ---

124. % CMatrix(idr(1),id-d)=1; 

     ---

125. % sc=[sc,id-d]; 

     ---

126. % else 

     ---

127. % tr=A(id,:); 

     ---

128. % [~,idc]=find(tr==1); 

     ---

129. % A(:,idc(1))=0; 

     ---

130. % nLines=nLines+1; 

     ---

131. % idx=[idx,idc(1)]; 

     ---

132. % CMatrix(id,idc(1))=1; 

     ---

133. % sr=[sr,id]; 

     ---

134. % end  

     ---

135. % else 

     ---

136. % break; 

     ---

137. % end 

     ---

138. % end 

     ---

139. % if nLines==d 

     ---

140. % if ismin 

     ---

141. % cost=sum(copyC(:).*CMatrix(:)); 

     ---

142. % else 

     ---

143. % cost=sum((max(copyC(:))-copyC(:)).*CMatrix(:)); 

     ---

144. % end 

     ---

145. % CMatrix=CMatrix(1:m,1:n); 

     ---

146. % break; 

     ---

147. % else 

     ---

148. % tempC=C; 

     ---

149. % tempC(idy,:)=[]; 

     ---

150. % tempC(:,idx)=[]; 

     ---

151. % minUnOverlap=min(tempC(:)); 

     ---

152. % C=C-minUnOverlap; 

     ---

153. % C(idy,:)=C(idy,:)+minUnOverlap; 

     ---

154. % C(:,idx)=C(:,idx)+minUnOverlap; 

     ---

155. % end 

     ---

156. % end 

     ---

157.  

     ---

158. end 

     ---

159.  

     ---

匈牙利算法的拓展

• workers數小於tasks數目
  此時能夠虛擬出若干個workers使個數相等，對於虛擬出的workers對於tasks的代價是相同的都是0.

• workers數目大於tasks數目
  能夠和上述相似，增長虛擬tasks，在係數矩陣中對應列元素都爲0

• 最大值問題
  找到係數矩陣中最大的元素，而後使用最大元素分別減去矩陣中元素，這樣最大化問題就化爲最小化問題，一樣可使用匈牙利算法。

測試

1. 如今使用代碼求解上述例題

1476016827034.jpg

2. 注意方法一中存在不足之處是必須找到僅含一個0元素的行或列，這並不科學，好比下面這個係數矩陣，就會出現死循環，因此能夠尋找含有最少0元素的行或列，這就和方法二相似了，因此我沒再實現。方法二能夠解決下面這個極端例子

咱們來使用方法二(將代碼中註釋掉的部分反註釋，方法一註釋掉)測試下這個極端的例子

1476017296572.jpg

3. 測試下面這個係數矩陣，匈牙利算法拓展1

1476017394268.jpg

結果

1476017434080.jpg

4. 對上個測試用例計算最大代價

1476017934004.jpg

綜上，相較於方法一，方法二可以處理各類情況。