Initializing help system before first use

ELS - Solving several model instances in parallel


Type: Lot sizing
Rating: 5 (difficult)
Description: This implementation (program: runelsd.* starting submodel: elsd.mos) extends the parallel version of the ELS model to solving of submodels distributed to various computing nodes, showing the following features:
  • parallel remote execution of submodels
  • communication between different models (for bound updates on the objective function)
  • sending and receiving events
  • stopping submodels
File(s): runelsd.c, runelsd.java, elsd.mos (submodel), readelsdem.mos (submodel)
Data file(s): els.dat


runelsd.c
/*******************************************************
   Mosel Example Problems 
   ======================

   file runelsd.c
   ``````````````
   Run several instances of the model elsd.mos in
   parallel and coordinate the solution update.

   Before running this program, you need to set up the array
   NodeNames with machine names/addresses of your local network.
   All nodes that are used need to have the same version of
   Xpress installed and suitably licensed, and the server 
   "xprmsrv" must have been started on these machines.
   
   All files are local to the root node, no write access is
   required at remote nodes.

  *** The model started by this program cannot be run with 
      a Community Licence for the provided data instance ***

   (c) 2012 Fair Isaac Corporation
       author: S. Heipcke, Apr 2012, rev. Jul. 2014
*******************************************************/

#include <stdio.h>
#include <stdlib.h>
#include <limits.h>
#include <string.h>
#include "xprd.h"
#include "bindrv.h"

#define DATAFILE "els.dat"
#define T 15                   /* Time periods */
#define P 4                    /* No. of products */
#define I 2                    /* No. of (remote) Mosel instances */
#define M 6                    /* No. of Mosel models */

#define NEWSOL 2               /* Identifier for "sol. found" event */
#define M1 1
#define M2 3

/* Setting up remote connections:
 * Use machine names within your local network, IP addresses, or 
 * empty string for the current node running this model */
static char *NodeNames[] = {"", "localhost"};

double **solprod;              /* Sol. values for var.s produce */
double **solsetup;             /* Sol. values for var.s setup */
double **DEMAND;               /* Demand per period */

static void *dummy_malloc(size_t s,void *ctx);

/***************** Reading result data ******************/
/**** Read a table of the form 'double ar[][]' ****/
/*
 static int assignTableVal(s_bindrvctx bdrv, double **ar) 
 {
   int i,j,ctrl;
   union { int i; double r; char b; char *str; BINDRV_LONG l;} val;

   bindrv_getctrl(bdrv,&(val.i));
   if (val.i==BINDRV_CTRL_OPENLST)
   {
     while(bindrv_nexttoken(bdrv)>=0) {
     bindrv_getctrl(bdrv,&(val.i));
     ctrl=val.i;
     if(ctrl==BINDRV_CTRL_CLOSELST) break;
     else
       if(ctrl==BINDRV_CTRL_OPENNDX)
       {
         bindrv_getint(bdrv,&(val.i));
         i=val.i;
         bindrv_getint(bdrv,&(val.i));
         j=val.i;
	 bindrv_getctrl(bdrv,&(val.i));
         if(val.i==BINDRV_CTRL_CLOSENDX)
         {
           bindrv_getreal(bdrv,&(val.r));
           ar[i-1][j-1]=val.r;
         }
         else
         {
           printf("Wrong file format. ')' expected.\n");
	   return 1;
         }
       }
       else
       {
         printf("Wrong file format. '(' expected.\n");
	 return 2;
       }
     }
   }
   else
   { 
     printf("Wrong file format. '[' expected.\n");
     return 3;
   }  
   return 0;  
 }
*/

/**** Fixed size version for reading 'double ar[][]' ****/
 static int assignFixedTableVal(s_bindrvctx bdrv, double **ar) 
 {
   int i,j,ctrl,j2;
   union { int i; double r; char b; char *str; BINDRV_LONG l;} val;

   bindrv_getctrl(bdrv,&(val.i));
   if (val.i==BINDRV_CTRL_OPENLST)
   {
     while(bindrv_nexttoken(bdrv)>=0) {
     bindrv_getctrl(bdrv,&(val.i));
     ctrl=val.i;
     if(ctrl==BINDRV_CTRL_CLOSELST) break;
     else
       if(ctrl==BINDRV_CTRL_OPENNDX)
       {
         bindrv_getint(bdrv,&(val.i));
         i=val.i;
         bindrv_getint(bdrv,&(val.i));
         j=val.i;
	 bindrv_getctrl(bdrv,&(val.i));
         if(val.i==BINDRV_CTRL_CLOSENDX)
         {
           for(j2=j;j2<=T;j2++) 
           {
             bindrv_getreal(bdrv,&(val.r));
             ar[i-1][j2-1]=val.r;
           }
         }
         else
         {
           printf("Wrong file format. ')' expected.\n");
	   return 1;
         }
       }
       else
       {
         printf("Wrong file format. '(' expected.\n");
	 return 2;
       }
     }
   }
   else
   {
     printf("Wrong file format. '[' expected.\n");
     return 3;
   }  
   return 0;  
 }
 
/**** Identify the label to read ****/
 static int readSol(char *filename) 
 {
   s_bindrvctx bdrv;
   FILE *f;
   union { int i; double r; char b; char *str; BINDRV_LONG l;} val;
   int res=0;

   f=fopen(filename,"r");                        /* Use Mosel bin reader */
   bdrv=bindrv_newreader((size_t (*)(void *,size_t,size_t,void*))fread,f);

 /* By default the library will allocate strings (returned by 'getstring') */
 /* using 'malloc'. The following routine allows to replace 'malloc' by some */
 /* user-defined memory allocation routine (that gets a context as an */
 /* extra parameter) */
   bindrv_setalloc(bdrv,dummy_malloc,NULL);

   while(bindrv_nexttoken(bdrv)>=0 && res==0) {
     bindrv_getctrl(bdrv,&(val.i));
     if(val.i==BINDRV_CTRL_LABEL)
     {
       bindrv_getstring(bdrv,&(val.str));
       printf("Reading %s\n", val.str);
       if (strcmp(val.str,"solprod")==0)
         res=assignFixedTableVal(bdrv, solprod);
       else if (strcmp(val.str,"solsetup")==0)
         res=assignFixedTableVal(bdrv, solsetup);
       else if (strcmp(val.str,"DEMAND")==0)
         res=assignFixedTableVal(bdrv, DEMAND);
       /* Unknow labels are simply ignored */
     }
     else
     {
       printf("Wrong file format\n");
       res=4;
     }  
   }

   bindrv_delete(bdrv);
   fclose(f);
   
   return res;
 }

/**** A dummy malloc function using a static region ****/
static void *dummy_malloc(size_t s,void *ctx)
{
 static char buf[1024];

 return (s>1024)?NULL:buf;
}

/****** Run a submodel to read demand data for pretty solution printing ******/
 static int readDem(XPRDcontext xprd,XPRDmosel moselInst)
 {
   int e;
   char params[200];
   XPRDmodel modDem = NULL;                /* Model */

                                           /* Compile the model */
   if((e=XPRDcompmod(moselInst, "", "rmt:readelsdem.mos", 
                                    "rmt:readelsdem.bim", ""))!=0)
   {
     printf("Compilation failed: %d\n", e);
     return 2;
   }
                                          /* Load model */
   modDem = XPRDloadmod(moselInst, "rmt:readelsdem.bim");
   if(modDem==NULL)
   {
     printf("Loading failed\n");
     return 3;
   }  

   sprintf(params, "DATAFILE=%s,T=%d,P=%d", DATAFILE, T, P); 
   if(XPRDrunmod(modDem, params)!=0)       /* Run the model */
   {
     printf("Failed to run model\n");
     return 4;
   }
   XPRDwaitevent(xprd,-1);                 /* Wait for termination event */

   readSol("Demand");                      /* Read the demand data */

   remove("readelsdem.bim");               /* Cleaning up temporary files */
   remove("Demand");

   return 0;
 }

/***************** Main ******************/

 int main(int argv,char *args[])
 {
  XPRDcontext xprd;
  XPRDmosel moselInst[I];                  /* Mosel instances */
  XPRDmodel modELS[M],evsender;            /* Models */
  char msg[256];
  char params[200];
  double objval,evvalue,sp;
  int algsol,algopt,p,e,i,m,t,evclass,senderid;

  DEMAND = (double **)malloc(P * sizeof(double *));
  for(p=0;p<P;p++) DEMAND[p] = (double *)malloc(T * sizeof(double));
  solprod = (double **)malloc(P * sizeof(double *));
  for(p=0;p<P;p++) solprod[p] = (double *)malloc(T * sizeof(double));
  solsetup = (double **)malloc(P * sizeof(double *));
  for(p=0;p<P;p++) solsetup[p] = (double *)malloc(T * sizeof(double));

  xprd=XPRDinit();                        /* Create an XPRD context */

  for(i=0;i<I;i++)                        /* Establish remote connections */
  {
    moselInst[i]=XPRDconnect(xprd, NodeNames[i], NULL, NULL, msg,sizeof(msg));       if(moselInst[i]==NULL)
    {
      printf("Connection failed: %s\n", msg);
      return 1;
    }  
  }
                                          /* Compile the model */
  if((e=XPRDcompmod(moselInst[0], "", "rmt:elsd.mos", "rmt:elsd.bim", ""))!=0)
  {
    printf("Compilation failed: %d\n", e);
    return 2;
  }
      
                                         /* Load models */
  modELS[M1-1] = XPRDloadmod(moselInst[0], "rmt:elsd.bim");
  modELS[M2-1] = XPRDloadmod(moselInst[1], "rmt:elsd.bim");
  if((modELS[M1-1]==NULL)||(modELS[M2-1]==NULL))
  {
    printf("Loading failed\n");
    return 3;
  }  

                                         /* Run the models */
  sprintf(params, "ALG=%d,DATAFILE=%s,T=%d,P=%d", M1, DATAFILE, T, P); 
  if(XPRDrunmod(modELS[M1-1], params)!=0)
  {
    printf("Failed to run model %d\n", M1);
    return 4;
  }
  sprintf(params, "ALG=%d,DATAFILE=%s,T=%d,P=%d", M2, DATAFILE, T, P); 
  if(XPRDrunmod(modELS[M2-1], params)!=0)
  {
    printf("Failed to run model %d\n", M2);
    return 4;
  }

  objval = 1.0e+20;  
  algsol = -1; algopt = -1;

  do
  {
    XPRDwaitevent(xprd,-1);        /* Wait for the next event */
    XPRDgetevent(xprd, &evsender, &evclass, &evvalue);    /* Get the event */
    senderid = (evsender==modELS[M1-1]?M1:M2);
    if (evclass==NEWSOL)           /* Test the event class */
    {
      if (evvalue<objval)          /* Value of the event (= obj. value) */
      {                            /* ID of model sending the event */
        algsol = senderid;
	objval = evvalue;
        printf("Improved solution %g found by model %d\n", objval, algsol);
	for(m=0;m<M;m++)
          if (m!=algsol && (m==M1-1 || m==M2-1))
            XPRDsendevent(modELS[m], NEWSOL, objval);
      }
      else
        printf("Solution %g found by model %d\n", evvalue, senderid);
    }
  } while (evclass!=XPRD_EVENT_END);        /* No model has finished */
                           
  algopt = senderid;               /* Retrieve ID of terminated model */
  for(m=0;m<M;m++)
    if (m==M1-1 || m==M2-1)
      XPRDstoprunmod(modELS[m]);   /* Stop all running models */

  while((XPRDgetstatus(modELS[M1-1])==XPRD_RT_RUNNING)||
        (XPRDgetstatus(modELS[M2-1])==XPRD_RT_RUNNING)||
        !XPRDqueueempty(xprd))
  {
    XPRDwaitevent(xprd,-1);        /* Wait for the next event */
    XPRDdropevent(xprd);           /* Ignore termination event */
  }

  if (algsol==-1)
  {
    printf("No solution available");
    return 1;
  }
  else                     
  {                  /* Retrieve the best solution from shared memory */
    sprintf(params, "sol_%d", algsol);
    if((e=readSol(params))!=0) 
    {
      printf("%d Could not read file %d\n", e, algsol);
      return 5;
    }
                     /* Read the demand data by running another Mosel model */
    if((e=readDem(xprd,moselInst[0]))!=0) 
    {
      printf("Could not read demand data - %d\n", e);
      return 6;
    }
  
  /* Solution printing */
    printf("Best solution found by model %d\n", algsol); 
    printf("Optimality proven by model %d\n", algopt); 
    printf("Objective value: %g\n", objval); 
    printf("Period  setup\n");
    for(p=0;p<P;p++) printf("      %d", (p+1));
    for(t=0;t<T;t++) 
    {	  
      sp=0;
      for(p=0;p<P;p++) sp+=solsetup[p][t];
      printf("\n %2d%8.0f    ", (t+1), sp);
      for(p=0;p<P;p++)
        printf("%3.0f (%1.0f)", solprod[p][t], DEMAND[p][t]);
    }
    printf("\n");
  }


/* Cleaning up temporary files */
  for(m=0;m<M;m++) 
    if (m==M1 || m==M2)
    {
      sprintf(params, "sol_%d", m);
      remove(params);
    }
  remove("elsd.bim");

  for(i=0;i<I;i++)                /* Terminate remote connections */
    XPRDdisconnect(moselInst[i]); 
  XPRDfinish(xprd);               /* Terminate XPRD */
  
  return 0;
}

runelsd.java
/*******************************************************
   Mosel Example Problems 
   ======================

   file runelsd.java
   `````````````````
   Run several instances of the model elsd.mos in
   parallel and coordinate the solution update.

   Before running this program, you need to set up the array
   NodeNames with machine names/addresses of your local network.
   All nodes that are used need to have the same version of
   Xpress installed and suitably licensed, and the server 
   "xprmsrv" must have been started on these machines.
   
   All files are local to the root node, no write access is
   required at remote nodes.

   *** The model started by this program cannot be run with 
       a Community Licence for the provided data instance ***

   (c) 2012 Fair Isaac Corporation
       author: S. Heipcke, Jan 2012, rev. Jul. 2014
*******************************************************/

// javac -cp "xprd.jar;bindrv.jar;." runelsd.java

import com.dashoptimization.*;
import java.lang.*;
import java.io.*;

public class runelsd
{
 static final String DATAFILE = "els.dat";
 static final int T = 15;                  // Time periods
 static final int P = 4;                   // No. of products
 static final int I = 2;                   // No. of (remote) Mosel instances
 static final int M = 6;                   // No. of Mosel models

 static final int NEWSOL = 2;              // Identifier for "sol. found" event
 static final int M1 = 1;
 static final int M2 = 3;

/* Setting up remote connections:
 * Use machine names within your local network, IP addresses, or 
 * empty string for the current node running this model */
 static final String NodeNames[] = {"", "localhost"};

 static double[][] solprod;               // Sol. values for var.s produce 
 static double[][] solsetup;              // Sol. values for var.s setup
 static double[][] DEMAND = new double[P][T];    // Demand per period


/***************** Reading result data ******************/
/**** Read a table of the form 'double ar[][]' ****/
/*
 static void assignTableVal(BinDrvReader bdrv, double[][] ar) throws IOException
 {
   int i,j,ctrl;

   if (bdrv.getControl()==BinDrvReader.CTRL_OPENLST)
   {
     while(bdrv.nextToken()>=0) {
     ctrl=bdrv.getControl();
     if(ctrl==BinDrvReader.CTRL_CLOSELST) break;
     else
       if(ctrl==BinDrvReader.CTRL_OPENNDX)
       {
         i=bdrv.getInt();
         j=bdrv.getInt();
// System.out.println(i + " " + j);
         if(bdrv.getControl()==BinDrvReader.CTRL_CLOSENDX)
         {
           ar[i-1][j-1]=bdrv.getReal();
// System.out.println(ar[i-1][j-1]);
         }
         else
           throw new IOException("Wrong file format. ')' expected.");
       }
       else
         throw new IOException("Wrong file format. '(' expected.");
     }
   }
   else
     throw new IOException("Wrong file format. '[' expected.");
 }
*/

/**** Fixed size version for reading 'double ar[][]' ****/
 static void assignFixedTableVal(BinDrvReader bdrv, double[][] ar) throws IOException
 {
   int i,j,ctrl,i2,j2;

   if (bdrv.getControl()==BinDrvReader.CTRL_OPENLST)
   {
     while(bdrv.nextToken()>=0) {
     ctrl=bdrv.getControl();
     if(ctrl==BinDrvReader.CTRL_CLOSELST) break;
     else
       if(ctrl==BinDrvReader.CTRL_OPENNDX)
       {
         i=bdrv.getInt();
         j=bdrv.getInt();
         if(bdrv.getControl()==BinDrvReader.CTRL_CLOSENDX)
         {
           for(j2=j;j2<=T;j2++) 
           {
             ar[i-1][j2-1]=bdrv.getReal();
// System.out.println(i + " " + j2 + " " + ar[i-1][j2-1]);
           }
         }
         else
           throw new IOException("Wrong file format. ')' expected.");
       }
       else
         throw new IOException("Wrong file format. '(' expected.");
     }
   }
   else
     throw new IOException("Wrong file format. '[' expected.");
 }
/**** Identify the label to read ****/
 static void readSol(String filename) throws IOException
 {
   BinDrvReader bdrv;
   FileInputStream f;
   String s;

   f=new FileInputStream(filename);
   bdrv=new BinDrvReader(f);              // Use Mosel bin reader

   while(bdrv.nextToken()>=0) {
     if(bdrv.getControl()==BinDrvReader.CTRL_LABEL)
     {
       s=bdrv.getString();
       System.out.println("Reading "+s);
       if (s.equals("solprod"))  assignFixedTableVal(bdrv, solprod);
       else if (s.equals("solsetup"))  assignFixedTableVal(bdrv, solsetup);
       else if (s.equals("DEMAND"))  assignFixedTableVal(bdrv, DEMAND);
       // Unknow labels are simply ignored
     }
     else
       throw new IOException("Wrong file format");
   }

   f.close();
 }

/****** Run a submodel to read demand data for pretty solution printing ******/
 static void readDem(XPRD xprd,XPRDMosel moselInst) throws IOException
 {
   XPRDModel modDem = null;                // Model

   try {                                   // Compile the model
     moselInst.compile("", "rmt:readelsdem.mos", "rmt:readelsdem.bim");
   }
   catch(Exception e) {
     System.out.println("Compilation failed: " + e);
     System.exit(2);
   }
   try {                                   // Load model
    modDem=moselInst.loadModel("rmt:readelsdem.bim");
   }
   catch(IOException e) {
     System.out.println("Loading failed: " + e);
     System.exit(3);
   }  
   modDem.execParams = "DATAFILE="+DATAFILE + ",T="+T + ",P="+P; 
   try {
     modDem.run();                         // Run the model
   }
   catch(Exception e)
   {
     System.out.println("Failed to run model: "+e);
     System.exit(4);
   }
   xprd.waitForEvent();                    // Wait for termination event

   readSol("Demand");                      // Read the demand data

   new File("readelsdem.bim").delete();    // Cleaning up temporary files
   new File("Demand").delete();
 }


/***************** Main ******************/

 public static void main(String[] args)
 {
  XPRD xprd = new XPRD();                  // Initialize XPRD
  XPRDMosel[] moselInst = new XPRDMosel[I];  // Mosel instances
  XPRDModel[] modELS = new XPRDModel[M];   // Models
  XPRDEvent Msg;                           // Messages sent by models 
  solprod = new double[P][T];
  solsetup = new double[P][T];
  double objval;
  int algsol,algopt,senderid=0;

  try {
    for(int i=0;i<I;i++)                   // Establish remote connections
      moselInst[i]=xprd.connect(NodeNames[i]);
  }
  catch(IOException e) {
    System.out.println("Connection failed: " + e);
    System.exit(1);
  }

  try {                                    // Compile the model
    moselInst[0].compile("", "rmt:elsd.mos", "rmt:elsd.bim");
  }
  catch(Exception e) {
    System.out.println("Compilation failed: " + e);
    System.exit(2);
  }
      
  try {                                    // Load models
    modELS[M1-1]=moselInst[0].loadModel("rmt:elsd.bim");
    modELS[M2-1]=moselInst[1].loadModel("rmt:elsd.bim");
  }
  catch(IOException e) {
    System.out.println("Loading failed: " + e);
    System.exit(3);
  }  

                                            // Run-time parameters
  modELS[M1-1].execParams = "ALG="+M1 + ",DATAFILE="+DATAFILE + 
	                    ",T="+T + ",P="+P;
  modELS[M2-1].execParams = "ALG="+M2 + ",DATAFILE="+DATAFILE + 
	                    ",T="+T + ",P="+P;
  try {
    modELS[M1-1].run();                     // Run the models
    modELS[M2-1].run();
  }
  catch(Exception e)
  {
   System.out.println("Failed to run model: "+e);
   System.exit(4);
  }

  objval = 1.0e+20;
  algsol = -1; algopt = -1;

  do
  {
    xprd.waitForEvent();           // Wait for the next event
    Msg = xprd.getNextEvent();     // Get the event
    if (Msg.eventClass==NEWSOL)    // Get the event class
    {                              // Identify the sending model
      senderid = Msg.sender.getNumber()==modELS[M1-1].getNumber()?M1:M2;
      if (Msg.value<objval)        // Value of the event (= obj. value)
      {                            
        algsol = senderid;
	objval = Msg.value;
        System.out.println("Improved solution " + objval + 
			   " found by model " + algsol);
	for(int m=0;m<M;m++)
          if (m!=algsol-1 && (m==M1-1 || m==M2-1))
            modELS[m].sendEvent(NEWSOL, objval);
      }
      else
        System.out.println("Solution " + Msg.value + 
			   " found by model " + senderid );
    }

  } while (Msg.eventClass!=XPRDEvent.EVENT_END);   // No model has finished
  
  algopt = senderid;              // Store ID of terminated model
  for(int m=0;m<M;m++)
    if (m==M1-1 || m==M2-1)
      modELS[m].stop();           // Stop all running models

  while((modELS[M1-1].getExecStatus()==XPRDModel.RT_RUNNING)||
        (modELS[M2-1].getExecStatus()==XPRDModel.RT_RUNNING)||
        !xprd.isQueueEmpty() )
  {
    xprd.waitForEvent();          // Wait for the next event
    xprd.dropNextEvent();         // Ignore termination event
  }

  if (algsol==-1)
  {
    System.out.println("No solution available");
    System.exit(1);
  }
  else                     
  {                  // Retrieve the best solution from shared memory
    try {
      readSol("sol_"+algsol); 
    } catch(IOException e) {
      System.out.println(e + "Could not read file " + algsol);
      System.exit(5);
    }
                     // Read the demand data by running another Mosel model
    try {
      readDem(xprd,moselInst[0]); 
    } catch(IOException e) {
      System.out.println(e+"Could not read demand data");
      System.exit(6);
    }
  
  // Solution printing
    System.out.println("Best solution found by model " + algsol); 
    System.out.println("Optimality proven by model " + algopt); 
    System.out.println("Objective value: " + objval); 
    System.out.print("Period  setup");
    for(int p=0;p<P;p++) System.out.print("      "+(p+1));
    for(int t=0;t<T;t++) 
    {	  
      double sp=0;
      for(int p=0;p<P;p++) sp+=solsetup[p][t];
      System.out.print("\n " + String.format("%2d",(t+1)) + 
                       String.format("%8.0f",sp) + "    ");
      for(int p=0;p<P;p++)

        System.out.print(String.format("%3.0f",solprod[p][t]) +
	                 " (" + String.format("%1.0f",DEMAND[p][t]) + ")");
    }
    System.out.println();
  }


// Cleaning up temporary files
  for(int m=0;m<M;m++) new File("sol_"+m).delete();
  new File("elsd.bim").delete();

  for(int i=0;i<I;i++)                   // Terminate remote connections
    moselInst[i].disconnect(); 
 }
}

elsd.mos
(!*******************************************************
   Mosel Example Problems
   ======================

   file elsd.mos
   `````````````
   Economic lot sizing, ELS, problem
   (Cut generation algorithm adding (l,S)-inequalities 
    in one or several rounds at the root node or in 
    tree nodes)
   
   -- Distributed computing version --

   *** Not intended to be run standalone - run from runelsd.*  ***
    
   ELS considers production planning over a horizon
   of T periods. In period t, t=1,...,T, there is a
   given demand DEMAND[p,t] that must be satisfied by
   production produce[p,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[p,t]. There is no inventory
   or stock-holding cost.
   
   (c) 2010 Fair Isaac Corporation
       author: S. Heipcke, May 2010, rev. Sep. 2014
  *******************************************************!)

model Els
 uses "mmxprs","mmsystem","mmjobs"

 parameters
  ALG = 0                              ! Default algorithm: no user cuts
  CUTDEPTH = 3                         ! Maximum tree depth for cut generation
  DATAFILE = "els4.dat"
  T = 60
  P = 4
  RMT = "rmt:"                         ! Files are on root node
 end-parameters 

 forward public function cb_node: boolean
 forward procedure tree_cut_gen
 forward public function cb_updatebnd: boolean
 forward public procedure cb_intsol
 
 declarations
  NEWSOL = 2                           ! "New solution" event class
  EPS = 1e-6                           ! Zero tolerance
  TIMES = 1..T                         ! Time periods
  PRODUCTS = 1..P                      ! Set of products

  DEMAND: array(PRODUCTS,TIMES) of integer  ! Demand per period
  SETUPCOST: array(TIMES) of integer        ! Setup cost per period
  PRODCOST: array(PRODUCTS,TIMES) of real   ! Production cost per period
  CAP: array(TIMES) of integer              ! Production capacity per period
  D: array(PRODUCTS,TIMES,TIMES) of integer ! Total demand in periods t1 - t2

  produce: array(PRODUCTS,TIMES) of mpvar   ! Production in period t
  setup: array(PRODUCTS,TIMES) of mpvar     ! Setup in period t

  solprod: array(PRODUCTS,TIMES) of real    ! Sol. values for var.s produce 
  solsetup: array(PRODUCTS,TIMES) of real   ! Sol. values for var.s setup
  starttime: real
 
  Msg: Event                           ! An event
 end-declarations

 initializations from RMT+DATAFILE
  DEMAND SETUPCOST PRODCOST CAP
 end-initializations

 forall(p in PRODUCTS,s,t in TIMES) D(p,s,t):= sum(k in s..t) DEMAND(p,k)

! Objective: minimize total cost
 MinCost:= sum(t in TIMES) (SETUPCOST(t) * sum(p in PRODUCTS) setup(p,t) + 
                            sum(p in PRODUCTS) PRODCOST(p,t) * produce(p,t) )

! 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
 forall(p in PRODUCTS, t in TIMES) 
  ProdSetup(p,t):= produce(p,t) <= D(p,t,getlast(TIMES)) * setup(p,t)

! Production in periods 0 to t must satisfy the total demand
! during this period of time
 forall(p in PRODUCTS,t in TIMES) 
   sum(s in 1..t) produce(p,s) >= sum (s in 1..t) DEMAND(p,s)

! Capacity limits
 forall(t in TIMES) sum(p in PRODUCTS) produce(p,t) <= CAP(t)

! Variables setup are 0/1
 forall(p in PRODUCTS, t in TIMES) setup(p,t) is_binary 

 setparam("zerotol", EPS/100)           ! Set Mosel comparison tolerance      
 starttime:=gettime

 setparam("XPRS_THREADS", 1)            ! No parallel threads for optimization

! Uncomment to get detailed MIP output
! setparam("XPRS_VERBOSE", true)
 setparam("XPRS_LPLOG", 0)
 setparam("XPRS_MIPLOG", -1000)

 writeln("**************ALG=",ALG,"***************")

 SEVERALROUNDS:=false; TOPONLY:=false

 case ALG of
  1: do
      setparam("XPRS_CUTSTRATEGY", 0)   ! No cuts
      setparam("XPRS_HEURSTRATEGY", 0)  ! No heuristics
     end-do
  2: do 
      setparam("XPRS_CUTSTRATEGY", 0)   ! No cuts
      setparam("XPRS_HEURSTRATEGY", 0)  ! No heuristics  
      setparam("XPRS_PRESOLVE", 0)      ! No presolve
     end-do
  3: tree_cut_gen                       ! User branch-and-cut (single round)
  4: do                                 ! User branch-and-cut (several rounds)
      tree_cut_gen
      SEVERALROUNDS:=true
     end-do
  5: do                                 ! User cut-and-branch (several rounds)
      tree_cut_gen
      SEVERALROUNDS:=true
      TOPONLY:=true
     end-do
 end-case

! Parallel setup
 setcallback(XPRS_CB_PRENODE, "cb_updatebnd") 
 setcallback(XPRS_CB_INTSOL, "cb_intsol")

! Solve the problem 
 minimize(MinCost) 


!*************************************************************************
!  Cut generation loop:                                 
!    get the solution values                                            
!    identify violated constraints and add them as cuts to the problem
!    re-solve the modified problem                 
!*************************************************************************
 public function cb_node:boolean
  declarations
   ncut:integer                         ! Counter for cuts
   cut: array(range) of linctr          ! Cuts
   cutid: array(range) of integer       ! Cut type identification
   type: array(range) of integer        ! Cut constraint type 
   objval,ds: real
  end-declarations

  depth:=getparam("XPRS_NODEDEPTH")

  if((TOPONLY and depth<1) or (not TOPONLY and depth<=CUTDEPTH)) then
   ncut:=0 

 ! Get the solution values
   forall(t in TIMES, p in PRODUCTS) do
     solprod(p,t):=getsol(produce(p,t))
     solsetup(p,t):=getsol(setup(p,t))
   end-do
  
 ! Search for violated constraints
   forall(p in PRODUCTS,l in TIMES) do
    ds:=0 
    forall(t in 1..l)
      if(solprod(p,t) < D(p,t,l)*solsetup(p,t) + EPS) then ds += solprod(p,t)
      else  ds += D(p,t,l)*solsetup(p,t)
      end-if
  
      ! Add the violated inequality: the minimum of the actual production
      ! produce(p,t) and the maximum potential production D(p,t,l)*setup(t) 
      ! in periods 1 to l must at least equal the total demand in periods 
      ! 1 to l.
      ! sum(t=1:l) min(produce(p,t), D(p,t,l)*setup(p,t)) >= D(p,1,l)
    
    if(ds < D(p,1,l) - EPS) then
      cut(ncut):= sum(t in 1..l) 
       if(solprod(p,t)<(D(p,t,l)*solsetup(p,t))+EPS, produce(p,t), 
          D(p,t,l)*setup(p,t)) - D(p,1,l)
      cutid(ncut):= 1
      type(ncut):= CT_GEQ
      ncut+=1
    end-if   
  end-do

  returned:=false                     ! Call this function once per node
   
 ! Add cuts to the problem
   if(ncut>0) then 
    addcuts(cutid, type, cut);  
    writeln("Model ", ALG, ": Cuts added : ", ncut, 
            " (depth ", depth, ", node ", getparam("XPRS_NODES"), 
            ", obj. ", getparam("XPRS_LPOBJVAL"), ")")
    if SEVERALROUNDS then
     returned:=true                   ! Repeat until no new cuts generated
    end-if 
   end-if
  end-if
 end-function

! ****Optimizer settings for using the cut manager****

 procedure tree_cut_gen
  setparam("XPRS_HEURSTRATEGY", 0)    ! Switch heuristics off
  setparam("XPRS_CUTSTRATEGY", 0)     ! Switch automatic cuts off
  setparam("XPRS_PRESOLVE", 0)        ! Switch presolve off
  setparam("XPRS_EXTRAROWS", 5000)    ! Reserve extra rows in matrix

  setcallback(XPRS_CB_CUTMGR, "cb_node")  ! Define the cut manager callback
 end-procedure

!*************************************************************************
!  Setup for parallel solving:                                 
!    check whether cutoff update required at every node                
!    store and communicate any new solution found                 
!*************************************************************************
! Update cutoff value
 public function cb_updatebnd: boolean 
  if not isqueueempty then
   repeat
    Msg:= getnextevent
   until isqueueempty
   newcutoff:= getvalue(Msg)
   cutoff:= getparam("XPRS_MIPABSCUTOFF")
   writeln("Model ", ALG, ": New cutoff: ", newcutoff, 
           " old: ", cutoff)
   if newcutoff<cutoff then
    setparam("XPRS_MIPABSCUTOFF", newcutoff)
   end-if
   if (newcutoff < getparam("XPRS_LPOBJVAL")) then
    returned:= true                    ! Node becomes infeasible
   end-if
  end-if
 end-function

! Store and communicate new solution
 public procedure cb_intsol
  objval:= getparam("XPRS_LPOBJVAL")   ! Retrieve current objective value
  cutoff:= getparam("XPRS_MIPABSCUTOFF")

  writeln("Model ", ALG, ": Solution: ", objval, " cutoff: ", cutoff)
  if(cutoff > objval) then
   setparam("XPRS_MIPABSCUTOFF", objval)
  end-if
  
 ! Get the solution values
  forall(t in TIMES, p in PRODUCTS) do
   solprod(p,t):=getsol(produce(p,t))
   solsetup(p,t):=getsol(setup(p,t))
  end-do
  
 ! Store the solution in shared memory  
  initializations to "bin:" + RMT + "sol_"+ALG
   solprod
   solsetup
  end-initializations
  
 ! Send "solution found" signal  
  send(NEWSOL, objval) 
 end-procedure
 
end-model



Economic Lot-Sizing (ELS)
=========================

A well-known class of valid inequalities for ELS are the
(l,S)-inequalities.  Let D(pq) denote the total demand in periods p 
to q and y(t) be a binary variable indicating whether there is any 
production in period t.  For each period l and each subset of periods S 
of 1 to l, the (l,S)-inequality is

    sum (t=1:l | t in S) x(t) + sum (t=1:l | t not in S) D(tl) * y(t)
        >= D(1l)

It says that actual production x(t) in periods included S plus maximum 
potential production D(tl)*y(t) in the remaining periods (those not in 
S) must at least equal total demand in periods 1 to l.  Note that in 
period t at most D(tl) production is required to meet demand up to 
period l.

Based on the observation that

    sum (t=1:l | t in S) x(t) + sum (t=1:l | t not in S) D(tl) * y(t)
        >= sum (t=1:l) min(x(t), D(tl) * y(t))
        >= D(1l)

it is easy to develop a separation algorithm and thus a cutting plane
algorithm based on these (l,S)-inequalities.

readelsdem.mos
(!*******************************************************
   Mosel Example Problems
   ======================

   file readelsdem.mos
   ```````````````````
   Read the ELS demand data and write them to memory.

   *** Not intended to be run standalone - run from runelsd.*  ***
   
   (c) 2012 Fair Isaac Corporation
       author: S. Heipcke, Feb 2012
  *******************************************************!)

model "Read Els Demand"
 uses "mmjobs"

 parameters
  DATAFILE = "els4.dat"
  T = 60
  P = 4
  RMT = "rmt:"                         ! Files are on root node
 end-parameters 
 
 declarations
  TIMES = 1..T                         ! Time periods
  PRODUCTS = 1..P                      ! Set of products

  DEMAND: array(PRODUCTS,TIMES) of integer  ! Demand per period
 end-declarations

 initializations from RMT+DATAFILE
  DEMAND
 end-initializations

 initializations to "bin:"+RMT+"Demand"
  DEMAND
 end-initializations
 
end-model


© 2001-2020 Fair Isaac Corporation. All rights reserved. This documentation is the property of Fair Isaac Corporation (“FICO”). Receipt or possession of this documentation does not convey rights to disclose, reproduce, make derivative works, use, or allow others to use it except solely for internal evaluation purposes to determine whether to purchase a license to the software described in this documentation, or as otherwise set forth in a written software license agreement between you and FICO (or a FICO affiliate). Use of this documentation and the software described in it must conform strictly to the foregoing permitted uses, and no other use is permitted.