/********************************************************
  Xpress-BCL Java Example Problems
  ================================

  file xbdlvriis2rep.java
  ```````````````````````
  Transportation problem (infeasible data).
  Repairing infeasibility.
  - Using Optimizer functions -

  (c) 2008 Fair Isaac Corporation
      author: S.Heipcke, Jan. 2008, rev. Mar. 2011
********************************************************/

import java.io.*;
import java.util.*;
import java.text.DecimalFormat;
import com.dashoptimization.*;

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() {
        XPRBexpr lobj, lc;
        int s,c,i;
        XPRBvar[][] x;
        XPRBctr[] CSupply, CDemand;
        int nrow, ncol;
        IntHolder iscode;
        double[] lrp, grp, lbp, ubp;

        try (XPRBprob p = new XPRBprob("Delivery"); /* Initialize BCL and create a new problem */
             XPRS xprs = new XPRS()) {              /* Initialize Xpress-Optimizer */

            /****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 */
    }
}
