Quadratic Assignment Problems: solution enumeration

In the last few posts we have given code for computing the Gilmore Lawler bound for quadratic assignment problems, and described two different branching techniques.   In this post we show how to enumerate solutions of a small QAP.  The possible solutions to a QAP range over permutations over unassigned facilities and locations.   If lower bounds do not increase after repeated branching, we may eventually assign nearly all the facilities to locations.  If there are only 2 or 3 remaining, we may as well just enumerate all the possible solutions and pick the best.   In practice, we hope that we will rarely have to do this, otherwise the branch-and-bound tree will be too large to process efficiently.  Here is a procedure that enumerates over all solutions of a size 3 QAP.  Iterating over permutations is a standard textbook problem, so maybe you can write cleaner code than I did here!

/// <summary>
/// Assume that all facilities have been assigned to locations,
/// except three. EnumerateSolutions enumerates all six remaining assignments
/// returns the best.
/// </summary>
private double EnumerateSolutions(Qap q, BranchAndBoundNode node, int[] pBest) {
double objBest;
int nPerm = perm3.GetLength(0);
int n = perm3.GetLength(1);
double[] obj = new double[nPerm];
int[] r = node.p.UnusedIndices().ToArray();
int[] c = node.pInv.UnusedIndices().ToArray();
int[] p = (int[])node.p.Clone();
Debug.Assert(r.Length == n, “Only call me for size “ + n + ” QAPs.”);

// Try all the possible permutations.
for (int i = 0; i < nPerm; i++) {
for (int j = 0; j < n; j++) {
p[r[j]] = c[perm3[i, j]];
}
obj[i] = q.Evaluate(p);
}

// Get the best permutation and objective.
int best = obj.ArgMin(out objBest);
p.CopyTo(pBest, 0);
for (int j = 0; j < n; j++) {
pBest[r[j]] = perm3[best, j];
}
return objBest;
}

If you combine the code from the following posts you will see that we almost have a working B&B solver:

• Linear assignment problem solver.
• QAP object.
• QAP Gilmore-Lawler bound.
• Branch-and-bound skeleton.
• Vector and matrix utilities.
• Branching.
• Strong and weak branching.

In my next post I will fill in a couple of the gaps, and then I will re-post all of the code.  Then we’ll think about what we can actually do with a QAP solver!

Quadratic Assignment Problems: strong and weak branching

Let’s pick up where we left off last time and write a score-based Branch delegate.  It will be used as the basis for both of our “real” branching functions.

 

private BranchingDecision BranchCore(BranchAndBoundNode node, double bound, double[][] S) {
RowColumnSums(S, node.Size);
double rowVal, colVal;
int rowBest = _rowSum.ArgMax(out rowVal);
int colBest = _colSum.ArgMax(out colVal);
if (rowVal > colVal) {
return new BranchingDecision(true, rowBest);
}
else {
return new BranchingDecision(false, colBest);
}
}

 

private void RowColumnSums(double[][] U, int size) {
if (_rowSum == null) {
_rowSum = new double[_qap.Size];
_colSum = new double[_qap.Size];
}
_rowSum.ConstantFill(Double.MinValue);
_colSum.ConstantFill(Double.MinValue);

for (int i = 0; i < size; i++) {
for (int j = 0; j < size; j++) {
_rowSum[i] += U[i][j];
_colSum[j] += U[i][j];
}
}
}

 

If we can create a reasonable score matrix and pass it into BranchCore, then we’ve got something.  Two common techniques for QAP are:

1. Weak branching. Let S_ij be the reduced cost U_ij that we got from computing the lower bound.  (U_ij is a lower bound on the amount that the bound on the subproblem will increase.  So large U_ij means the subproblem will be a lot easier.)
2. Strong branching. Let S_ij be the bound for subproblem ij.  The bigger the bound, the easier the subproblem.

Weak branching is trivial to implement for the Gilmore-Lawler (or GLB) bound.  As I mentioned a long time ago, the reduced cost matrix essentially comes for free.

 

private BranchingDecision WeakBranch(BranchAndBoundNode node, double bound, double[][] U) {
return BranchCore(node, bound, U);
}

For strong branching we need to compute bounds for each subproblem.  We overwrite the U matrix and pass it into BranchCore.  Since the U matrix is used to fathom subproblems later on in the code, we subtract the parent bound from the bound of each subproblem.

 

private BranchingDecision StrongBranch(BranchAndBoundNode node, double bound, double[][] U) {
Qap r = new Qap(node.Size – 1);
double[][] Ur = MatrixUtilities.NewMatrix(node.Size – 1, node.Size – 1);
for (int i = 0; i < node.Size; i++) {
for (int j = 0; j < node.Size; j++) {
Qap.Reduce(node.Qap, r, i, j);
U[i][j] = GLB.Bound(r, Ur) – bound;
}
}
return BranchCore(node, bound, U);
}

 

Strong branching in general produces better guesses that weak branching.  However, it comes at a computational cost.  Therefore I want to introduce one last complication…er…feature.  Since strong branching is smarter but more costly, we should only use it where it is most needed.  I want to have a set of rules that governs when to use strong or weak branching.  Branching decisions are most important at the top of the tree, or on nodes where the bound sucks.  We can measure the suckiness of the bound for a node by computing the relative gap:  the ratio of the gap between the bound and the incumbent solution, and the corresponding gap at the root.   So our branching rule data structure will have Depth and Gap properties, as well as properties that store the Branch delegate information.

 

/// <summary>A description of a branching rule: when and how.
/// </summary>
public class BranchingRule {
/// <summary>The ID for this rule.
/// </summary>
public int Index { get; set; }
/// <summary>Maximum depth for which the rule applies.
/// </summary>
public int Depth { get; set; }
/// <summary>Maximum relative gap for which the rule applies.
/// </summary>
public double Gap { get; set; }
/// <summary>Strong or Weak.
/// </summary>
public BranchType Type { get; set; }
/// <summary>StrongBranch() or WeakBranch().
/// </summary>
internal Branch Rule { get; set; }

public BranchingRule(int index, double gap, int depth, BranchType type) {
Index = index;
Depth = depth;
Gap = gap;
Type = type;
}

public override string ToString() {
return “Depth = “ + Depth + “, Gap = “ + Gap + “, Type = “ + Type.ToString();
}
}

We will have an array of BranchingRule called _branchingRule.  Given a node, we will select the first rule in the array whose depth is larger than the node depth, and whose gap is larger than the node gap.

 

private Branch SelectBranchingRule(BranchAndBoundNode node) {
double gap = _results.RelativeGap(node.LowerBound);
BranchingRule criterion = _branchingRules.First(entry => (node.Level <= entry.Depth) && (gap >= entry.Gap));
node.RuleIndex = criterion.Index;
return criterion.Rule;
}

 

I am cheating here because I have not defined the local variable _results.  We’ll get to that in a few posts.  For now, just assume that RelativeGap is defined as I described above: (Objective – bound) / (Objective – RootBound).

At long last, that covers branching.  There are many variations, but they are essentially variations of strong or weak branching.  At this point we only have two more main topics to cover: what to do when you get to the bottom of the tree, and how do you capture and update results.