/********************************************************
 * Xpress-BCL Java Example Problems
 * ================================
 *
 * file xbdlvriis2rep.java
 * ```````````````````````
 * Transportation problem (infeasible data).
 * Repairing infeasibility.
 * - Using Optimizer functions -
 *
 * (c) 2008-2024 Fair Isaac Corporation
 * author: S.Heipcke, Jan. 2008, rev. Mar. 2011
 ********************************************************/

import com.dashoptimization.*;
import java.io.*;
import java.util.*;

public class xbdlvriis2rep {
  static final int NSupp = 10; /* Number of suppliers */
  static final int NCust = 7; /* Number of customers */
  static final int MaxArcs = 100; /* Max. num. of non-zero cost values */

  static final String VANFILE = System.getProperty("XPRBDATA") + "/delivery/ifvan.dat";
  static final String COSTFILE = System.getProperty("XPRBDATA") + "/delivery/cost.dat";

  /****DATA****/
  /* Supplier:                    London  Luton  B'ham Bristl  Derby Stckpt */
  static final double SUPPLY[] = {
    140.0,
    200.0,
    50.0,
    10.0,
    400.0,
    200.0,
    /* Supplier:                     York  Derby Soton Scnthp */
    20.0,
    90.0,
    30.0,
    12.0
  };

  /* Customer:                     London Livpol Doncst   York   Hull  Manchr */
  static final double DEMAND[] = {
    1230.3,
    560.4,
    117.1,
    592.8,
    310.0,
    1247.0,
    /* Customer:                     Shffld */
    86.0
  };

  static double[][] COST; /* Cost per supplier-customer pair */
  static double[][] IFVAN; /* Non-zero if route uses vans instead
                                        of lorries */
  static final double VANCAP = 40.0; /* Capacity on routes that use vans */

  static XPRSprob op;

  /***********************************************************************/

  static void modDelivery() {
    try (XPRBprob p = new XPRBprob("Delivery"); /* Initialize BCL and create a new problem */
        XPRBexprContext context =
            new XPRBexprContext(); /* Release XPRBexpr instances at end of block. */
        XPRS xprs = new XPRS()) {
        /* Initialize Xpress-Optimizer */
      XPRBexpr lobj, lc;
      int s, c, i;
      XPRBvar[][] x;
      XPRBctr[] CSupply, CDemand;
      int nrow, ncol;
      IntHolder iscode;
      double[] lrp, grp, lbp, ubp;

      /****VARIABLES****/
      x = new XPRBvar[NSupp][NCust];
      for (s = 0; s < NSupp; s++)
        for (c = 0; c < NCust; c++) x[s][c] = p.newVar("x_s" + s + "_" + c);

      /****OBJECTIVE****/
      lobj = new XPRBexpr();
      for (s = 0; s < NSupp; s++) /* Objective: Minimize total cost */
        for (c = 0; c < NCust; c++) lobj.add(x[s][c].mul(COST[s][c]));
      p.setObj(lobj); /* Set objective function */

      /****CONSTRAINTS****/
      CDemand = new XPRBctr[NCust];
      for (c = 0; c < 4; c++) {
          /* Satisfy demand of each customer */
        lc = new XPRBexpr();
        for (s = 0; s < 5; s++) lc.add(x[s][c]);
        CDemand[c] = p.newCtr("Demand", lc.gEql(DEMAND[c]));
      }
      for (c = 4; c < NCust; c++) {
          /* Satisfy demand of each customer */
        lc = new XPRBexpr();
        for (s = 5; s < NSupp; s++) lc.add(x[s][c]);
        CDemand[c] = p.newCtr("Demand", lc.gEql(DEMAND[c]));
      }

      CSupply = new XPRBctr[NSupp];
      for (s = 0; s < 5; s++) {
          /* Keep within supply at each supplier*/
        lc = new XPRBexpr();
        for (c = 0; c < 4; c++) lc.add(x[s][c]);
        CSupply[s] = p.newCtr("Supply", lc.lEql(SUPPLY[s]));
      }
      for (s = 5; s < NSupp; s++) {
          /* Keep within supply at each supplier*/
        lc = new XPRBexpr();
        for (c = 4; c < NCust; c++) lc.add(x[s][c]);
        CSupply[s] = p.newCtr("Supply", lc.lEql(SUPPLY[s]));
      }

      /****BOUNDS****/
      for (s = 0; s < NSupp; s++)
        for (c = 0; c < NCust; c++) if (IFVAN[s][c] != 0) x[s][c].setUB(VANCAP);

      /****SOLVING + OUTPUT****/
      p.setSense(XPRB.MINIM); /* Set objective sense to minimization */
      p.lpOptimize(""); /* Solve the LP-problem */

      System.out.println("LP status: " + p.getLPStat());
      if (p.getLPStat() == XPRB.LP_OPTIMAL)
        System.out.println("Objective: " + p.getObjVal()); /* Get objective value */
      else if (p.getLPStat() == XPRB.LP_INFEAS) {
        op = p.getXPRSprob(); /* Retrieve the Optimizer problem */

        /**** Trying to fix infeasibilities ****/
        /*
          lrp: (affects = and <= rows)
          ax - aux_var  = b
          ax - aux_var <= b
          grp: (affects = and >= rows)
          ax + aux_var  = b
          ax + aux_var >= b
          lbp:
          x_i + aux_var >= l
          ubp:
          x_i - aux_var <= u
        */

        /**** Simplified infeasibility repair:
         * specifying preferences per constraint/bound type ****/

        System.out.println("\n**** Repair infeasibility:");
        iscode = new IntHolder();
        op.repairInfeas(iscode, 'c', 'o', ' ', 10, 9, 0, 20, 0.001);
        printSolution(p, iscode.value, x);

        /**** Weighted infeasibility repair:
         * specifying preferences for every constraint/bound separately ****/

        ncol = op.getIntAttrib(XPRS.ORIGINALCOLS);
        nrow = op.getIntAttrib(XPRS.ORIGINALROWS);
        lrp = new double[nrow];
        grp = new double[nrow];
        lbp = new double[ncol];
        ubp = new double[ncol];

        /* Relax bounds due to van capacity */
        /* Repairweightedinfeas for upper bounds of concerned flow variables */

        for (s = 0; s < NSupp; s++)
          for (c = 0; c < NCust; c++) if (IFVAN[s][c] != 0) ubp[x[s][c].getColNum()] = 20;

        System.out.println("\n**** Relax van capacity:");
        op.repairWeightedInfeas(iscode, lrp, grp, lbp, ubp, 'd', 0.001, "");
        printSolution(p, iscode.value, x);

        /* Relax supply limits (may buy in additional quantities) */
        /* Repairinfeas for 'less or equal' side of Supply constraints */

        for (s = 0; s < NSupp; s++) lrp[CSupply[s].getRowNum()] = 10;

        System.out.println("\n**** Relax supply limits:");
        op.repairWeightedInfeas(iscode, lrp, grp, lbp, ubp, 'd', 0.001, "");
        printSolution2(op, iscode.value, x);

        /* Relax demand constraints (may not satisfy all customers) */
        /* Repairinfeas for 'greater or equal' side of Demand constraints */

        for (c = 0; c < NCust; c++) grp[CDemand[c].getRowNum()] = 9;

        System.out.println("\n**** Relax demand constraints:");
        op.repairWeightedInfeas(iscode, lrp, grp, lbp, ubp, 'd', 0.001, "");
        printSolution2(op, iscode.value, x);
      }
    }
  }

  /***********************************************************************/

  /**** Initialize the stream tokenizer ****/
  static StreamTokenizer initST(FileReader file) {
    StreamTokenizer st = null;

    st = new StreamTokenizer(file); /* Initialize the stream tokenizer */
    st.commentChar('!'); /* Use the character '!' for comments */
    st.eolIsSignificant(true); /* Return end-of-line character */
    st.ordinaryChar(','); /* Use ',' as separator */
    st.parseNumbers(); /* Read numbers as numbers (not strings)*/
    return st;
  }

  /**** Read multi line data file ****/
  static void readMultiLineFile(String filename, double[][] data) throws IOException {
    FileReader datafile = null;
    StreamTokenizer st;
    int i = 0, j = 0;

    datafile = new FileReader(filename);
    st = initST(datafile);
    do {
      do {
        st.nextToken();
      } while (st.ttype == st.TT_EOL); /* Skip empty lines */
      while (st.ttype == st.TT_NUMBER) {
        data[j][i++] = st.nval;
        if (st.nextToken() != ',') {
          j++;
          i = 0;
          break;
        }
        st.nextToken();
      }
    } while (st.ttype == st.TT_EOL);
    datafile.close();
  }

  /**** Read data from files ****/
  static void readData() throws IOException {
    int s, c;

    COST = new double[NSupp][NCust];
    IFVAN = new double[NSupp][NCust];

    /* Initialize data tables to 0: in the van data file some entries that are
     * zero are simply left out (that is, the lines are incomplete) */
    for (s = 0; s < NSupp; s++)
      for (c = 0; c < NCust; c++) {
        COST[s][c] = 0;
        IFVAN[s][c] = 0;
      }
    /* Read the demand data file */
    readMultiLineFile(COSTFILE, COST);

    /* Read the van data file */
    readMultiLineFile(VANFILE, IFVAN);
  }

  /**** Print out the solution values ****/
  static void printSolution(XPRBprob p, int scode, XPRBvar[][] x) {
    int s, c;
    double sup, dem;
    java.lang.String rstat[] = {
      "relaxed optimum found",
      "relaxed problem infeasible",
      "relaxed problem unbounded",
      "solution nonoptimal for original objective",
      "error",
      "numerical instability"
    };

    System.out.println("Status: " + rstat[scode]);
    if (scode == 0) {
      p.sync(XPRB.XPRS_SOL);
      for (s = 0; s < NSupp; s++)
        for (c = 0; c < NCust; c++)
          if (x[s][c].getSol() > 0.01)
            System.out.print(x[s][c].getName() + ":" + x[s][c].getSol() + " ");
      System.out.println();
      System.out.println("Violations:");
      for (c = 0; c < NCust; c++) {
        sup = 0;
        for (s = 0; s < NSupp; s++) sup += x[s][c].getSol();
        if (sup < DEMAND[c]) System.out.println("  Customer " + c + ": " + (DEMAND[c] - sup));
      }
      for (s = 0; s < NSupp; s++) {
        dem = 0;
        for (c = 0; c < NCust; c++) dem += x[s][c].getSol();
        if (dem > SUPPLY[s]) System.out.println("  Supplier " + s + ": " + (dem - SUPPLY[s]));
      }
      for (s = 0; s < NSupp; s++)
        for (c = 0; c < NCust; c++)
          if (IFVAN[s][c] != 0 && VANCAP < x[s][c].getSol())
            System.out.println("  Van " + s + "-" + c + ": " + (x[s][c].getSol() - VANCAP));
    }
  }

  static void printSolution2(XPRSprob op, int scode, XPRBvar[][] x) {
    int s, c, ncol;
    double sup, dem;
    java.lang.String rstat[] = {
      "relaxed optimum found",
      "relaxed problem infeasible",
      "relaxed problem unbounded",
      "solution nonoptimal for original objective",
      "error",
      "numerical instability"
    };
    double[] sol;

    System.out.println("Status: " + rstat[scode]);
    if (scode == 0) {
      ncol = op.getIntAttrib(XPRS.ORIGINALCOLS);
      sol = new double[ncol];
      op.getLpSol(sol, null, null, null);

      for (s = 0; s < NSupp; s++)
        for (c = 0; c < NCust; c++)
          if (x[s][c].getColNum() >= 0 && sol[x[s][c].getColNum()] > 0.01)
            System.out.print(x[s][c].getName() + ":" + sol[x[s][c].getColNum()] + " ");
      System.out.println();
      System.out.println("Violations:");
      for (c = 0; c < NCust; c++) {
        sup = 0;
        for (s = 0; s < NSupp; s++) if (x[s][c].getColNum() >= 0) sup += sol[x[s][c].getColNum()];
        if (sup < DEMAND[c]) System.out.println("  Customer " + c + ": " + (DEMAND[c] - sup));
      }
      for (s = 0; s < NSupp; s++) {
        dem = 0;
        for (c = 0; c < NCust; c++) if (x[s][c].getColNum() >= 0) dem += sol[x[s][c].getColNum()];
        if (dem > SUPPLY[s]) System.out.println("  Supplier " + s + ": " + (dem - SUPPLY[s]));
      }
      for (s = 0; s < NSupp; s++)
        for (c = 0; c < NCust; c++)
          if (IFVAN[s][c] != 0 && x[s][c].getColNum() >= 0 && VANCAP < sol[x[s][c].getColNum()])
            System.out.println("  Van " + s + "-" + c + ": " + (sol[x[s][c].getColNum()] - VANCAP));
    }
  }

  /***********************************************************************/

  public static void main(String[] args) {
    try {
      readData(); /* Data input from file */
    } catch (IOException e) {
      System.err.println(e.getMessage());
      System.exit(1);
    }
    modDelivery(); /* Problem formulation and solution */
  }
}
