/********************************************************
  BCL Example Problems
  ====================

  file xbcoco3.c
  ``````````````
  Coco Problem Phase 3.
  Introduce time periods and inventory.

  (c) 2008-2023 Fair Isaac Corporation
      author: S.Heipcke, Jan. 2000, rev. Mar. 2011
********************************************************/

#include <stdio.h>
#include "xprb.h"

#define NP 2           /* Number of products (p) */
#define NF 2           /* Factories (f) */
#define NR 2           /* Raw materials (r) */
#define NT 4           /* Time periods (t) */

#define REVFILE     XPRBDATAPATH "/coco/revt.dat"
#define CMAKEFILE   XPRBDATAPATH "/coco/cmake.dat"
#define CBUYFILE    XPRBDATAPATH "/coco/cbuyt.dat"
#define REQFILE     XPRBDATAPATH "/coco/req.dat"
#define MXSELLFILE  XPRBDATAPATH "/coco/maxsellt.dat"
#define MXMAKEFILE  XPRBDATAPATH "/coco/mxmake.dat"
#define PSTOCK0FILE XPRBDATAPATH "/coco/pstock0.dat"
#define RSTOCK0FILE XPRBDATAPATH "/coco/rstock0.dat"

/****TABLES****/
double REV[NP][NT];     /* Unit selling price of product p in period t */
double CMAK[NP][NF];    /* Unit cost to make product p at factory f */
double CBUY[NR][NT];    /* Unit cost to buy raw material r in period t */
double REQ[NP][NR];     /* Requirement by unit of prod. p for raw material r */
double MXSELL[NP][NT];  /* Max. amount of p that can be sold in period t */
double MXMAKE[NF];      /* Max. amount factory f can make over all products */
double PSTOCK0[NP][NF]; /* Initial product p stock level at factory f */
double RSTOCK0[NR][NF]; /* Initial raw material r stock level at factory f*/

/****DATA****/
double CPSTOCK = 2;   /* Unit cost to store any product p */
double CRSTOCK = 1;   /* Unit cost to store any raw material r */
double MXRSTOCK = 300;  /* Max. amount of r that can be stored each f and t */

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

void modcoco3(XPRBprob prob)
{
 XPRBvar make[NP][NF][NT], sell[NP][NF][NT], pstock[NP][NF][NT+1],
         buy[NR][NF][NT], rstock[NR][NF][NT+1];
 XPRBctr ctr;
 int p,f,r,t;

/****VARIABLES****/
 for(p=0;p<NP;p++)
  for(f=0;f<NF;f++)
  {
   for(t=0;t<NT;t++)
   {
    make[p][f][t]=XPRBnewvar(prob,XPRB_PL, XPRBnewname("make_p%d_f%d",p+1,f+1),
        0, XPRB_INFINITY);
        /* Amount of prod. p to make at factory f in period t */
    sell[p][f][t]=XPRBnewvar(prob,XPRB_PL, XPRBnewname("sell_p%d_f%d",p+1,f+1),
        0, XPRB_INFINITY);
        /* Amount of prod. p sold from factory f in period t */
   }
   for(t=0;t<NT+1;t++)
    pstock[p][f][t]=XPRBnewvar(prob,XPRB_PL, XPRBnewname("pstock_p%d_f%d",p+1,f+1),
        0, XPRB_INFINITY);
        /* Stock level of prod. p at factory f at start of period t */
  }
 for(r=0;r<NR;r++)
  for(f=0;f<NF;f++)
  {
   for(t=0;t<NT;t++)
    buy[r][f][t]=XPRBnewvar(prob,XPRB_PL, XPRBnewname("buy_r%d_f%d",r+1,f+1),
        0, XPRB_INFINITY);
        /* Amount of raw material r bought for factory f in period t */
   for(t=0;t<NT+1;t++)
    rstock[r][f][t]=XPRBnewvar(prob,XPRB_PL, XPRBnewname("rstock_r%d_f%d",r+1,f+1),
        0, XPRB_INFINITY);
        /* Stock level of raw mat. r at factory f at start of per. t */
  }

/****OBJECTIVE****/
 ctr = XPRBnewctr(prob,"OBJ",XPRB_N);
 for(f=0;f<NF;f++)                     /* Objective: maximize total profit */
 {
  for(p=0;p<NP;p++)
  {
   for(t=0;t<NT;t++)
   {
    XPRBaddterm(ctr, sell[p][f][t], REV[p][t]);
    XPRBaddterm(ctr, make[p][f][t], -1*CMAK[p][f]);
   }
   for(t=1;t<NT+1;t++)
    XPRBaddterm(ctr, pstock[p][f][t], -1*CPSTOCK);
  }
  for(r=0;r<NR;r++)
  {
   for(t=0;t<NT;t++)
    XPRBaddterm(ctr, buy[r][f][t], -1*CBUY[r][t]);
   for(t=1;t<NT+1;t++)
    XPRBaddterm(ctr, rstock[r][f][t], -1*CRSTOCK);
  }
 }
 XPRBsetobj(prob,ctr);                       /* Select objective function */

/****CONSTRAINTS****/
 for(p=0;p<NP;p++)
  for(f=0;f<NF;f++)
   for(t=0;t<NT;t++)
   {
    ctr=XPRBnewctr(prob,"PBal",XPRB_E);      /* Product stock balance */
    XPRBaddterm(ctr, pstock[p][f][t], 1);
    XPRBaddterm(ctr, make[p][f][t], 1);
    XPRBaddterm(ctr, sell[p][f][t], -1);
    XPRBaddterm(ctr, pstock[p][f][t+1], -1);
   }

 for(r=0;r<NR;r++)
  for(f=0;f<NF;f++)
   for(t=0;t<NT;t++)
   {
    ctr=XPRBnewctr(prob,"RBal",XPRB_E);      /* Raw material stock balance */
    XPRBaddterm(ctr, rstock[r][f][t], 1);
    XPRBaddterm(ctr, buy[r][f][t], 1);
    XPRBaddterm(ctr, rstock[r][f][t+1], -1);
    for(p=0;p<NP;p++)
     XPRBaddterm(ctr, make[p][f][t], -1*REQ[p][r]);
   }

 for(p=0;p<NP;p++)
  for(t=0;t<NT;t++)
  {         /* Limit on the amount of product p to be sold */
   ctr=XPRBnewctr(prob,"MxSell",XPRB_L);
   for(f=0;f<NF;f++)
    XPRBaddterm(ctr, sell[p][f][t], 1);
   XPRBaddterm(ctr, NULL, MXSELL[p][t]);
  }

 for(f=0;f<NF;f++)
  for(t=0;t<NT;t++)
  {
   ctr=XPRBnewctr(prob,"MxMake",XPRB_L);     /* Capacity limit at factory f */
   for(p=0;p<NP;p++)
    XPRBaddterm(ctr, make[p][f][t], 1);
   XPRBaddterm(ctr, NULL, MXMAKE[f]);
  }

 for(f=0;f<NF;f++)
  for(t=1;t<NT+1;t++)
  {
   ctr=XPRBnewctr(prob,"MxRStock",XPRB_L);   /* Raw material stock limit */
   for(r=0;r<NR;r++)
    XPRBaddterm(ctr, rstock[r][f][t], 1);
   XPRBaddterm(ctr, NULL, MXRSTOCK);
  }

/****BOUNDS****/
 for(p=0;p<NP;p++)
  for(f=0;f<NF;f++)
   XPRBfixvar(pstock[p][f][0], PSTOCK0[p][f]);    /* Initial product levels */

 for(r=0;r<NR;r++)
  for(f=0;f<NF;f++)
   XPRBfixvar(rstock[r][f][0], RSTOCK0[r][f]);    /* Initial raw mat. levels */

/****SOLVING + OUTPUT****/
 XPRBsetsense(prob,XPRB_MAXIM);              /* Choose sense of optimization */
 XPRBlpoptimize(prob, "");                   /* Solve the LP-problem */
 printf("Objective: %g\n", XPRBgetobjval(prob));      /* Get objective value */

        /* Uncomment to print out the solution values */
/* for(p=0;p<NP;p++)
  for(f=0;f<NF;f++)
  {
   for(t=0;t<NT;t++)
    printf("%s:%g %s:%g  ", XPRBgetvarname(make[p][f][t]),
      XPRBgetsol(make[p][f][t]), XPRBgetvarname(sell[p][f][t]),
      XPRBgetsol(sell[p][f][t]));
   for(t=0;t<NT+1;t++)
    printf("%s:%g ", XPRBgetvarname(pstock[p][f][t]), XPRBgetsol(pstock[p][f][t]));
  }
 printf("\n");

 for(r=0;r<NR;r++)
  for(f=0;f<NF;f++)
  {
   for(t=0;t<NT;t++)
    printf("%s:%g ", XPRBgetvarname(buy[r][f][t]), XPRBgetsol(buy[r][f][t]));
   for(t=0;t<NT+1;t++)
    printf("%s:%g ", XPRBgetvarname(rstock[r][f][t]), XPRBgetsol(rstock[r][f][t]));
  }
 printf("\n");
*/
}

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

    /**** Read data from files ****/
void readdata(void)
{
 FILE *datafile;
 int p,r;
        /* Read the revenu data file */
 datafile=fopen(REVFILE,"r");
 for(p=0;p<NP;p++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", REV[p], NT);
 fclose(datafile);

        /* Read the production cost data file */
 datafile=fopen(CMAKEFILE,"r");
 for(p=0;p<NP;p++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", CMAK[p], NF);
 fclose(datafile);

        /* Read the raw material cost data file */
 datafile=fopen(CBUYFILE,"r");
 for(r=0;r<NR;r++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", CBUY[r], NT);
 fclose(datafile);

        /* Read the resource requirement data file */
 datafile=fopen(REQFILE,"r");
 for(p=0;p<NP;p++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", REQ[p], NR);
 fclose(datafile);

        /* Read the max. sales quantities data file */
 datafile=fopen(MXSELLFILE,"r");
 for(p=0;p<NP;p++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", MXSELL[p], NT);
 fclose(datafile);

        /* Read the production capacities data file */
 datafile=fopen(MXMAKEFILE,"r");
 XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", MXMAKE, NF);
 fclose(datafile);

        /* Read the product stock data file */
 datafile=fopen(PSTOCK0FILE,"r");
 for(p=0;p<NP;p++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", PSTOCK0[p], NF);
 fclose(datafile);

        /* Read the raw material stock data file */
 datafile=fopen(RSTOCK0FILE,"r");
 for(r=0;r<NR;r++)
  XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", RSTOCK0[r], NF);
 fclose(datafile);
}

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

int main(int argc, char **argv)
{
 XPRBprob prob;

 prob=XPRBnewprob("Coco3");   /* Initialize a new problem in BCL */
 readdata();                  /* Data input from file */
 modcoco3(prob);              /* Model and solve the problem */

 return 0;
}