/********************************************************
 * Xpress-BCL Java Example Problems
 * ================================
 *
 * file xbcutstk.java
 * ``````````````````
 * Cutting stock problem, solved by column (= cutting
 * pattern) generation heuristic looping over the
 * root node.
 *
 * (c) 2008-2024 Fair Isaac Corporation
 * author: S.Heipcke, 2001, rev. Mar. 2014
 ********************************************************/

import com.dashoptimization.*;

public class xbcutstk {
  static final int NWIDTHS = 5;
  static final int MAXWIDTH = 94;

  static final double EPS = 1e-6;
  static final int MAXCOL = 10;

  /****DATA****/
  static final double[] WIDTH = {17, 21, 22.5, 24, 29.5}; /* Possible widths */

  static final int[] DEMAND = {150, 96, 48, 108, 227}; /* Demand per width */
  static int[][] PATTERNS; /* (Basic) cutting patterns */

  static XPRBvar[] pat; /* Rolls per pattern */
  static XPRBctr[] dem; /* Demand constraints */
  static XPRBctr cobj; /* Objective function */

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

  static void modCutStock(XPRBprob p) throws XPRSexception {
    int i, j;
    XPRBexpr le;

    PATTERNS = new int[NWIDTHS][NWIDTHS];
    for (j = 0; j < NWIDTHS; j++) PATTERNS[j][j] = (int) Math.floor((double) MAXWIDTH / WIDTH[j]);

    /****VARIABLES****/
    pat = new XPRBvar[NWIDTHS + MAXCOL];
    for (j = 0; j < NWIDTHS; j++)
      pat[j] =
          p.newVar(
              "pat_" + (j + 1), XPRB.UI, 0, (int) Math.ceil((double) DEMAND[j] / PATTERNS[j][j]));

    /****OBJECTIVE****/
    le = new XPRBexpr();
    for (j = 0; j < NWIDTHS; j++) le.add(pat[j]);
    cobj = p.newCtr("OBJ", le);
    p.setObj(cobj); /* Minimize total number of rolls */

    /****CONSTRAINTS****/
    dem = new XPRBctr[NWIDTHS];
    for (i = 0; i < NWIDTHS; i++) {
        /* Satisfy the demand per width */
      le = new XPRBexpr();
      for (j = 0; j < NWIDTHS; j++) le.add(pat[j].mul(PATTERNS[i][j]));
      dem[i] = p.newCtr("Demand", le.gEql(DEMAND[i]));
    }
  }

  /**************************************************************************/
  /*  Column generation loop at the top node:                               */
  /*    solve the LP and save the basis                                     */
  /*    get the solution values                                             */
  /*    generate new column(s) (=cutting pattern)                           */
  /*    load the modified problem and load the saved basis                  */
  /**************************************************************************/
  static void solveCutStock(XPRB bcl, XPRBprob p) throws XPRSexception {
    double objval = 0; /* Objective value */
    int i, j;
    int starttime;
    int npatt, npass; /* Counters for columns and passes */
    double[] solpat; /* Solution values for variables pat */
    double[] dualdem; /* Dual values of demand constraints */
    XPRBbasis basis;
    double dw, z;
    int[] x;

    x = new int[NWIDTHS];
    solpat = new double[NWIDTHS + MAXCOL];
    dualdem = new double[NWIDTHS];

    starttime = XPRB.getTime();
    npatt = NWIDTHS;

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

      /* Get the solution values: */
      for (j = 0; j < npatt; j++) solpat[j] = pat[j].getSol();
      for (i = 0; i < NWIDTHS; i++) dualdem[i] = dem[i].getDual();

      /* Solve integer knapsack problem  z = min{cx : ax<=r, x in Z^n}
      with r=MAXWIDTH, n=NWIDTHS */
      z = knapsack(bcl, NWIDTHS, dualdem, WIDTH, (double) MAXWIDTH, DEMAND, x);
      /* Force garbage collection to clean up subproblem: */
      /*   System.gc();
           System.runFinalization();
      */
      System.out.print(
          "(" + (XPRB.getTime() - starttime) / 1000.0 + " sec) Pass " + (npass + 1) + ": ");

      if (z < 1 + EPS) {
        System.out.println("no profitable column found.");
        System.out.println();
        basis = null; /* No need to keep the basis any longer */
        break;
      } else {
        /* Print the new pattern: */
        System.out.println("new pattern found with marginal cost " + (z - 1));
        System.out.print("   Widths distribution: ");
        dw = 0;
        for (i = 0; i < NWIDTHS; i++) {
          System.out.print(WIDTH[i] + ":" + x[i] + "  ");
          dw += WIDTH[i] * x[i];
        }
        System.out.println("Total width: " + dw);

        /* Create a new variable for this pattern: */
        pat[npatt] = p.newVar("pat_" + (npatt + 1), XPRB.UI);

        cobj.add(pat[npatt]); /* Add new var. to the objective */
        dw = 0;
        for (i = 0; i < NWIDTHS; i++) /* Add new var. to demand constraints*/
          if (x[i] > EPS) {
            dem[i].add(pat[npatt].mul(x[i]));
            if ((int) Math.ceil((double) DEMAND[i] / x[i]) > dw)
              dw = (int) Math.ceil((double) DEMAND[i] / x[i]);
          }
        pat[npatt].setUB(dw); /* Change the upper bound on the new var.*/
        npatt++;

        p.loadMat(); /* Reload the problem */
        p.loadBasis(basis); /* Load the saved basis */
        basis = null; /* No need to keep the basis any longer */
      }
    }

    p.mipOptimize(""); /* Solve the MIP */

    System.out.println(
        "("
            + (XPRB.getTime() - starttime) / 1000.0
            + " sec) Optimal solution: "
            + p.getObjVal()
            + " rolls, "
            + npatt
            + " patterns");
    System.out.print("   Rolls per pattern: ");
    for (i = 0; i < npatt; i++) System.out.print(pat[i].getSol() + ", ");
    System.out.println();
  }

  /**************************************************************************/
  /* Integer Knapsack Algorithm for solving the integer knapsack problem    */
  /*    z = max{cx : ax <= R, x <= d, x in Z^N}                             */
  /*                                                                        */
  /* Input data:                                                            */
  /*   N:        Number of item types                                       */
  /*   c[i]:     Unit profit of item type i, i=1..n                         */
  /*   a[i]:     Unit resource use of item type i, i=1..n                   */
  /*   R:        Total resource available                                   */
  /*   d[i]:     Demand for item type i, i=1..n                             */
  /* Return values:                                                         */
  /*   xbest[i]: Number of items of type i used in optimal solution         */
  /*   zbest:    Value of optimal solution                                  */
  /**************************************************************************/
  static double knapsack(XPRB bcl, int N, double[] c, double[] a, double R, int[] d, int[] xbest) {
    /*
      System.out.println("Solving z = max{cx : ax <= b; x in Z^n}");
      System.out.print("   c   =");
      for(j = 0; j < N; j++)  System.out.print(" " + c[j] + ",");
      System.out.println();
      System.out.print("   a   =");
      for(j = 0; j < N; j++)  System.out.print(" " + a[j] + ",");
      System.out.println();
      System.out.print("   c/a =");
      for (j = 0; j < N; j++)  System.out.print(" " + c[j]/a[j] + ",");
      System.out.println();
      System.out.print("   b   = " + R);
    */

    try (XPRBprob pk = bcl.newProb("Knapsack")) {
      XPRBvar[] x = new XPRBvar[N];
      for (int j = 0; j < N; j++) x[j] = pk.newVar("x", XPRB.UI, 0, d[j]);

      XPRBexpr klobj = new XPRBexpr();
      for (int j = 0; j < N; j++) klobj.add(x[j].mul(c[j]));
      pk.setObj(klobj);
      XPRBexpr knap = new XPRBexpr();
      for (int j = 0; j < N; j++) knap.add(x[j].mul(a[j]));
      pk.newCtr("KnapXPRBctr", knap.lEql(R));
      pk.setSense(XPRB.MAXIM);
      pk.mipOptimize("");

      double zbest = pk.getObjVal();
      for (int j = 0; j < N; j++) xbest[j] = (int) Math.floor(x[j].getSol() + 0.5);
      /*
        System.out.println("   z = " + zbest);
        System.out.print("   x =");
        for(j=0; j<N; j++) System.out.print(" " + xbest[j] + ",");
        System.out.println();
      */

      return (zbest);
    }
  }

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

  public static void main(String[] args) throws XPRSexception {
    try (XPRB bcl = new XPRB(); /* Initialize BCL */
        XPRBprob p = bcl.newProb("Els"); /* Create a new problem */
        XPRBexprContext context =
            new XPRBexprContext(); /* Release XPRBexpr instances at end of block. */
        XPRS xprs = new XPRS()) {
        /* Initialize Xpress-Optimizer */
      modCutStock(p); /* Model the problem */
      solveCutStock(bcl, p); /* Solve the problem */
    }
  }
}