AnsweredAssumed Answered

TSP Subtour-Elimination Callback Class

Question asked by behnam on Sep 7, 2017
Latest reply on Sep 8, 2017 by Susanne.Heip

I am new to Xpress and trying to implement a simple TSP instance using BCL Java to get myself up to speed. I have done this in Cplex before using the cut callbacks, which work fine for adding the subtour elimination constraints to the B&B tree during the solve. However, the following code snippet does not quite work in Xpress. It is a simple 13-node TSP where root node's integer solution is as follows:

x[0][2] = 1
x[0][7] = 1
x[1][8] = 1
x[1][11] = 1
x[2][3] = 1
x[3][7] = 1
x[4][11] = 1
x[4][12] = 1
x[5][9] = 1
x[5][10] = 1
x[6][8] = 1
x[6][12] = 1
x[9][10] = 1

 

Obviously, we need to add a subtour elimination constraint, such as:

 

CUT(0):x[0][2] + x[0][7] + x[2][3] + x[3][7] <= 3

 

However, the solver declare the problem infeasible after adding this valid cut.

 

 *** Search completed ***     Time:     0 Nodes:          1
Problem is integer infeasible
Number of integer feasible solutions found is 0
Best bound is  1.00000E+40

 

What am I missing?

 

 

########################### Code ######################


import java.util.HashSet;
import java.util.Set;
import com.dashoptimization.*;

 

/**
 * Simple TSP implementation using callback classes of BCL Subtour elimination
 * constraints are added at each node of the branch and bound tree.
 *
 */
public class TSP {
    static final double EPS = 1e-6; /* Epsilon */
    static int n = 13; /* Number of cities */
    static int[][] D = { { 0, 2451, 713, 1018, 1631, 1374, 2408, 213, 2571, 875, 1420, 2145, 1972 },
            { 2451, 0, 1745, 1524, 831, 1240, 959, 2596, 403, 1589, 1374, 357, 579 },
            { 713, 1745, 0, 355, 920, 803, 1737, 851, 1858, 262, 940, 1453, 1260 },
            { 1018, 1524, 355, 0, 700, 862, 1395, 1123, 1584, 466, 1056, 1280, 987 },
            { 1631, 831, 920, 700, 0, 663, 1021, 1769, 949, 796, 879, 586, 371 },
            { 1374, 1240, 803, 862, 663, 0, 1681, 1551, 1765, 547, 225, 887, 999 },
            { 2408, 959, 1737, 1395, 1021, 1681, 0, 2493, 678, 1724, 1891, 1114, 701 },
            { 213, 2596, 851, 1123, 1769, 1551, 2493, 0, 2699, 1038, 1605, 2300, 2099 },
            { 2571, 403, 1858, 1584, 949, 1765, 678, 2699, 0, 1744, 1645, 653, 600 },
            { 875, 1589, 262, 466, 796, 547, 1724, 1038, 1744, 0, 679, 1272, 1162 },
            { 1420, 1374, 940, 1056, 879, 225, 1891, 1605, 1645, 679, 0, 1017, 1200 },
            { 2145, 357, 1453, 1280, 586, 887, 1114, 2300, 653, 1272, 1017, 0, 504 },
            { 1972, 579, 1260, 987, 371, 999, 701, 2099, 600, 1162, 1200, 504, 0 } }; /* Distance matrix */

 

    static XPRBvar[][] x; /* binary variable for each edge */
    static XPRB bcl;
    static XPRBprob p;

 

    static class SubtourEliminationCallback implements XPRSoptNodeListener {
        @Override
        public void XPRSoptNodeEvent(XPRSprob oprob, Object data, IntHolder arg2) {
            double[][] solX; /* Solution values */
            XPRBprob bprob = (XPRBprob) data;

 

            try {
                /* Get the solution values */
                bprob.beginCB(oprob);
                bprob.loadMat();
                bprob.sync(XPRB.XPRS_SOL);
                solX = new double[n][n];
                for (int i = 0; i < n; i++) {
                    for (int j = 0; j < n; j++) {
                        solX[i][j] = 0;
                        if (i < j) {
                            solX[i][j] = x[i][j].getSol();
                            // just for fun - print the solution
                            if (solX[i][j] > EPS) {
                                System.out.println("x[" + i + "][" + j + "] = 1");
                            }
                        }
                    }
                }
                /* Search for violated constraints: */
                Set<Edge> tour = findsubtour(solX);
                /*
                 * If the tour that includes node 0 (depot) does not visit all ndoes, we need to
                 * add a cut
                 */
                if (tour.size() < n) {
                    XPRBexpr le = new XPRBexpr();
                    for (Edge e : tour) {
                        le.add(x[e.origin][e.destination]);
                    }
                    XPRBcut[] c = new XPRBcut[1];
                    c[0] = bprob.newCut(le.lEql(-tour.size() + 1));
                    // just for fun - print the cut
                    c[0].print();
                    bprob.addCuts(c);
                }
                bprob.endCB();
            } catch (XPRSprobException e) {
                System.out.println("Error  " + e.getCause() + ": " + e.getMessage());
            }
        }
    }

 

    static void buildTSP() {

 

        // only define arcs from smaller node index to higher
        x = new XPRBvar[n][n];
        for (int i = 0; i < n - 1; i++) {
            for (int j = i + 1; j < n; j++) {
                x[i][j] = p.newVar("x[" + i + "][" + j + "]", XPRB.BV);
            }
        }

 

        // objective function
        XPRBexpr cobj = new XPRBexpr();
        for (int i = 0; i < n - 1; i++) {
            for (int j = i + 1; j < n; j++) {
                cobj.addTerm(D[i][j], x[i][j]);
            }
        }
        p.setObj(cobj);

 

        // Degree-2 constraints
        for (int i = 0; i < n; i++) {
            XPRBexpr le = new XPRBexpr();
            for (int j = 0; j < i; j++) {
                le.add(x[j][i]);
            }
            for (int j = i + 1; j < n; j++) {
                le.add(x[i][j]);
            }
            p.newCtr("Degree 2 Constraints for node " + i, le.eql(2));
        }
    }

 


    public static void main(String[] args) throws XPRSexception {

 

        bcl = new XPRB(); /* Initialize BCL */
        p = bcl.newProb("TSP"); /* Create a new problem */
        p.setCutMode(1); /* Switch presolve off, enable cut manager */
        p.setSense(XPRB.MINIM);
        
        buildTSP(); /* Basic problem definition */
        
        XPRS.init(); /* Initialize Xpress-Optimizer */
        XPRSprob oprob = p.getXPRSprob(); /* Get Optimizer problem */
        oprob.setIntControl(XPRS.HEURSTRATEGY, 0);
        oprob.setIntControl(XPRS.LPLOG, 0);
        oprob.setIntControl(XPRS.MIPLOG, 2);
        oprob.setIntControl(XPRS.CUTSTRATEGY, 0); /* Disable automatic cuts */
        oprob.setIntControl(XPRS.PRESOLVE, 0); /* Switch presolve off */
        
        // Instantiate subtour elimination cut callback
        SubtourEliminationCallback cb = new SubtourEliminationCallback();
        oprob.addOptNodeListener(cb, p);

 

        try {
            p.mipOptimize(""); /* Solve the problem as MIP */
            int statmip = p.getMIPStat(); /* Get the MIP problem status */

 

            if ((statmip == XPRB.MIP_SOLUTION) || (statmip == XPRB.MIP_OPTIMAL))
            /* An integer solution has been found */
            {
                System.out.println("Objective: " + p.getObjVal());
                /* Print out the solution: */
                for (int i = 0; i < n - 1; i++) {
                    for (int j = i + 1; j < n; j++) {
                        if (x[i][j].getSol() > EPS) {
                            System.out.println("x[" + i + "][" + j + "] = 1");
                        }
                    }
                }
            }
        } catch (XPRSexception e) {
            System.err.println(e.getMessage());
            System.exit(1);
        }

 

        p.finalize(); /* Delete the problem */
        p = null;
    }

 

    /**
     * Auxiliary method that returns the subtour that includes node 0 (depot)
     * @param sol
     * @return
     */
    protected static Set<Edge> findsubtour(double[][] sol) {
        int i, j;
        Set<Edge> subtourEdgeSet = new HashSet<Edge>();
        Set<Integer> subtourNodes = new HashSet<Integer>();
        Set<Integer> unprocessedNodes = new HashSet<Integer>();
        subtourNodes.add(0);
        unprocessedNodes.add(0);

 

        Set<Integer> newlyDiscoveredNodes;
        while (true) {
            newlyDiscoveredNodes = new HashSet<Integer>();
            for (Integer node : unprocessedNodes) {
                i = node.intValue();
                for (j = 0; j < i; j++) {
                    if (sol[j][i] > 0.99) {
                        subtourEdgeSet.add(new Edge(j, i));
                        if (!subtourNodes.contains(j)) {
                            newlyDiscoveredNodes.add(j);
                            subtourNodes.add(j);
                        }
                    }
                }
                for (j = i + 1; j < sol.length; j++) {
                    if (sol[i][j] > 0.99) {
                        subtourEdgeSet.add(new Edge(i, j));
                        if (!subtourNodes.contains(j)) {
                            newlyDiscoveredNodes.add(j);
                            subtourNodes.add(j);
                        }

 

                    }
                }
            }
            // if no new node was found, break
            if (newlyDiscoveredNodes.size() == 0) {
                break;
            } else {
                unprocessedNodes = new HashSet<Integer>();
                unprocessedNodes.addAll(newlyDiscoveredNodes);
            }
        }

 

        return subtourEdgeSet;
    }
}

Outcomes