Examples with UniformIntDist umontreal.ssj.probdist.UniformIntDist used on opensource projects

Example 1 with UniformIntDist

use of umontreal.ssj.probdist.UniformIntDist in project Stochastic-Inventory by RobinChen121.

the class GetPmfMulti method getPmf.

public double[][] getPmf(int t) {
    stepSize = 1;
    if (t > 4)
        stepSize = 4;
    if (distributionGeneral[0][t] instanceof NormalDist) {
        NormalDist distribution1 = new NormalDist(distributionGeneral[0][t].getMean(), distributionGeneral[0][t].getStandardDeviation());
        NormalDist distribution2 = new NormalDist(distributionGeneral[1][t].getMean(), distributionGeneral[1][t].getStandardDeviation());
        double[] supportLB = new double[2];
        double[] supportUB = new double[2];
        supportLB[0] = (int) distribution1.inverseF(1 - truncationQuantile);
        supportUB[0] = (int) distribution1.inverseF(truncationQuantile);
        supportLB[1] = (int) distribution2.inverseF(1 - truncationQuantile);
        supportUB[1] = (int) distribution2.inverseF(truncationQuantile);
        int demandLength1 = (int) ((supportUB[0] - supportLB[0] + 1) / stepSize);
        int demandLength2 = (int) ((supportUB[1] - supportLB[1] + 1) / stepSize);
        double[][] pmf = new double[demandLength1 * demandLength2][3];
        int index = 0;
        double probilitySum = (2 * truncationQuantile - 1) * (2 * truncationQuantile - 1);
        for (int i = 0; i < demandLength1; i++) for (int j = 0; j < demandLength2; j++) {
            pmf[index][0] = supportLB[0] + i * stepSize;
            pmf[index][1] = supportLB[1] + j * stepSize;
            pmf[index][2] = (distribution1.cdf(pmf[index][0] + 0.5 * stepSize) - distribution1.cdf(pmf[index][0] - 0.5 * stepSize)) * (distribution2.cdf(pmf[index][1] + 0.5 * stepSize) - distribution2.cdf(pmf[index][1] - 0.5 * stepSize)) / probilitySum;
            index++;
        }
        return pmf;
    }
    if (distributionGeneral[0][t] instanceof GammaDist) {
        double mean1 = distributionGeneral[0][t].getMean();
        double mean2 = distributionGeneral[1][t].getMean();
        double variance1 = distributionGeneral[0][t].getVariance();
        double variance2 = distributionGeneral[1][t].getVariance();
        double scale1 = mean1 / variance1;
        double shape1 = mean1 * scale1;
        double scale2 = mean2 / variance2;
        double shape2 = mean2 * scale2;
        GammaDist distribution1 = new GammaDist(shape1, scale1);
        GammaDist distribution2 = new GammaDist(shape2, scale2);
        double[] supportLB = new double[2];
        double[] supportUB = new double[2];
        supportLB[0] = (int) distribution1.inverseF(1 - truncationQuantile);
        supportUB[0] = (int) distribution1.inverseF(truncationQuantile);
        supportLB[1] = (int) distribution2.inverseF(1 - truncationQuantile);
        supportUB[1] = (int) distribution2.inverseF(truncationQuantile);
        int demandLength1 = (int) ((supportUB[0] - supportLB[0] + 1) / stepSize);
        int demandLength2 = (int) ((supportUB[1] - supportLB[1] + 1) / stepSize);
        double[][] pmf = new double[demandLength1 * demandLength2][3];
        int index = 0;
        double probilitySum = (2 * truncationQuantile - 1) * (2 * truncationQuantile - 1);
        for (int i = 0; i < demandLength1; i++) for (int j = 0; j < demandLength2; j++) {
            pmf[index][0] = supportLB[0] + i * stepSize;
            pmf[index][1] = supportLB[1] + j * stepSize;
            pmf[index][2] = (distribution1.cdf(pmf[index][0] + 0.5 * stepSize) - distribution1.cdf(pmf[index][0] - 0.5 * stepSize)) * (distribution2.cdf(pmf[index][1] + 0.5 * stepSize) - distribution2.cdf(pmf[index][1] - 0.5 * stepSize)) / probilitySum;
            index++;
        }
        return pmf;
    }
    if (distributionGeneral[0][t] instanceof PoissonDist) {
        PoissonDist distribution1 = new PoissonDist(distributionGeneral[0][t].getMean());
        PoissonDist distribution2 = new PoissonDist(distributionGeneral[1][t].getMean());
        double[] supportLB = new double[2];
        double[] supportUB = new double[2];
        supportLB[0] = (int) distribution1.inverseF(1 - truncationQuantile);
        supportUB[0] = (int) distribution1.inverseF(truncationQuantile);
        supportLB[1] = (int) distribution2.inverseF(1 - truncationQuantile);
        supportUB[1] = (int) distribution2.inverseF(truncationQuantile);
        int demandLength1 = (int) ((supportUB[0] - supportLB[0] + 1) / stepSize);
        int demandLength2 = (int) ((supportUB[1] - supportLB[1] + 1) / stepSize);
        double[][] pmf = new double[demandLength1 * demandLength2][3];
        int index = 0;
        double probilitySum = (2 * truncationQuantile - 1) * (2 * truncationQuantile - 1);
        for (int i = 0; i < demandLength1; i++) for (int j = 0; j < demandLength2; j++) {
            pmf[index][0] = supportLB[0] + i * stepSize;
            pmf[index][1] = supportLB[1] + j * stepSize;
            pmf[index][2] = distribution1.prob(i) * distribution2.prob(j) / probilitySum;
            index++;
        }
        return pmf;
    }
    if (distributionGeneral[0][t] instanceof UniformIntDist) {
        UniformIntDist distribution1 = (UniformIntDist) distributionGeneral[0][t];
        UniformIntDist distribution2 = (UniformIntDist) distributionGeneral[1][t];
        int demandLength1 = distribution1.getJ() - distribution1.getI() + 1;
        int demandLength2 = distribution2.getJ() - distribution2.getI() + 1;
        double[][] pmf = new double[demandLength1 * demandLength2][3];
        int index = 0;
        for (int i = distribution1.getXinf(); i <= distribution1.getXsup(); i++) for (int j = distribution2.getXinf(); j <= distribution2.getXsup(); j++) {
            pmf[index][0] = i;
            pmf[index][1] = j;
            pmf[index][2] = distribution1.prob(i) * distribution2.prob(j);
            index++;
        }
        return pmf;
    }
    return null;
}

Example 2 with UniformIntDist

use of umontreal.ssj.probdist.UniformIntDist in project Stochastic-Inventory by RobinChen121.

the class GetPmf method getpmf.

/**
 * @return
 */
public double[][][] getpmf() {
    int T = distributions.length;
    double[] supportLB = new double[T];
    double[] supportUB = new double[T];
    for (int i = 0; i < T; i++) {
        supportLB[i] = (int) distributions[i].inverseF(1 - truncationQuantile);
        if (distributions[0] instanceof DiscreteDistributionInt)
            supportLB[i] = 0;
        supportUB[i] = (int) distributions[i].inverseF(truncationQuantile);
    }
    double[][][] pmf = new double[T][][];
    if (distributions[0] instanceof UniformIntDist) {
        for (int i = 0; i < T; i++) {
            UniformIntDist distribution = (UniformIntDist) distributions[0];
            int demandLength = distribution.getJ() - distribution.getI() + 1;
            pmf[i] = new double[demandLength][];
            int index = 0;
            for (int j = distribution.getXinf(); j <= distribution.getXsup(); j++) {
                pmf[i][index] = new double[2];
                pmf[i][index][0] = j;
                pmf[i][index][1] = distribution.prob(j);
                index++;
            }
        }
        return pmf;
    }
    for (int i = 0; i < T; i++) {
        int demandLength = (int) ((supportUB[i] - supportLB[i] + 1) / stepSize);
        pmf[i] = new double[demandLength][];
        // demand values are all integers
        for (int j = 0; j < demandLength; j++) {
            pmf[i][j] = new double[2];
            pmf[i][j][0] = supportLB[i] + j * stepSize;
            if (distributions[0] instanceof DiscreteDistributionInt || distributions[0] instanceof PoissonDist) {
                // may be something wrong, poisson[] can't be (casted) delivered
                // but the results are correct
                double probilitySum = distributions[i].cdf(supportUB[i]) - distributions[i].cdf(supportLB[i] - 1);
                pmf[i][j][1] = ((DiscreteDistributionInt) distributions[i]).prob(j) / probilitySum;
            } else {
                double probilitySum = distributions[i].cdf(supportUB[i] + 0.5 * stepSize) - distributions[i].cdf(supportLB[i] - 0.5 * stepSize);
                pmf[i][j][1] = (distributions[i].cdf(pmf[i][j][0] + 0.5 * stepSize) - distributions[i].cdf(pmf[i][j][0] - 0.5 * stepSize)) / probilitySum;
            }
        }
    }
    return pmf;
}

Example 3 with UniformIntDist

use of umontreal.ssj.probdist.UniformIntDist in project Stochastic-Inventory by RobinChen121.

the class MultiItemCashXW method main.

public static void main(String[] args) {
    double[] price = { 2, 10 };
    // higher margin vs lower margin
    double[] variCost = { 1, 1 };
    double depositeRate = 0;
    // initial cash
    double iniCash = 10;
    // initial inventory
    int iniInventory1 = 0;
    int iniInventory2 = 0;
    // gamma distribution:mean demand is shape / beta and variance is shape / beta^2
    // beta = 1 / scale
    // shape = demand * beta
    // variance = demand / beta
    // gamma in ssj: alpha is alpha, and lambda is beta(beta)
    // horizon length
    int T = 4;
    double[] meanDemands = new double[] { 10, 3 };
    // higher average demand vs lower average demand
    double[][] demand = new double[2][T];
    // higher variance vs lower variance
    double[] beta = { 10, 1 };
    double d1 = meanDemands[0];
    double d2 = meanDemands[1];
    double v1 = variCost[0];
    double v2 = variCost[1];
    double p1 = price[0];
    double p2 = price[1];
    for (int t = 0; t < T; t++) {
        demand[0][t] = d1;
        demand[1][t] = d2;
    }
    double[] salPrice = Arrays.stream(variCost).map(a -> a * 0.5).toArray();
    // number of products
    int m = demand.length;
    // may affect poisson results
    double truncationQuantile = 0.9999;
    double stepSize = 1;
    double minCashState = 0;
    double maxCashState = 10000;
    int minInventoryState = 0;
    int maxInventoryState = 200;
    int Qbound = 20;
    double discountFactor = 1;
    // get demand possibilities for each period
    Distribution[][] distributions = new GammaDist[m][T];
    // Distribution[][] distributions =  new UniformIntDist[m][T];
    for (int i = 0; i < m; i++) for (int t = 0; t < T; t++) {
        distributions[i][t] = new GammaDist(demand[i][t] * beta[i], beta[i]);
    // distributions[i][t] = new UniformIntDist((int)(demand[i][t] * 0.6), (int)(demand[i][t] * 1.4));
    // distributions[i][t] = new PoissonDist(demand[i][t]);
    // distributions[i][t]= new NormalDist(demand[i][t], 0.1 * demand[i][t]);
    }
    GetPmfMulti PmfMulti = new GetPmfMulti(distributions, truncationQuantile, stepSize);
    // build action list (y1, y2) for pai(x1, x2, w)
    // CashStateMulti are states (x1, x2, w)
    Function<CashStateMulti, ArrayList<double[]>> buildActionListPai = s -> {
        ArrayList<double[]> actions = new ArrayList<>();
        double Ybound = Qbound;
        for (double i = 0; i < Ybound; i = i + 1) for (double j = 0; j < Ybound; j = j + 1) {
            double[] thisActions = { i, j };
            actions.add(thisActions);
        }
        return actions;
    };
    // build action list (y1, y2) for V(x1, x2, w)
    Function<CashStateMulti, ArrayList<double[]>> buildActionListV = s -> {
        ArrayList<double[]> actions = new ArrayList<>();
        int miny1 = (int) s.getIniInventory1();
        int miny2 = (int) s.getIniInventory2();
        double iniR = s.getIniCash() + v1 * s.getIniInventory1() + v2 * s.getIniInventory2();
        for (double i = miny1; i < miny1 + Qbound; i = i + 1) for (double j = miny2; j < miny2 + Qbound; j = j + 1) {
            if (v1 * i + v2 * j < iniR + 0.1) {
                double[] thisActions = { i, j };
                actions.add(thisActions);
            }
        }
        return actions;
    };
    BoundaryFuncton<CashStateMulti, Double> boundFinalCash = (IniState) -> {
        return IniState.getIniCash() + salPrice[0] * IniState.getIniInventory1() + salPrice[1] * IniState.getIniInventory2();
    };
    // State Transition Function
    // revise
    StateTransitionFunction<CashStateMulti, double[], double[], CashStateMulti> stateTransition = (IniState, actions, RandomDemands) -> {
        double x1 = IniState.getIniInventory1();
        double x2 = IniState.getIniInventory2();
        double y1 = actions[0];
        double y2 = actions[1];
        double endInventory1 = y1 - RandomDemands[0];
        endInventory1 = Math.max(0, endInventory1);
        double endInventory2 = y2 - RandomDemands[1];
        endInventory2 = Math.max(0, endInventory2);
        double revenue1 = p1 * Math.min(y1, RandomDemands[0]);
        double revenue2 = p2 * Math.min(y2, RandomDemands[1]);
        double nextW = revenue1 + revenue2 + (1 + depositeRate) * (IniState.getIniCash() - v1 * (y1 - x1) - // revise
        v2 * (y2 - x2));
        // rounding the states to one decimal 10.0
        endInventory1 = Math.round(endInventory1 * 10) / 10;
        endInventory2 = Math.round(endInventory2 * 10) / 10;
        nextW = Math.round(nextW * 10) / 10;
        nextW = nextW > maxCashState ? maxCashState : nextW;
        nextW = nextW < minCashState ? minCashState : nextW;
        endInventory1 = endInventory1 > maxInventoryState ? maxInventoryState : endInventory1;
        endInventory2 = endInventory2 < minInventoryState ? minInventoryState : endInventory2;
        return new CashStateMulti(IniState.getPeriod() + 1, endInventory1, endInventory2, nextW);
    };
    /**
     *****************************************************************
     * Solve
     */
    CashRecursionV2 recursion = new CashRecursionV2(discountFactor, PmfMulti, buildActionListV, buildActionListPai, stateTransition, boundFinalCash, T, variCost);
    int period = 1;
    CashStateMulti iniState = new CashStateMulti(period, iniInventory1, iniInventory2, iniCash);
    long currTime = System.currentTimeMillis();
    double finalValue = recursion.getExpectedValueV(iniState);
    System.out.println("final optimal cash  is " + finalValue);
    System.out.println("optimal order quantity in the first priod is :  y1 = " + recursion.getAction(iniState)[0] + ", y2 = " + recursion.getAction(iniState)[1]);
    double time = (System.currentTimeMillis() - currTime) / 1000.0;
    System.out.println("running time is " + time + "s");
    double[] optY = recursion.getYStar(iniState);
    System.out.println("optimal order quantity y* in the first priod is : " + Arrays.toString(optY));
    /**
     *****************************************************************
     * Simulate
     *
     * basically, this simulation is testing for Theorem 1:
     * optimal ordering decisions depend on y*(R)
     */
    int sampleNum = 10000;
    currTime = System.currentTimeMillis();
    CashSimulationY simulation = new CashSimulationY(sampleNum, distributions, discountFactor, recursion, stateTransition);
    double simFinalValue = simulation.simulateSDPGivenSamplNum2(iniState, variCost);
    double gap = (finalValue - simFinalValue) / finalValue;
    System.out.printf("optimality gap for this policy y* is %.2f%%\n", gap * 100);
    time = (System.currentTimeMillis() - currTime) / 1000.0;
    System.out.println("running time is " + time + "s");
    /**
     *****************************************************************
     * Compute a1* and a2*
     *
     * and simulate their results to test Theorem 2
     */
    stepSize = 0.1;
    double[][][] pmf1 = new GetPmf(distributions[0], truncationQuantile, stepSize).getpmf();
    Distribution[] distributions1 = distributions[0];
    double[][][] pmf2 = new GetPmf(distributions[1], truncationQuantile, stepSize).getpmf();
    Distribution[] distributions2 = distributions[1];
    RecursionG recursionG1 = new RecursionG(pmf1, distributions1, price[0], variCost[0], 0, salPrice[0]);
    RecursionG recursionG2 = new RecursionG(pmf2, distributions2, price[1], variCost[1], 0, salPrice[1]);
    double[] opta1 = recursionG1.getOptY();
    double[] opta2 = recursionG2.getOptY();
    System.out.println("---------------");
    System.out.println("a1* in each period:");
    DecimalFormat df = new DecimalFormat("0.00");
    Arrays.stream(opta1).forEach(e -> System.out.print(df.format(e) + " "));
    System.out.println("");
    System.out.println("a2* in each period:");
    Arrays.stream(opta2).forEach(e -> System.out.print(df.format(e) + " "));
    // double simFinalValue2 = simulation.simulateSDPGivenSamplNuma1a22(iniState, variCost, opta1, opta2);
    // double gap2 = (simFinalValue2 - finalValue) / finalValue;
    // System.out.printf("optimality gap for this policy a* is %.2f%%\n", gap2 * 100);
    double[] mean = new double[] { demand[0][0], demand[1][0] };
    double[] variance = new double[] { demand[0][0] / beta[0], demand[1][0] / beta[1] };
    double[][] optTable = recursion.getOptTableDetail(mean, variance, price, opta1, opta2);
    double[] gaps = new double[] { gap };
    WriteToExcel wr = new WriteToExcel();
    String fileName = "run" + ".xls";
    String headString = "meanD1" + "\t" + "meanD2" + "\t" + "variance1" + "\t" + "variance2" + "\t" + "period" + "\t" + "x1" + "\t" + "x2" + "\t" + "w" + "\t" + "p1" + "\t" + "p2" + "\t" + "c1" + "\t" + "c2" + "\t" + "R" + "\t" + "y1*" + "\t" + "y2*" + "\t" + "cashSituation" + "\t" + "alpha" + "\t" + "yHead1" + "\t" + "yHead2" + "\t" + "a1*" + "\t" + "a2*" + "\t" + "Theorem1Gap";
    wr.writeArrayToExcel2(optTable, fileName, headString, gaps);
}

Example 4 with UniformIntDist

use of umontreal.ssj.probdist.UniformIntDist in project Stochastic-Inventory by RobinChen121.

the class MultiItemCashXWTesting method main.

public static void main(String[] args) {
    double[] price = { 2, 10 };
    // higher margin vs lower margin
    double[] variCost = { 1, 2 };
    double depositeRate = 0;
    // initial cash
    double iniCash = 10;
    // initial inventory
    int iniInventory1 = 0;
    int iniInventory2 = 0;
    // gamma distribution:mean demand is shape / beta and variance is shape / beta^2
    // beta = 1 / scale
    // shape = demand * beta
    // variance = demand / beta
    // gamma in ssj: alpha is alpha, and lambda is beta(beta)
    // horizon length
    int T = 4;
    double[] meanDemands = new double[] { 10, 3 };
    // read parameter settings from excel files
    ReadExcel re = new ReadExcel();
    double[][] paraSettings = re.readExcelXLSX("Numerical experiments-2021-02-06.xlsx", 2);
    for (int runTime = 1; runTime < 2; runTime = runTime + 1) {
        price = new double[] { paraSettings[runTime][2], paraSettings[runTime][8] };
        variCost = new double[] { paraSettings[runTime][1], paraSettings[runTime][7] };
        double[] beta = new double[] { paraSettings[runTime][6], paraSettings[runTime][12] };
        // a = new int[] {(int)paraSettings[runTime][5], (int)paraSettings[runTime][11]}; // integer for uniform
        meanDemands = new double[] { paraSettings[runTime][3], paraSettings[runTime][9] };
        // higher average demand vs lower average demand
        double[][] demand = new double[2][T];
        double d1 = meanDemands[0];
        double d2 = meanDemands[1];
        double v1 = variCost[0];
        double v2 = variCost[1];
        double p1 = price[0];
        double p2 = price[1];
        for (int t = 0; t < T; t++) {
            demand[0][t] = d1;
            demand[1][t] = d2;
        }
        double[] salPrice = Arrays.stream(variCost).map(a -> a * 0.5).toArray();
        // number of products
        int m = demand.length;
        // may affect poisson results
        double truncationQuantile = 0.9999;
        double stepSize = 1;
        double minCashState = 0;
        double maxCashState = 10000;
        int minInventoryState = 0;
        int maxInventoryState = 200;
        int Qbound = 20;
        double discountFactor = 1;
        // get demand possibilities for each period
        // Distribution[][] distributions =  new GammaDist[m][T];
        // Distribution[][] distributions =  new PoissonDist[m][T];
        // Distribution[][] distributions =  new NormalDist[m][T];
        Distribution[][] distributions = new UniformIntDist[m][T];
        for (int i = 0; i < m; i++) for (int t = 0; t < T; t++) {
            // distributions[i][t] = new GammaDist(demand[i][t]* beta[i], beta[i]);
            distributions[i][t] = new UniformIntDist((int) (demand[i][t] * 0.2), (int) (demand[i][t] * 1.8));
        // distributions[i][t] = new PoissonDist(demand[i][t]);
        // distributions[i][t]= new NormalDist(demand[i][t], 0.1 * demand[i][t]);
        // distributions[i][t] = new UniformIntDist(a[i], b[i]);
        }
        GetPmfMulti PmfMulti = new GetPmfMulti(distributions, truncationQuantile, stepSize);
        // build action list (y1, y2) for pai(x1, x2, w)
        // CashStateMulti are states (x1, x2, w)
        Function<CashStateMulti, ArrayList<double[]>> buildActionListPai = s -> {
            ArrayList<double[]> actions = new ArrayList<>();
            double Ybound = Qbound;
            for (double i = 0; i < Ybound; i = i + 1) for (double j = 0; j < Ybound; j = j + 1) {
                double[] thisActions = { i, j };
                actions.add(thisActions);
            }
            return actions;
        };
        // build action list (y1, y2) for V(x1, x2, w), no use in the recursion
        Function<CashStateMulti, ArrayList<double[]>> buildActionListV = s -> {
            ArrayList<double[]> actions = new ArrayList<>();
            int miny1 = (int) s.getIniInventory1();
            int miny2 = (int) s.getIniInventory2();
            double iniR = s.getIniCash() + v1 * s.getIniInventory1() + v2 * s.getIniInventory2();
            for (double i = miny1; i < miny1 + Qbound; i = i + 1) for (double j = miny2; j < miny2 + Qbound; j = j + 1) {
                if (v1 * i + v2 * j < iniR + 0.1) {
                    double[] thisActions = { i, j };
                    actions.add(thisActions);
                }
            }
            return actions;
        };
        BoundaryFuncton<CashStateMulti, Double> boundFinalCash = (IniState) -> {
            return IniState.getIniCash() + salPrice[0] * IniState.getIniInventory1() + salPrice[1] * IniState.getIniInventory2();
        };
        // State Transition Function
        // revise
        StateTransitionFunction<CashStateMulti, double[], double[], CashStateMulti> stateTransition = (IniState, actions, RandomDemands) -> {
            double x1 = IniState.getIniInventory1();
            double x2 = IniState.getIniInventory2();
            double y1 = actions[0];
            double y2 = actions[1];
            double endInventory1 = y1 - RandomDemands[0];
            endInventory1 = Math.max(0, endInventory1);
            double endInventory2 = y2 - RandomDemands[1];
            endInventory2 = Math.max(0, endInventory2);
            double revenue1 = p1 * Math.min(y1, RandomDemands[0]);
            double revenue2 = p2 * Math.min(y2, RandomDemands[1]);
            double nextW = revenue1 + revenue2 + (1 + depositeRate) * (IniState.getIniCash() - v1 * (y1 - x1) - // revise
            v2 * (y2 - x2));
            // rounding the states to one decimal 10.0
            endInventory1 = Math.round(endInventory1 * 1) / 1;
            // very slow when decimal
            endInventory2 = Math.round(endInventory2 * 1) / 1;
            nextW = Math.round(nextW * 1) / 1;
            nextW = nextW > maxCashState ? maxCashState : nextW;
            nextW = nextW < minCashState ? minCashState : nextW;
            endInventory1 = endInventory1 > maxInventoryState ? maxInventoryState : endInventory1;
            endInventory2 = endInventory2 < minInventoryState ? minInventoryState : endInventory2;
            return new CashStateMulti(IniState.getPeriod() + 1, endInventory1, endInventory2, nextW);
        };
        /**
         *****************************************************************
         * Solve
         */
        CashRecursionV2 recursion = new CashRecursionV2(discountFactor, PmfMulti, buildActionListV, buildActionListPai, stateTransition, boundFinalCash, T, variCost);
        int period = 1;
        CashStateMulti iniState = new CashStateMulti(period, iniInventory1, iniInventory2, iniCash);
        long currTime = System.currentTimeMillis();
        double finalValue = recursion.getExpectedValueV(iniState);
        System.out.println("final optimal cash  is " + finalValue);
        System.out.println("optimal order quantity in the first priod is :  y1 = " + recursion.getAction(iniState)[0] + ", y2 = " + recursion.getAction(iniState)[1]);
        double time = (System.currentTimeMillis() - currTime) / 1000.0;
        System.out.println("running time is " + time + "s");
        double[] optY = recursion.getYStar(iniState);
        System.out.println("optimal order quantity y* in the first priod is : " + Arrays.toString(optY));
        /**
         *****************************************************************
         * Simulate
         *
         * basically, this simulation is testing for Theorem 1:
         * optimal ordering decisions depend on y*(R)
         */
        int sampleNum = 10000;
        currTime = System.currentTimeMillis();
        CashSimulationY simulation = new CashSimulationY(sampleNum, distributions, discountFactor, recursion, stateTransition);
        double simFinalValue = simulation.simulateSDPGivenSamplNum2(iniState, variCost);
        double gap = (finalValue - simFinalValue) / finalValue;
        System.out.printf("optimality gap for this policy y* is %.2f%%\n", gap * 100);
        time = (System.currentTimeMillis() - currTime) / 1000.0;
        System.out.println("running time is " + time + "s");
        /**
         *****************************************************************
         * Compute a1* and a2*
         *
         * and simulate their results to test Theorem 2
         */
        // can be changed
        stepSize = 1;
        double[][][] pmf1 = new GetPmf(distributions[0], truncationQuantile, stepSize).getpmf();
        Distribution[] distributions1 = distributions[0];
        double[][][] pmf2 = new GetPmf(distributions[1], truncationQuantile, stepSize).getpmf();
        Distribution[] distributions2 = distributions[1];
        RecursionG recursionG1 = new RecursionG(pmf1, distributions1, price[0], variCost[0], 0, salPrice[0]);
        RecursionG recursionG2 = new RecursionG(pmf2, distributions2, price[1], variCost[1], 0, salPrice[1]);
        double[] opta1 = recursionG1.getOptY();
        double[] opta2 = recursionG2.getOptY();
        System.out.println("---------------");
        System.out.println("a1* in each period:");
        DecimalFormat df = new DecimalFormat("0.00");
        Arrays.stream(opta1).forEach(e -> System.out.print(df.format(e) + " "));
        System.out.println("");
        System.out.println("a2* in each period:");
        Arrays.stream(opta2).forEach(e -> System.out.print(df.format(e) + " "));
        // double simFinalValue2 = simulation.simulateSDPGivenSamplNuma1a22(iniState, variCost, opta1, opta2);
        // double gap2 = (simFinalValue2 - finalValue) / finalValue;
        // System.out.printf("optimality gap for this policy a* is %.2f%%\n", gap2 * 100);
        System.out.println();
        System.out.println("*****************************************");
        double[] mean = new double[] { demand[0][0], demand[1][0] };
        double[] variance = new double[] { demand[0][0] / beta[0], demand[1][0] / beta[1] };
        double[][] optTable = recursion.getOptTableDetail(mean, variance, price, opta1, opta2);
        double[] gaps = new double[] { gap };
        WriteToExcel wr = new WriteToExcel();
        String fileName = "run" + (int) paraSettings[runTime][0] + ".xls";
        String headString = "meanD1" + "\t" + "meanD2" + "\t" + "variance1" + "\t" + "variance2" + "\t" + "period" + "\t" + "x1" + "\t" + "x2" + "\t" + "w" + "\t" + "p1" + "\t" + "p2" + "\t" + "c1" + "\t" + "c2" + "\t" + "R" + "\t" + "y1*" + "\t" + "y2*" + "\t" + "cashSituation" + "\t" + "alpha" + "\t" + "yHead1" + "\t" + "yHead2" + "\t" + "a1*" + "\t" + "a2*" + "\t" + "Theorem1Gap";
        wr.writeArrayToExcel2(optTable, fileName, headString, gaps);
    }
}

Example 5 with UniformIntDist

use of umontreal.ssj.probdist.UniformIntDist in project Stochastic-Inventory by RobinChen121.

the class MultiItemYRTesting method main.

public static void main(String[] args) {
    double[] price = { 2, 10 };
    // higher margin vs lower margin
    double[] variCost = { 1, 2 };
    double depositebeta = 0;
    // initial cash
    double iniCash = 10;
    // initial inventory
    int iniInventory1 = 0;
    int iniInventory2 = 0;
    // gamma distribution:mean demand is shape / beta and variance is shape / beta^2
    // beta = 1 / scale
    // shape = demand * beta
    // variance = demand / beta
    // gamma in ssj: alpha is alpha, and lambda is beta(beta)
    // horizon length
    int T = 4;
    double[] meanDemands = new double[2];
    // higher average demand vs lower average demand
    double[][] demand = new double[2][T];
    // higher variance vs lower variance
    double[] beta = { 10, 1 };
    int[] b = new int[2];
    int[] a = new int[2];
    // read parameter settings from excel files
    ReadExcel re = new ReadExcel();
    double[][] paraSettings = re.readExcelXLSX("Numerical experiments-settings.xlsx", 2);
    for (int runTime = 3; runTime < 16; runTime = runTime + 2) {
        price = new double[] { paraSettings[runTime][2], paraSettings[runTime][8] };
        variCost = new double[] { paraSettings[runTime][1], paraSettings[runTime][7] };
        b = new int[] { (int) paraSettings[runTime][6], (int) paraSettings[runTime][12] };
        a = new int[] { (int) paraSettings[runTime][5], (int) paraSettings[runTime][11] };
        meanDemands = new double[] { paraSettings[runTime][3], paraSettings[runTime][9] };
        double d1 = meanDemands[0];
        double d2 = meanDemands[1];
        double v1 = variCost[0];
        double v2 = variCost[1];
        double p1 = price[0];
        double p2 = price[1];
        for (int t = 0; t < T; t++) {
            demand[0][t] = d1;
            demand[1][t] = d2;
        }
        // unit salvage value is half of the unit variable cost
        double[] salPrice = Arrays.stream(variCost).map(s -> s * 0.5).toArray();
        // number of products
        int m = demand.length;
        // may affect poisson results
        double truncationQuantile = 0.9999;
        int stepSize = 1;
        double minCashState = 0;
        double maxCashState = 10000;
        int minInventoryState = 0;
        int maxInventoryState = 200;
        int Qbound1 = 20;
        int Qbound2 = 10;
        double discountFactor = 1;
        // get demand possibilities for each period
        Distribution[][] distributions = new GammaDist[m][T];
        // Distribution[][] distributions =  new NormalDist[m][T];
        for (int i = 0; i < m; i++) for (int t = 0; t < T; t++) {
            // distributions[i][t] = new UniformIntDist(a[i], b[i]);
            distributions[i][t] = new GammaDist(demand[i][t] * beta[i], beta[i]);
        // distributions[i][t] = new PoissonDist(demand[i][t]);
        // distributions[i][t]= new NormalDist(demand[i][t], 0.1 * demand[i][t]);
        }
        // build action list (y1, y2) for V(x1, x2, R)
        Function<CashStateMultiYR, ArrayList<double[]>> buildActionListPai = s -> {
            ArrayList<double[]> actions = new ArrayList<>();
            for (double i = 0; i < Qbound1; i++) for (double j = 0; j < Qbound2; j++) {
                double[] thisActions = { i, j };
                actions.add(thisActions);
            }
            return actions;
        };
        // build action list (y1, y2) for Pai(x1, x2, R)
        Function<CashStateMulti, ArrayList<double[]>> buildActionListV = s -> {
            ArrayList<double[]> actions = new ArrayList<>();
            int miny1 = (int) s.getIniInventory1();
            int miny2 = (int) s.getIniInventory2();
            double iniR = s.getIniCash() + v1 * s.getIniInventory1() + v2 * s.getIniInventory2();
            for (int i = miny1; i < miny1 + Qbound1; i++) for (int j = miny2; j < miny2 + Qbound2; j++) {
                if (v1 * i + v2 * j < iniR + 0.1) {
                    double[] thisActions = { i, j };
                    actions.add(thisActions);
                }
            }
            return actions;
        };
        BoundaryFuncton<CashStateMulti, Double> boundFinalCash = (IniState) -> {
            return IniState.getIniCash() + salPrice[0] * IniState.getIniInventory1() + salPrice[1] * IniState.getIniInventory2();
        };
        // State Transition Function
        StateTransitionFunctionV<CashStateMultiYR, double[], CashStateMulti> stateTransition = (IniState, RandomDemands) -> {
            double endInventory1 = IniState.getIniInventory1() - RandomDemands[0];
            endInventory1 = Math.max(0, endInventory1);
            double endInventory2 = IniState.getIniInventory2() - RandomDemands[1];
            endInventory2 = Math.max(0, endInventory2);
            double revenue1 = p1 * Math.min(IniState.getIniInventory1(), RandomDemands[0]);
            double revenue2 = p2 * Math.min(IniState.getIniInventory2(), RandomDemands[1]);
            double nextW = revenue1 + revenue2 + (1 + depositebeta) * (IniState.getIniR() - v1 * IniState.getIniInventory1() - // revise
            v2 * IniState.getIniInventory2());
            endInventory1 = Math.round(endInventory1);
            endInventory2 = Math.round(endInventory2);
            nextW = Math.round(nextW);
            nextW = nextW > maxCashState ? maxCashState : nextW;
            nextW = nextW < minCashState ? minCashState : nextW;
            endInventory1 = endInventory1 > maxInventoryState ? maxInventoryState : endInventory1;
            endInventory2 = endInventory2 < minInventoryState ? minInventoryState : endInventory2;
            return new CashStateMulti(IniState.getPeriod() + 1, endInventory1, endInventory2, nextW);
        };
        GetPmfMulti PmfMulti = new GetPmfMulti(distributions, truncationQuantile, stepSize);
        /**
         *****************************************************************
         * Solve
         */
        CashRecursionV recursion = new CashRecursionV(discountFactor, PmfMulti, buildActionListV, buildActionListPai, stateTransition, boundFinalCash, T, variCost);
        int period = 1;
        CashStateMulti iniState = new CashStateMulti(period, iniInventory1, iniInventory2, iniCash);
        long currTime = System.currentTimeMillis();
        double finalValue = recursion.getExpectedValueV(iniState);
        System.out.println("final optimal cash  is " + finalValue);
        System.out.println("optimal order quantity in the first priod is :  y1 = " + recursion.getAction(iniState)[0] + ", y2 = " + recursion.getAction(iniState)[1]);
        double time = (System.currentTimeMillis() - currTime) / 1000.0;
        System.out.println("running time is " + time + "s");
        CashStateR iniState2 = new CashStateR(period, iniCash);
        double[] optY = recursion.getYStar(iniState2);
        System.out.println("optimal order quantity y* in the first priod is : " + Arrays.toString(optY));
        double[] mean = new double[] { demand[0][0], demand[1][0] };
        double[] variance = new double[] { demand[0][0] / beta[0], demand[1][0] / beta[1] };
        /**
         *****************************************************************
         * Simulate
         *
         * basically, this simulation is testing for Theorem 1:
         * optimal ordering decisions depend on y*(R)
         */
        int sampleNum = 10000;
        currTime = System.currentTimeMillis();
        CashSimulationY simulation = new CashSimulationY(sampleNum, distributions, discountFactor, recursion, stateTransition);
        double simFinalValue = simulation.simulateSDPGivenSamplNum(iniState, variCost);
        double gap = (simFinalValue - finalValue) / finalValue;
        System.out.printf("optimality gap for this policy y* is %.2f%%\n", gap * 100);
        time = (System.currentTimeMillis() - currTime) / 1000.0;
        System.out.println("running time is " + time + "s");
        /**
         *****************************************************************
         * Compute a1* and a2*
         *
         * and simulate their results to test Theorem 2
         */
        double[][][] pmf1 = new GetPmf(distributions[0], truncationQuantile, stepSize).getpmf();
        Distribution[] distributions1 = distributions[0];
        double[][][] pmf2 = new GetPmf(distributions[1], truncationQuantile, stepSize).getpmf();
        Distribution[] distributions2 = distributions[1];
        RecursionG recursionG1 = new RecursionG(pmf1, distributions1, price[0], variCost[0], 0, salPrice[0]);
        RecursionG recursionG2 = new RecursionG(pmf2, distributions2, price[1], variCost[1], 0, salPrice[1]);
        double[] opta1 = recursionG1.getOptY();
        double[] opta2 = recursionG2.getOptY();
        System.out.println("a1* in each period:");
        DecimalFormat df = new DecimalFormat("0.00");
        Arrays.stream(opta1).forEach(e -> System.out.print(df.format(e) + " "));
        System.out.println("");
        System.out.println("a2* in each period:");
        Arrays.stream(opta2).forEach(e -> System.out.print(df.format(e) + " "));
        double simFinalValue2 = simulation.simulateSDPGivenSamplNuma1a2(iniState, variCost, opta1, opta2);
        double gap2 = (simFinalValue2 - finalValue) / finalValue;
        System.out.printf("optimality gap for this policy a* is %.2f%%\n", gap2 * 100);
        double[][] optTable = recursion.getOptTableDetail2(mean, variance, price, opta1, opta2);
        double[] gaps = new double[] { gap, gap2 };
        WriteToExcel wr = new WriteToExcel();
        String fileName = "run" + (int) paraSettings[runTime][0] + ".xls";
        String headString = "meanD1" + "\t" + "meanD2" + "\t" + "variance1" + "\t" + "variance2" + "\t" + "period" + "\t" + "x1" + "\t" + "x2" + "\t" + "w" + "\t" + "p1" + "\t" + "p2" + "\t" + "c1" + "\t" + "c2" + "\t" + "R" + "\t" + "y1*" + "\t" + "y2*" + "\t" + "cashSituation" + "\t" + "alpha" + "\t" + "yHead1" + "\t" + "yHead2" + "\t" + "a1*" + "\t" + "a2*" + "\t" + "Theorem1Gap" + "Theorem2Gap";
        wr.writeArrayToExcel2(optTable, fileName, headString, gaps);
    // System.out.println("alpha in the first period: " + optTable[0][10]);
    // System.out.println("*******************************");
    }
}