Initializing help system before first use

Els - An economic lot-sizing problem solved by cut-and-branch and branch-and-cut heuristics


Type: Lot-sizing
Rating: 5 (difficult)
Description: The version 'xbels' of this example shows how to implement cut-and-branch (= cut generation at the root node of the MIP search) and 'xbelsc' implements a branch-and-cut (= cut generation at the MIP search tree nodes) algorithm using the cut manager.
File(s): xbels.cs, xbelsc.cs


xbels.cs
/********************************************************/
/*  Xpress-BCL C# Example Problems                      */
/*  ==============================                      */
/*                                                      */
/*  file xbels.cs                                       */
/*  `````````````                                       */
/*  Example for the use of Xpress-BCL                   */
/*  (Economic lot sizing, ELS, problem, solved by       */
/*   adding(l,S)-inequalities) in several rounds        */
/*   looping over the root node)                        */
/*                                                      */
/*  ELS considers production planning over a horizon    */
/*  of T periods. In period t, t=1,...,T, there is a    */
/*  given demand DEMAND[t] that must be satisfied by    */
/*  production prod[t] in period t and by inventory     */
/*  carried over from previous periods. There is a      */
/*  set-up up cost SETUPCOST[t] associated with         */
/*  production in period t.  The unit production cost   */
/*  in period t is PRODCOST[t]. There is no inventory   */
/*  or stock-holding cost.                              */
/*                                                      */
/*  (c) 2008 Fair Isaac Corporation                     */
/*      authors: S.Heipcke, D.Brett.                    */
/********************************************************/

using System;
using System.Text;
using System.IO;
using Optimizer;
using BCL;


namespace Examples
{

    public class TestAdvEls
    {
        const double EPS = 1e-6;

        const int T = 6;                       /* Number of time periods */

        /****DATA****/
        int[] DEMAND    = { 1, 3, 5, 3, 4, 2}; /* Demand per period */
        int[] SETUPCOST = {17,16,11, 6, 9, 6}; /* Setup cost per period */
        int[] PRODCOST  = { 5, 3, 2, 1, 3, 1}; /* Production cost per period */
        int[,] D = new int[T,T];               /* Total demand,periods t1-t2 */

        XPRBvar[] prod = new XPRBvar[T];       /* Production in period t */
        XPRBvar[] setup = new XPRBvar[T];      /* Setup in period t */

        XPRBprob p = new XPRBprob("Els");      /* Initialize a new BCL prob */

        /*********************************************************************/

        public void modEls()
        {
            int s,t,k;
            XPRBexpr cobj,le; 

            for(s=0;s<T;s++)
                for(t=0;t<T;t++)
                    for(k=s;k<=t;k++)
                        D[s,t] += DEMAND[k];

            /****VARIABLES****/
            for(t=0;t<T;t++)
            {
                prod[t]=p.newVar("prod" + (t+1));
                setup[t]=p.newVar("setup" + (t+1), BCLconstant.XPRB_BV);   
            }

            /****OBJECTIVE****/
            cobj = new XPRBexpr();
            for(t=0;t<T;t++)                       /* Minimize total cost */
                cobj += SETUPCOST[t]*setup[t] + PRODCOST[t]*prod[t];
            p.setObj(p.newCtr("Objective", cobj));

            /****CONSTRAINTS****/
            /* Production in period t must not exceed the total demand for the
            remaining periods; if there is production during t then there
            is a setup in t */
            for(t=0;t<T;t++)
                p.newCtr("Production", prod[t] <= D[t,T-1]*setup[t]); 

            /* Production in periods 0 to t must satisfy the total demand
            during this period of time */
            for(t=0;t<T;t++)
            {
                le = new XPRBexpr(0);
                for(s=0;s<=t;s++) le += prod[s];
                p.newCtr("Demand", le >= D[0,t]);
            }

        }

        /*********************************************************************/
        /*  Cut generation loop at the top node:                             */
        /*    solve the LP and save the basis                                */
        /*    get the solution values                                        */
        /*    identify and set up violated constraints                       */
        /*    load the modified problem and load the saved basis             */
        /*********************************************************************/
        public void solveEls()
        {
            double objval;               /* Objective value */
            int t,l;
            int starttime;
            int ncut, npass, npcut;      /* Counters for cuts and passes */
             
             /* Solution values for var.s prod & setup */
            double[] solprod = new double[T];
            double[] solsetup = new double[T];  
            
            double ds;
            XPRBbasis basis;
            XPRBexpr le;
            XPRSprob xprsp = p.getXPRSprob();

            starttime=XPRB.getTime();

            /* Disable automatic cuts - we use our own */
            xprsp.CutStrategy = 0;           
            
            /* Switch presolve off */
            xprsp.Presolve = 0;
            
            ncut = npass = 0;

            do
            {
                npass++;
                npcut = 0;
                p.lpOptimize("p");              /* Solve the LP */
                basis = p.saveBasis();      /* Save the current basis */
                objval = p.getObjVal();     /* Get the objective value */

                /* Get the solution values: */
                for(t=0;t<T;t++)
                {
                    solprod[t]=prod[t].getSol();
                    solsetup[t]=setup[t].getSol();
                }

                /* Search for violated constraints: */
                for(l=0;l<T;l++)
                {
                    for (ds=0.0, t=0; t<=l; t++)
                    {
                        if(solprod[t] < D[t,l]*solsetup[t] + EPS)  
                            ds += solprod[t];
                        else  
                            ds += D[t,l]*solsetup[t];
                    }

                    /* Add the violated inequality: the minimum of the actual
                    production prod[t] and the maximum potential production 
                    D[t][l]*setup[t] in periods 0 to l must at least equal the
                    total demand in periods 0 to l.
                    sum(t=1:l) min(prod[t], D[t][l]*setup[t]) >= D[0][l]
                    */
                    if(ds < D[0,l] - EPS) 
                    {
                        le = new XPRBexpr(0);
                        for(t=0;t<=l;t++) 
                        {
                            if (solprod[t] < D[t,l]*solsetup[t] + EPS)
                                le += prod[t];
                            else
                                le += D[t,l]*setup[t];
                        }
                        p.newCtr("cut" + (ncut+1), le >= D[0,l]);
                        ncut++; 
                        npcut++;
                    }
                }

                System.Console.WriteLine("Pass " + npass + " (" + 
                	(XPRB.getTime()-starttime)/1000.0 + 
                	" sec), objective value " + objval + ", cuts added: " +
                	npcut + " (total " + ncut + ")");

                if(npcut==0)
                    System.Console.Write("Optimal integer solution found:\n");
                else
                { 
                    p.loadMat();         /* Reload the problem */
                    p.loadBasis(basis);  /* Load the saved basis */
                    basis.reset();       /* No need to keep basis any longer */
                }
            } while(npcut>0);

            /* Print out the solution: */
            for(t=0;t<T;t++)
                System.Console.WriteLine("Period " + (t+1) + ": prod " + 
                	prod[t].getSol() + " (demand: " + DEMAND[t] + 
                	", cost: " + PRODCOST[t] + "), setup " + 
                	setup[t].getSol() + " (cost: " + SETUPCOST[t] + ")"); 
 
        }

        /*********************************************************************/

        public static void Main()
        {
           XPRB.init();
           TestAdvEls TestInstance = new TestAdvEls();

            TestInstance.modEls();       /* Model the problem */
            TestInstance.solveEls();     /* Solve the problem */

            return;
        } 

    }

}

xbelsc.cs
/********************************************************
  Xpress-BCL C# Example Problems
  ==============================

  file xbelsc.cs 
  ``````````````
  Economic lot sizing, ELS, problem, solved by adding
  (l,S)-inequalities) in a branch-and-cut heuristic 
  (using the cut manager).
  
  ELS considers production planning over a horizon
  of T periods. In period t, t=1,...,T, there is a
  given demand DEMAND[t] that must be satisfied by
  production prod[t] in period t and by inventory
  carried over from previous periods. There is a
  set-up up cost SETUPCOST[t] associated with
  production in period t. The unit production cost
  in period t is PRODCOST[t]. There is no inventory
  or stock-holding cost.

*** This model cannot be run with a Community Licence ***

  (c) 2008 Fair Isaac Corporation
      authors: S.Heipcke, D.Brett.
********************************************************/

using System;
using System.Text;
using System.IO;
using Optimizer;
using BCL;


namespace Examples
{

    public struct XPRBprobStruct
    {
        public XPRBprob prob;
        public double tol;
    };

    public class TestAdvElsCallback
    {
        const int T = 6;                             /* Number of time periods */

        /****DATA****/
        int[] DEMAND    = { 1, 3, 5, 3, 4, 2};  /* Demand per period */
        int[] SETUPCOST = {17,16,11, 6, 9, 6};  /* Setup cost per period */
        int[] PRODCOST  = { 5, 3, 2, 1, 3, 1};  /* Production cost per period */
        int[,] D = new int[T, T];               /* Total demand in periods t1 - t2 */

        XPRBvar[] prod = new XPRBvar[T];                        /* Production in period t */
        XPRBvar[] setup = new XPRBvar[T];                       /* Setup in period t */


        XPRBprob p = new XPRBprob("ElsC");                      /* Initialize a new problem in BCL */

        /***********************************************************************/

        public void modEls()
        {
            int s, t, k;
            XPRBexpr cobj, le;

            for (s = 0; s < T; s++)
                for (t = 0; t < T; t++)
                    for (k = s; k <= t; k++)
                        D[s, t] += DEMAND[k];

            /****VARIABLES****/
            for (t = 0; t < T; t++)
            {
                prod[t] = p.newVar("prod" + (t + 1));
                setup[t] = p.newVar("setup" + (t + 1), BCLconstant.XPRB_BV);
            }

            /****OBJECTIVE****/
            cobj = new XPRBexpr();
            for (t = 0; t < T; t++)                       /* Minimize total cost */
                cobj += SETUPCOST[t] * setup[t] + PRODCOST[t] * prod[t];
            p.setObj(p.newCtr(cobj));

            /****CONSTRAINTS****/
            /* Production in period t must not exceed the total demand for the
            remaining periods; if there is production during t then there
            is a setup in t */
            for (t = 0; t < T; t++)
                p.newCtr("Production", prod[t] <= D[t, T - 1] * setup[t]);

            /* Production in periods 0 to t must satisfy the total demand
            during this period of time */
            for (t = 0; t < T; t++)
            {
                le = new XPRBexpr(0);
                for (s = 0; s <= t; s++) le += prod[s];
                p.newCtr("Demand", le >= D[0, t]);
            }

        }

        /**************************************************************************/
        /*  Cut generation loop at the tree node:                                 */
        /*    get the solution values                                             */
        /*    identify and set up violated constraints                            */
        /*    add cuts to the matrix                                              */
        /**************************************************************************/
        public int cbNode(XPRSprob xprsp, object mobj)
        {
            XPRBprobStruct mo;
            double objval;                  /* Objective value */
            int t, l;
            int ncut;                       /* Counters for cuts */
            double[] solprod = new double[T];
            double[] solsetup = new double[T]; /* Solution values for var.s prod & setup */
            double ds;
            int depth, node;
            XPRBcut[] cut = new XPRBcut[T];
            XPRBexpr le;

            mo = (XPRBprobStruct)mobj;

            ncut = 0;

            depth = xprsp.NodeDepth;
            node = xprsp.Nodes;

            /* Get the solution values */
            mo.prob.beginCB(xprsp);
            mo.prob.sync(BCLconstant.XPRB_XPRS_SOL);
            for (t = 0; t < T; t++)
            {
                solprod[t] = prod[t].getSol();
                solsetup[t] = setup[t].getSol();
            }

            /* Search for violated constraints: */
            for (l = 0; l < T; l++)
            {
                for (ds = 0.0, t = 0; t <= l; t++)
                {
                    if (solprod[t] < D[t, l] * solsetup[t] + mo.tol) ds += solprod[t];
                    else ds += D[t, l] * solsetup[t];
                }

                /* Add the violated inequality: the minimum of the actual production
                prod[t] and the maximum potential production D[t][l]*setup[t] 
                in periods 0 to l must at least equal the total demand in periods 
                0 to l.
                sum(t=1:l) min(prod[t], D[t][l]*setup[t]) >= D[0][l]
                */
                if (ds < D[0, l] - mo.tol)
                {
                    le = new XPRBexpr(0);
                    for (t = 0; t <= l; t++)
                    {
                        if (solprod[t] < D[t, l] * solsetup[t] + mo.tol)
                            le += prod[t];
                        else
                            le += D[t, l] * setup[t];
                    }
                    cut[ncut] = mo.prob.newCut(le >= D[0, l]);
                    ncut++;
                }
            }

            /* Add cuts to the problem */
            if (ncut > 0)
            {
                mo.prob.addCuts(cut, ncut);
                objval = xprsp.LPObjVal;
                System.Console.Write("Cuts added : " + ncut + " (depth " + depth + ", node ");
                System.Console.WriteLine(node + "), obj. " + objval);
            }
            mo.prob.endCB();

            return 0;
        }

        /***********************************************************************/
        public void treeCutGen()
        {
            IntPtr oprob;
            XPRBprobStruct mo;
            double feastol;
            int starttime, t;
            CutmgrCallback del = new CutmgrCallback(cbNode);
            XPRSprob xprsp = p.getXPRSprob();

            starttime = XPRB.getTime();

            xprsp.LPLog = 0;
            xprsp.MIPLog = 3;

            xprsp.CutStrategy = 0;
            xprsp.Presolve = 0;
            xprsp.ExtraRows = 5000;

            feastol = xprsp.FeasTol;
            feastol *= 10;

            mo.prob = p;
            mo.tol = feastol;
            p.setCutMode(1);
            xprsp.AddCutmgrCallback(del, (object)mo);

            p.mipOptimize();                                    /* Solve the MIP */
            System.Console.Write("(" + (XPRB.getTime() - starttime) / 1000.0 + " sec) Global status ");
            System.Console.WriteLine(p.getMIPStat() + ", best solution: " + p.getObjVal());
            for (t = 0; t < T; t++)
            {
                System.Console.Write("Period " + (t + 1) + ": prod " + prod[t].getSol() + " (demand: ");
                System.Console.Write(DEMAND[t] + ", cost: " + PRODCOST[t] + "), setup ");
                System.Console.WriteLine(setup[t].getSol() + " (cost: " + SETUPCOST[t] + ")");
            }
        }

        /***********************************************************************/

        public static void Main()
        {
            XPRB.init();
            TestAdvElsCallback TestInstance = new TestAdvElsCallback();

            TestInstance.modEls();                    /* Model the problem */
            TestInstance.treeCutGen();                /* Solve the problem */

            return;
        }

    }

}