Search code examples
optimizationgraph-theorymathematical-optimizationpascalbipartite

LAPJVsp produces infeasible results during augmenting row reduction


In R. Jonker and A. Volgenant, A Shortest Augmenting Path Algorithm for Dense and Sparse Linear Assignment Problems (doi: 10.1007/BF02278710), the authors show that an implementation of their LAPJV algorithm, adapted to sparse graphs and referred to as LAPJVsp, performs well on various problems.

A Pascal implementation of LAPJVsp is currently available here. The augmenting row reduction step of the algorithm is mostly unchanged, and differs from the code provided in the published paper only by using a compressed sparse row representation of the biadjacency matrix of the graph, whose row indices, column indices, and weights are referred to as first, kk, and cc respectively:

tel:=0;
repeat
  h:=1; l0:=l; l:=0;
  while h<=l0 do
  begin
     i:=free[h]; h:=h+1; v0:=inf; vj:=inf;
     for t:=first[i] to first[i+1]-1 do
     begin
        j:=kk[t]; dj:=cc[t]-v[j];
        if dj<vj then if dj>=v0 then begin vj:=dj; j1:=j end
        else begin vj:=v0; v0:=dj; j1:=j0; j0:=j end;
     end;
     i0:=y[j0]; u[i]:=vj;
     if v0<vj then v[j0]:=v[j0]-vj+v0
              else if i0>0 then begin j0:=j1; i0:=y[j0] end;
     x[i]:=j0; y[j0]:=i;
     if i0>0 then
        if v0<vj then begin h:=h-1; free[h]:=i0 end
                 else begin l:=l+1; free[l]:=i0 end
  end;
  tel:=tel+1
until tel=2;

Now, for sparse inputs, there's no guarantee that feasible solutions to the linear assignment problem even exist, and in many cases, this infeasibility can be caught along the way, but the row reduction step given above ignores this issue, and in some cases outputs infeasible results, which then survive all the way to the end of the function.

For example, consider the bipartite graph whose biadjacency matrix is [[1, 1, 1], [*, *, 1], [*, *, 1]], in which * denotes missing edges. This can be run through the Pascal implementation as follows:

n:=3;
first[1]:=1; first[2]:=4; first[3]:=5; first[4]:=6;
kk[1]:=1; kk[2]:=2; kk[3]:=3; kk[4]:=3; kk[5]:=3;
cc[1]:=1; cc[2]:=1; cc[3]:=1; cc[4]:=1; cc[5]:=1;
zlap:=lap(x, y, u, v);

After the row reduction, the column assignments x become [1, 3, 2], and this also ends up being the end result of the function, even though the solution is clearly infeasible as there is no edge from the third row to the second column.

This leads to several questions: Is this a well-known bug? Can the algorithm be trusted to provide correct results assuming a feasible solution exists, so that one can just assume feasibility and argue that this isn't a bug? Can the row reduction step be salvaged and either provide correct results, or detect infeasibility?


Solution

  • This happens because the index of the column with the second-highest reduced cost for a given row, that is, j1, is never unset between different rows. We get rid of the error by explicitly unsetting the indices:

    ...
    tel:=0;
    repeat
       h:=1; l0:=l; l:=0;
       while h<=l0 do
       begin
          j0:=0; j1:=0;  {<-- This is new}
          i:=free[h]; h:=h+1; v0:=inf; vj:=inf;
    ...
    

    At this point, it is possible (as is indeed the case for the example given in the original post) that no value of dj for a given row is less than inf, meaning that j0 remains zero up until augmentation. Explicitly checking for this provides an infeasibility test.