A Graph Theoretic Approach to Non-Anticipativity Constraint Generation in Multistage Stochastic Programs with Incomplete Scenario Sets

TL;DR

This paper introduces a polynomially scalable algorithm for generating minimal non-anticipativity constraints in multistage stochastic programming with incomplete scenario sets, accommodating gradual uncertainty realizations.

Contribution

The paper presents SNAC, a novel algorithm that relaxes previous limitations, enabling efficient NAC generation for complex scenario structures in MSSPs.

Findings

SNAC scales polynomially with the number of scenarios.

It handles both endogenous and exogenous uncertainties.

The approach relaxes the need for full scenario sets.

Abstract

We propose an algorithm for generating a minimum-cardinality set of non-anticipativity constraints (NAC) for scenario-based multistage-stochastic programming (MSSP) problems with both endogenous and exogenous uncertainties which allow for gradual realizations. Recently several authors have considered approaches to generate the minimum cardinality NAC set for MSSPs for various scenario set structures. However, these approaches have been limited to uncertain parameters where the realizations occur instantaneously or the full set of scenarios is required. The proposed algorithm, referred to as Sample Non-Anticipativity Constraint algorithm (SNAC) relaxes this requirement. We show that as long as the number of uncertain parameters and parameter values are kept constant, the algorithm scales polynomially in the number of scenarios.

Tables

Table 1. Table 1: The event sets 𝒞 𝒞 \mathcal{C} for the manufacturing example.

Order	$𝒞$
0	{}
1	$\left\{e_{P1}^1\right\},\left\{e_{P2}^1\right\}$
2	$\{e_{P2}^1,e_{P2}^2\},\{e_{P1}^1,e_{P2}^1\},\{e_{P1}^1,e_{P1}^2\}$
3	$\{e_{P2}^1,e_{P2}^2,e_{P2}^3\},\{e_{P1}^1,e_{P2}^1,e_{P2}^2\},\{e_{P1}^1,e_{P1}^2,e_{P2}^1\},\{e_{P1}^1,e_{P1}^2,e_{P1}^3\}$
4	$\{e_{P1}^1,e_{P1}^2,e_{P2}^1,e_{P2}^2\},\{e_{P1}^1,e_{P1}^2,e_{P1}^3,e_{P2}^1\}\{e_{P1}^1,e_{P2}^1,e_{P2}^2,e_{P2}^3\}$
5	$\{e_{P1}^1,e_{P1}^2,e_{P2}^1,e_{P2}^2,e_{P2}^3\},\{e_{P1}^1,e_{P1}^2,e_{P1}^3,e_{P2}^1,e_{P2}^2\}$
6	$\{e_{P1}^1,e_{P1}^2,e_{P1}^3,e_{P2}^1,e_{P2}^2,e_{P2}^3\}$

Table 2. Table 2: Summary of problem sizes for case studies solving sample MSSPs based on the manufacturing example. All values are averages or maximums over 30 instances. Average and maximum for numbers of pairs of scenarios and actual NACs are reported for SNAC algorithm. Columns Full report the number of NACs and scenarios pairs in the full problem (does not change with each random instance).

			# of NACs			# of scenario pairs
Products	Outcomes	Scenarios	Average	Max	Full	Average	Max	Full
2	4	12	793.6	840	3192	16	17	66
2	10	24	4593.6	4944	39792	31.6	34	276
3	4	6	1351.2	1998	2988	6.7	10	15
3	5	24	10869.6	13008	72936	40.9	49	276
4	3	12	1730.7	2128	5328	21	26	66
4	4	128	39196	42032	975872	322.4	346	8128
4	5	24	9370.7	11136	44256	58	69	276
5	4	64	37350	42170	302720	246.9	279	2016
5	4	1024	581120	581120	78571520	3840	3840	523776

Table 3. Table 3: Summary of solution times for case studies solving sample MSSPs based on the manufacturing example. All times are reported in seconds and are averages over 30 instances. SNAC Running Time refers to the running time of the algorithm, Model Generation is the time required to generate the NACs (either all or based on minimum pairs), and CPLEX represents the time used by CPLEX to obtain the optimal solution (either full problem or the one produced by SNAC).

			SNAC Running	Model Generation		CPLEX
Products	Outcomes	Scenarios	Time	SNAC	Full	SNAC	Full
2	4	12	0.321	0.123	0.215	0.142	0.234
2	10	24	1.597	0.546	2.317	0.648	2.337
3	4	6	0.707	0.205	0.266	0.245	0.284
3	5	24	1.295	1.101	4.2	1.261	5.147
4	3	12	0.759	0.22	0.384	0.223	0.376
4	4	128	8.785	5.095	53.451	3.814	68.479
4	5	24	10.892	1.46	2.986	0.972	4.646
5	4	64	10.362	3.04	17.2	6.213	23.313
5	4	1024	5140.657	229.350	3545.935	43.481	4331.662

Table 4. Table 4: Parameters of the two-drug case study. Clinical trial plan for a 15-month planning horizon divided into five equal time periods.

Drug	Duration		Probability of Success		Cost ($M)		Resource 1 (Max = 2)		Resource 2 (Max = 3)		revmax	${\gamma }^L$	${\gamma }^D$
	PI	PII	PI	PII	PI	PII	PI	PII	PI	PII
D1	2	4	0.3	0.5	10	90	1	1	1	2	3100	19.2	44
D2	2	3	0.4	0.6	10	80	1	2	1	1	3250	19.6	56

Table 5. Table 5: Parameters of the three-drug case study. Clinical Trial plan for a 36-month planning horizon divided into 12 equal time periods.

Drug	Duration			Probability of Success			Trial Cost ($M)			Resource 1 (Max = 2)			Resource 1 (Max = 3)			revmax	${\gamma }^L$	${\gamma }^D$
	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII
D1	2	4	4	0.3	0.5	0.8	10	90	220	1	1	2	1	2	3	3100	19.2	44
D2	2	3	5	0.4	0.6	0.8	10	80	200	1	2	2	1	1	3	3250	19.6	56
D3	2	3	4	0.3	0.6	0.9	10	90	180	1	1	2	1	1	3	3300	20	52

Table 6. Table 6: Parameters of the four-drug case study. Clinical trial plan for an 18-month planning horizon divided into six equal time periods.

Drug	Duration			Probability of Success			Trial Cost ($M)			Resource 1 (Max = 4)			Resource 1 (Max = 3)			revmax	${\gamma }^L$	${\gamma }^D$
	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII
D1	1	1	3	0.3	0.5	0.8	10	90	220	1	1	2	1	2	3	3100	19.2	22
D2	1	2	2	0.4	0.6	0.8	10	80	200	1	2	2	1	1	3	3250	19.6	28
D3	1	1	3	0.3	0.6	0.9	10	90	180	1	1	2	1	1	3	3300	20	26
D4	1	2	2	0.4	0.6	0.8	10	100	170	1	1	2	1	2	3	3000	19.4	24

Table 7. Table 7: Parameters of the five-drug case study. Clinical trial plan for an 18-month planning horizon divided into six equal time periods.

Drug	Duration			Probability of Success			Trial Cost ($M)			Resource 1 (Max = 4)			Resource 1 (Max = 3)			revmax	${\gamma }^L$	${\gamma }^D$
	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII
D1	1	1	3	0.3	0.5	0.8	10	90	220	1	1	2	1	2	3	3100	19.2	22
D2	1	2	2	0.4	0.6	0.8	10	80	200	1	2	2	1	1	3	3250	19.6	28
D3	1	1	3	0.3	0.6	0.9	10	90	180	1	1	2	1	1	3	3300	20	26
D4	1	2	2	0.4	0.6	0.8	10	100	170	1	1	2	1	2	3	3000	19.4	24
D5	1	2	3	0.35	0.5	0.9	10	70	210	1	1	2	1	1	3	3150	19.6	24

Table 8. Table 8: Parameters for a six-drug case study. Clinical trial plan for an 18-month planning horizon divided into six equal time periods.

Drug	Duration			Probability of Success			Trial Cost ($M)			Resource 1 (Max = 4)			Resource 1 (Max = 3)			revmax	${\gamma }^L$	${\gamma }^D$
	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII	PI	PII	PIII
D1	1	1	3	0.3	0.5	0.8	10	90	220	1	1	2	1	2	3	3100	19.2	22
D2	1	2	2	0.4	0.6	0.8	10	80	200	1	2	2	1	1	3	3250	19.6	28
D3	1	1	3	0.3	0.6	0.9	10	90	180	1	1	2	1	1	3	3300	20	26
D4	1	2	2	0.4	0.6	0.8	10	100	170	1	1	2	1	2	3	3000	19.4	24
D5	1	2	3	0.35	0.5	0.9	10	70	210	1	1	2	1	1	3	3150	19.6	24
D6	1	2	3	0.45	0.45	0.8	10	85	195	1	2	2	2	1	3	3050	19	25

Equations

$$[Z_t^{rs}]\bigvee \left[\begin{array}{c}\neg Z_t^{rs}\\ {\mathbf{y}}_{tsc}={\mathbf{y}}_{trc}\\ {\mathbf{b}}_{tsc}={\mathbf{b}}_{trc}\end{array}\right]\forall r,s\in \mathbb{S},r\ne s,\forall t\in \mathcal{T},\forall c\in C$$

$$Z_t^{rs}=F(t,\mathbf{b},\mathbf{y}),\forall r,s\in \mathbb{S},\forall t\in \mathcal{T}$$

$$\left[\begin{array}{c}{\mathbf{y}}_{tsc}={\mathbf{y}}_{trc}\\ {\mathbf{b}}_{tsc}={\mathbf{b}}_{trc},\end{array}\right],\forall t\le \tau (r,s),\forall c\in C$$

$$[Z_t^{rs}]\bigvee \left[\begin{array}{c}\neg Z_t^{rs}\\ {\mathbf{y}}_{tsc}={\mathbf{y}}_{trc}\\ {\mathbf{b}}_{tsc}={\mathbf{b}}_{trc}\end{array}\right]\forall r,s\in \mathbb{S},r\ne s,\forall t\le \tau (r,s),\forall c\in C.$$

$$\Pi (c):S=\bigcup _j{S}_j^c,$$

$$\min$$

$$V_{i,j,t,s}=V_{i,j,t-1,s}+X_{i,j,t-{\tau }_{i,j},s},\forall i,j,t,s$$

$$Z_{i,1,1,s}=1-X_{i,1,t,s},\forall i,s$$

$$Z_{i,1,t,s}=Z_{i,1,t-1,s}-X_{i,j,t,s},\forall i,t>1,s$$

$$Z_{i,j,t,s}=Z_{i,j,t-1,s}+X_{i,j-1,t-{\tau }_{i,j-1},s}-X_{i,j,t,s},\forall i,j>1,t,s$$

$$\sum _t{X}_{i,j,t,s}\le 1,\forall i,j,s$$

$$\sum _{tt\le t}{X}_{i,j,tt,s}\le V_{i,j-1,t,s},\forall i,j>1,t,s$$

$$\sum _i\sum _j\sum _{tt>t-{\tau }_{i,j}}^{t\le t}{\rho }_{i,j,r}{X}_{i,j,tt,s}\le {{\rho }_r}^{max},\forall r,t,s$$

$$X_{i,1,1,s}=X_{i,1,1,1},\forall i,s$$

$${-V}_{i^{s,sp},j^{s,sp},t,s}\le X_{i,j,t,s}-X_{i,j,t,sp}\le V_{i^{s,sp},j^{s,sp},t,s}\forall i,j,\left(s,sp\right)\in B,t>1$$

$${Cst}_s=\sum _{i,j,t}{cd}_t{C}_{i,j}{X}_{i,j,t,s},\forall s$$

$${Rv}_s=\sum _i\sum _t{rev}_i^{max}{X}_{i,PII,t,s}-{{\gamma }_i}^D\left(Z_{i,PII,t,s}+Z_{i,PIII,t,s}\right)-{{\gamma }_i}^L(t+{\tau }_{i,PIII})X_{i,PIII,t,s},\forall s$$

$${FRev}_s=\sum _i\sum _j{rev}_{i,j}^{open}{f}_{i,j}{Z}_{i,j,|T|,s}+\sum _i\sum _{j\in \{PI,PII\}}\sum _{t>|T|-{\tau }_{i,j}}{rev}_{i,j,t}^{run}{f}_{i,j+1}{X}_{i,j,t,s},\forall s$$

$${{rev}_{i,j}}^{open}={{rev}_i}^{max}-{{\gamma }_i}^L\left(\left|T\right|+\sum _{jj\ge j}{\tau }_{i,jj}\right)$$

$${{rev}_{i,j}}^{run}={{rev}_i}^{max}-{{\gamma }_i}^L\left(t+\sum _{jj\ge j}{\tau }_{i,jj}\right)$$

$$\displaystyle f_{i,j}=0.9\bigg{[}\frac{rev_{i}^{max}-\gamma^{L}_{i}|T|-\sum_{j^{\prime}\geq j}C_{i,j}}{rev_{i}^{max}-\gamma_{i}^{L}|T|}\bigg{]}$$

$$X_{i,j,t,s},V_{i,j,t,s},Z_{i,j,t,s}\in \{0,1\}$$

Full text

1 Introduction

Stochastic programming is a scenario-based approach for optimization under uncertainty (Birge and Louveaux, 2011). Two-stage stochastic programs have a first stage where decisions are made and a second stage where uncertainty is realized and recourse actions can be taken. In contrast, multistage stochastic programs (MSSPs) have more than one stage where decisions are made and one or more periods where uncertainty is realized and recourse actions are taken. Note that the number of time periods in the planning horizon does not necessarily correspond to the number of stages in the problem.

The scenarios in a stochastic program are generated based on the outcomes of the uncertain parameters. A common approach for scenario generation in problems where uncertain parameters can be considered independent is to take the Cartesian product of the outcomes of the uncertain parameters, as described by Apap and Grossmann (2017). Uncertainty in stochastic programming can be classified based on how and when it is realized. Realization of exogenous, or decision-independent uncertain parameters occurs at a known time independent of the values of the decision variables (Birge and Louveaux, 2011). In operations research, the market conditions such as demand are often assumed to be independent of any investment decisions and are treated as exogenous uncertain parameters. With endogenous, or decision-dependent uncertainty, the timing of the realization or the value of the realized parameter is impacted by the value of the decision variables (Jonsbråten et al., 1998).

Another way to classify uncertainty is based on how the uncertainty is realized. Uncertainty in stochastic programs can either be realized instantaneously or gradually along the planning horizon (Solak, 2007). In the case of exogenous uncertainty, whether it is realized instantaneously or gradually has little to no effect on the multistage stochastic programming formulation, which is not the case for endogenous parameters. A well-known example of this type of uncertainty realization can be seen in the planning and scheduling of pharmaceutical clinical trials (Colvin and Maravelias, 2008). In the problem, the goal is to find the clinical trial schedule which provides the highest expected net present value given uncertainty in the outcome of clinical trials. The authors model the clinical trial outcome uncertainty for each potential new product using a single uncertain parameter. Uncertainty is only realized if a clinical trial fails or if all clinical trials are successfully completed. This leads to time periods in the planning horizon where a product may have successfully completed one or more clinical trials (i.e., no failures) but has not successfully completed all of them.

The formulation of an MSSP requires the equality of decision variables in different scenarios to ensure that the realized values of the uncertain parameters are not anticipated. Depending on the type of uncertainty in the problem, this can be accomplished either implicitly or explicitly. When the timing of the realization of uncertainty is known (i.e., exogenous uncertainty), anticipation of the underlying value of the uncertain parameters can be prevented by defining a variable used to relate the decision variables in two different scenarios. Using the same decision variable in multiple scenarios ensures that the decision variable values are identical in scenarios where uncertainty has not been realized. Alternatively, a set of equality constraints, called non-anticipativity constraints (NACs), which set the values of the decision variables equal to each other at time periods before uncertainty is realized can be used. If the time period of differentiation is not known (i.e., endogenous uncertainty), a set of NACs must be used as there is no way to implicitly incorporate the realization of uncertainty.

In this paper, we introduce and describe in detail a graph theoretic approach that generates a minimum cardinality NAC set for MSSPs with an incomplete scenario set where the realization of uncertainty occurs gradually.

While in principle it is possible to reformulate the problem as an instantaneous realization model by introducing additional uncertain parameters (and consequently apply some of the existing methods), our approach has the benefit of operating directly on the natural original formulation.

The remainder of the paper is organized as follows. In Section 1.1 we review existing relevant literature. Next, we present the notation and give the general NAC formulation for MSSPs with gradual realization of uncertainty in Section 2. In Section 3, we propose the algorithm, which generates a minimum cardinality NAC set. We also present an illustrative example for a small-scale problem (Section 3.2) and report results of a computational case study (Section 4). Conclusions and future directions are given in Section 5.

1.1 Review of the Literature on NAC Generation

Multistage stochastic programming problems that include endogenous uncertain parameters is a relatively new research area and has received far less attention in literature compared to its counterpart under exogenous uncertainties (Apap and Grossmann, 2017). To the best of our knowledge, Jonsbråten et al. (1998) were the first to introduce endogenous uncertainty by extending size selection problem, which was defined originally by Jorjani et al. (1999). Since then there has been an increasing number of studies applying endogenous uncertainty in gas-field development problem (Goel and Grossmann, 2004, 2006), in synthesis of process networks (Tarhan and Grossmann, 2008), in open pit mining scheduling problem by (Boland et al., 2008), in a pharmaceutical clinical trial planning problem (Colvin and Maravelias, 2008), in offshore oilfield development planning (Gupta and Grossmann, 2014), in vehicle routing problem (Hooshmand Khaligh and MirHassani, 2016), in various portfolio optimization of smart grid technologies (Giannelos et al., 2017), in artificial lift infrastructure planning problem (Zeng and Cremaschi, 2017), and in a new technology investment planning (NTIP) problem by (Christian and Cremaschi, 2018b).

However, application of MSSPs is limited due to curse of dimensionality. Rapid growth of the number of scenarios gives rise to a large number of NACs requiring a large number of additional indicator variables and logical constraints. Recent literature in the area focuses on the solution approaches to address the space and time complexity of MSSP models, especially due to the challenges associated with the large number of NACs.

Lagrangean relaxation approaches have received a significant attention. Here, NACs are removed or relaxed partially, which decomposes the MSSP into smaller problem of scenario groups. Gupta and Grossmann (2014) proposed a Lagrangean decomposition algorithm to obtain dual bounds for MSSP models under endogenous uncertainties. The authors decomposed scenarios into different scenario groups, only dualized the initial NACs linking two scenario groups, and kept initial and conditional NACs for the scenario pairs in the same scenario group. They applied this new Lagrangean decomposition algorithm to process-network-synthesis problems and reached 1% optimality gap within 30 iterations. Goel et al. (2006) presented a branch-and-bound algorithm that generates upper-bounds using the solution of the Lagrangean dual problem with relaxed NACs. The lower bounds were generated heuristically based on the solution of the Lagrangean dual problem. The results suggested that the branch-and-bound algorithm achieved significantly better solutions and tighter optimality gaps than the heuristic presented in Goel and Grossmann (2004). Tarhan et al. (2013) presented an approach that improves dual bounds for solving nonlinear MSSPs with decision-dependent uncertainty. The algorithm solves the Lagrangean relaxation of the dual problem where NACs have been removed to obtain the dual bound. The primal bound is generated by locating a feasible solution using a rolling-horizon approach.

Another approach utilizes scenario decomposition or feasible solution space decompositions. Goel and Grossmann (2004) considered an approximation approach, which searched a sub-space of the feasible region to find a “good” solution. They used this approach for solving a gas-field development problem, where the size and initial deliverability of reserves were endogenous uncertainties realized immediately after the drilling decisions were made. The sub-space is obtained by removing the scenario dependency of investment decisions for uncertain fields, yielding a more constrained version of the original problem. Then, this constrained version is solved to optimality. The results revealed that the solutions obtained to the constrained formulation were significant improvements compared to the deterministic solutions. Solak et al. (2010) developed a sample average approximation (SAA) based algorithm to generate candidate solutions for the R&D project portfolio optimization problem, where return levels for projects were endogenous uncertain parameters. Apap and Grossmann (2017) developed a sequential scenario decomposition (SSD) approach for MSSP models under both endogenous and exogenous uncertainties. The heuristic first solves endogenous sub-problems by selecting one scenario from each exogenous scenario group. Each sub-problem includes first-period NACs and endogenous NACs, but they do not include exogenous NACs. The approach solves a series of sub-problems with endogenous uncertainties to determine binary investment decisions, fixes the binary investment decisions to satisfy the first-time period endogenous and exogenous NACs, and solves the resulting model to obtain a feasible solution. Apap and Grossmann (2017) applied SSD to the process synthesis problem (Tarhan and Grossmann, 2008; Goel and Grossmann, 2006) and the offshore oilfield planning problem (Tarhan et al., 2009), and obtained high-quality feasible solutions.

In terms of branch and bound approaches, Tarhan and Grossmann (2008) proposed a solution strategy that used a duality based branch-and-bound algorithm for a MSSP model for solving process network synthesis problem, where the process yields are endogenous uncertain parameters with gradual realization as process operating decisions are made. The strategy was able to solve a 16 scenario example to within 3% optimality. Tarhan et al. (2009) explored an offshore oil planning problem with endogenous uncertainty in the initial maximum oil flowrate, and considered a single oil field with non-linear reservoir model and a gradual realization of uncertainty. Using a duality-based branch-and-bound, the authors located solutions that were up to 22% better than solutions obtained using an expected value approach. They, however, noted that the solution times for the MSSP were “rather long”. Christian and Cremaschi (2018a) developed a branch-and-bound algorithm combining a heuristic primal bounding approach, Knapsack problem based Decomposition Algorithm (KDA) (Christian and Cremaschi, 2015) and two dual-bounding approaches to solve large-scale MSSPs with endogenous uncertainty. The dual bounds were generated using a modified approach based on progressive hedging (PH), and by solving the individual scenario problems. The primal bounds were generated using the KDA. The branch-and-bound algorithm successfully addressed the space complexity of MSSP problems, but required considerable time to converge.

Colvin and Maravelias (2008) introduced a MSSP model for planning clinical trials in pharmaceutical R&D pipelines. The authors explored a rolling-horizon approximation approach, which yielded tight feasible solutions for large instances of the problem (Colvin and Maravelias, 2009), where the planning horizon was divided into a finite number of subsets. A relaxed MSSP is generated for the first subset by removing all inequality NACs for the stages beyond the first subset. The solution of the relaxed MSSP is implemented for the first subset, and related uncertainty is realized. Then, the process is repeated for each subset until the end of the planning horizon is reached. The authors were able to successfully solve cases with more than 1000 scenarios. Another approach, a branch-and-cut algorithm (Colvin and Maravelias, 2010a), initially adds a percentage of NACs and then iteratively adds constraints based on violations. The authors concluded that the algorithm reduced the number of NACs that should be included in the problem formulation significantly and that this method would be advantageous for any problem where the majority of constraints are NACs.

Vayanos et al. (2011) introduced a novel approach that uses decision rules for approximating MSSP problems under endogenous uncertainties. The approach employs piecewise constant and piecewise linear functions of uncertainties and can handle uncertain parameters with both discrete and continuous distributions. The authors applied the approach to solve oil infrastructure planning problem, and the approach yielded improved primal bounds with limited time budget.

In order to address the computational complexity caused by NACs and large size of scenarios, several authors considered approaches that remove the scenario and NAC structure to obtain good feasible solutions with limited computational effort. Mercier and Van Hentenryck (2008) proposed a multi-step anticipatory algorithm, which uses a sample average approximation approach to generate a Markov Decision Process (MDP). The MDP is then solved, and the greedy solution is returned. The algorithm was tested using 12 instances of the pharmaceutical R&D pipeline management problem. Solution qualities for the algorithm were 10% better than the dynamic programming equivalent for all instances. The authors concluded that the algorithm was, nevertheless, computationally expensive when applied to the pharmaceutical R&D pipeline management problems. Christian and Cremaschi (2015) proposed a knapsack-problem based heuristic (KDA) for pharmaceutical R&D pipeline management problems under endogenous uncertainties. The KDA generates and solves a series of knapsack problems based on realized outcomes and has been shown to obtain high-quality feasible solutions for problems with up to 1,048,576 scenarios. Zeng and Cremaschi (2017) developed a relaxed knapsack-problem based decomposition heuristic (RKDA) for the same problem to generate a valid dual bound. The relative gaps between KDA (primal bound) and RKDA (dual bound) were within 4.66% for problems with over 1 million scenarios. Christian and Cremaschi (2017) further extended the KDA with different sub-problem generation approaches: every time period generation scheme (ETP) and at each realization generation scheme (AER), and further tested the KDA for various instances of the problem. Zeng et al. (2018) developed a generalized knapsack-problem based decomposition algorithm (GKDA) to obtain feasible solutions for MSSPs under both endogenous and exogenous uncertainties. The GKDA is a scenario-free decomposition framework, which incorporates the initial and conditional NACs implicitly and addresses both space and time complexities of solving MSSP problems. The GDKA was applied to four planning problems and yielded feasible solutions with 0% to 20% relative gap for all cases.

In general, all approaches developed to date would realize significant computational benefits from reductions in numbers of NACs and scenarios as evidenced by heuristic approaches that employ scenario-free structures. A full set of NACs can be generated by relating each scenario to every other scenario. This produces a very large number of constraints which are often redundant. Several authors have presented properties of MSSPs, which reduce the number of redundant NACs. Applicability of a lot of these properties is limited to MSSPs where the scenario set is considered complete (i.e., given by the Cartesian product) (Goel and Grossmann, 2006).

Often situations arise where the set of scenarios considered in the MSSP is given by a sub-set of the full set of scenarios. For instance, every scenario in the full set generated by the Cartesian product may not be realizable. Alternatively, a sub-set of the full scenario set may be used in cases where solving the MSSP with the full scenario set is computationally prohibitive. Properties used to generate NACs that are applicable for the full set of scenarios do not necessarily apply to sub-sets of scenarios. Thus, recent advances have focused on the generation of NAC sets for sets of scenarios which are not generated by the Cartesian product. Hooshmand and MirHassani (2016) developed two approaches; the first one used a mixed integer linear programming formulation (MILP) to find the minimum cardinality NAC set, and the second one employed graph theory, which added and removed NACs until there were no violations of non-anticipativity. Apap and Grossmann (2017) extended the graph-theory-based approach to handle both endogenous and exogenous uncertainty. Boland et al. (2016) developed a series of proofs which define the conditions required to omit NACs from the problem formulation. In all these works, the authors used the concept of a differentiator set, i.e., a set of sets, which store information on the realization of uncertainty required for two scenarios to be distinguishable. The differentiator set is sufficient if the realization of the uncertain parameter is instantaneous, implying that the underlying value is revealed using a single indicator. If the realization of uncertainty is gradual there may be several indicators which result in partial revelation of uncertainty. The approaches developed to date cannot be applied to generate NACs of incomplete scenario sets for this type of uncertainty unless the underlying model of uncertain parameters is altered first.

2 Problem Formulations and Background

2.1 MSSPs with endogenous and exogenous uncertainty

The presentation of the general form of the MSSP is adapted from Apap and Grossmann (2017), Hooshmand and MirHassani (2016) and (Goel and Grossmann, 2006). Suppose that we have an MSSP where the planning horizon is defined as $\mathcal{T}=\{1,2,3,..,T\}$. We first define a set of sources of endogenous uncertainties, $i\in\mathbb{I}$, and a vector of uncertain parameters, $\theta_{i}$, associated with uncertainty source $i$. The realizable values of the uncertain parameter, $\theta_{i}$, are represented using a finite set $\Theta_{i}=\{\theta_{i}^{1},\theta_{i}^{2},\ldots,\theta_{i}^{\mathcal{R}_{i}}\}$. Further, we define $\mathbb{J}$ as the set of exogenous uncertain parameters. We denote by $\xi_{t}$ a vector of exogenous realizations defined for parameters $j\in\mathbb{J}$ and realized at time $t$. Each exogenous uncertain parameter is assumed to have a finite set of possible realizations, $\Xi_{jt}$.

The full scenario set is constructed using the Cartesian product, namely, the scenario set is defined as, $\mathbb{S}:=\prod_{i\in I}\Theta_{i}\times\prod_{(j,t)\in\mathbb{J}\times\mathcal{T}}\Xi_{jt}$. We concentrate on the cases where the actual set of scenarios is a subset of the full Cartesian product, and hence we also define the corresponding set $S\subset\mathbb{S}$.

Further, we are primarily interested in the case when uncertain parameters allow for gradual realization. Intuitively, it means that as time progresses, depending on the decisions made, the set of possible values for each uncertain parameter shrinks, until the true values are eventually realized. This then means that for each endogenous parameter $i$ there exists a subset of $2^{\Theta_{i}}$ which defines the permissible states of realized information about it. We will refer to transitions between these states as events. If all subsets of $2^{\Theta_{i}}$ are permissible, for each realization $\theta_{i}^{k}$ we can define an event $e_{i}^{k}$, indicating that $\theta_{i}^{k}$ has either been observed as the true realization of $\theta_{i}$ or eliminated. In general, events correspond to the transitions between information states and can be decoded in different ways, for example, see the manufacturing example described below, where events correspond to completion of processing stages. We can then define a set $\mathcal{C}$ as the set of permissible combinations of events, corresponding to possible trajectories of uncertainty realization. In the case of instantaneous uncertainty realization (or exogenous uncertainty) each parameter only has one event associated with it, indicating that with a single event the underlying value of the uncertain parameter is realized.

At every time period the decision maker has to select values for decision and the corresponding recourse variables. Observe that both should depend on scenario (i.e., the realization of uncertainty) and the set of realized events (i.e., the information available at the time of making the decision). We will then denote as $\boldsymbol{b}_{tsc}$ and $\boldsymbol{y}_{tsc}$ the decision and recourse variables at time $t\in\mathcal{T}$ under scenario $s\in S$ given realized events $c\in C$. Note that if uncertainties are realized instantaneously, there is no need to index on $c$. On the other hand, in the most general case of $C$ equal to the power-set over all realizations, the number of required variables can be prohibitively large.

NACs must be enforced in the optimization formulation of MSSP, in order to ensure that the decisions do not take into account information that is revealed in the future. We write the NACs in a general form as in (1).

[TABLE]

Disjunction (1a) provides the NACs, given that $Z_{t}^{rs}$ is defined as the indicator variable for whether scenarios $s$ and $r$ are distinguishable at time $t$, and is defined in equation (1b). Function $F$ is problem specific, and in (1b) $\boldsymbol{b}$ and $\boldsymbol{y}$ refer to the whole collection of decision and recourse variables.

Observe that NACs in (1) are written for all pairs of scenarios and all time periods. Depending on what the difference between the scenarios in a specific pair is, it may be possible to reduce the number of such constraints. Let us define $\tau(r,s)$ as the latest time period when scenarios $r,s$ are indistinguishable in exogenous uncertainty, i.e., $\tau(r,s)=\max⁡\{t\,|\,\xi_{t^{\prime}}^{s}=\xi_{t^{\prime}}^{r}\,\,\forall t^{\prime}\leq t\}$. Then if two scenarios differ in exogenous uncertainty only, then (1a) reduces to

[TABLE]

and there is no need for disjunction, or equivalently, $Z_{t}^{rs}=0$ for $t\leq\tau(r,s)$, and $Z_{t}^{rs}=1$ otherwise. If two scenarios differ in both exogenous and endogenous uncertainty, then (1a) can be written as

[TABLE]

Note that this construction can also be viewed from the point of view of Markov decision processes (MDP) with discrete state-space. At any point in time, the pair $r$ and $c$, i.e., scenario and realized event set, determines the information available to the decision maker on the realization of uncertainty, and hence corresponds to the states (i.e., contains the information necessary to make decisions). NACs along with indexing of the variables with $c$ and $r$ corresponds to the idea that an optimal solution to an MDP is a policy, i.e., a function prescribing actions to states. In the MDP framework NACs are usually not explicitly enforced, since the problem generally is not explicitly written in its entirety, and instead specialized algorithms are used.

In the next section we illustrate these constructions with a simple example, which we will then use to simplify the discussion of the proposed algorithm and as the basis for the numerical experiments.

2.2 An Illustrative Example

In order to illuminate the function of the endogenous NACs, specifically endogenous NACs with gradual realizations, here we present a simple demonstration, which we will refer to as the manufacturing example. We will then use this example throughout the paper to illustrate the developed algorithm. The example closely follows the description initially given in Colvin and Maravelias (2008).

Suppose that a company produces two different products, $p\in\{P1,P2\}$. Both products are required to complete three ordered processing stages using the same equipment. A product may develop a defect during any of the processing steps. If a product develops a defect, it does not complete the next step. Whether or not a product completes all processing stages without a defect is not known with certainty. We can represent this uncertainty associated with product $p$ using an uncertain parameter $\theta_{p}$. The possible outcomes for the uncertain parameter $\theta_{p}$ are $\Theta_{p}=\{\omega_{p}^{1},\omega_{p}^{2},\omega_{p}^{3},\omega_{p}^{4}\}$ where $\omega_{p}^{1},\omega_{p}^{2}$, and $\omega_{p}^{3}$ represent a defect developed on product $p$ in processing stages one, two, and three, respectively. The outcome $\omega_{p}^{4}$ represents the case where the product $p$ was successfully manufactured with no defects. In this example, the uncertainty is realized gradually, completing a processing stage will reveal if a defect is developed in that stage (a partial realization of uncertainty) but does not necessarily reveal whether or not the product will successfully complete all processing stages (complete realization of uncertainty). To have complete realization of uncertainty, all processing stages must be completed.

To illustrate how the non-anticipativity can be enforced for this example with the disjunctions (1a), we begin by defining a binary variable, $\chi_{p,sg,t}^{s}$, which takes a value of 1 indicating that processing stage $sg\in[SG]$ of product $p\in P$ is started at time period $t\in\mathcal{T}$ for scenario $s\in\mathbb{S}.$ Naturally it follows that we define a binary variable, $\zeta_{p,sg,t}^{s}$, which indicates whether the processing stage $sg$ has been completed for product $p$ at time $t$ in scenario $s$. For simplicity, in our example, we assume that each processing stage has a known duration for completion. We then can define intuitive constraint $\zeta_{p,sg,t}^{s}=\zeta_{p,sg,t-1}^{s}+\chi_{p,sg,t-\delta_{p,sg}}$. Depending on the objective of the problem, the model may contain additional constraints or variables. These variables and constraints do not impact the development of NACs, as such, we omit them from this discussion. Prior to completing any processing stages, all scenarios are indistinguishable. This implies that the values for the decision variable, $\chi_{p,sg,t}^{s}$, are identical until processing stages are completed.

For every scenario pair, $r,s\in\mathbb{S}=\times_{p\in P}\{\theta_{p}\}$, there exists one or more differentiating event. In the case of gradually realized uncertainty, the differentiating event is an event which causes a partial realization of an uncertain parameter. For example, consider two scenarios with the uncertain parameter realizations $(\omega_{P1}^{1},\omega_{P2}^{3})$ and $(\omega_{P1}^{2},\omega_{P2}^{4})$. The first scenario represents the case where product $P1$ develops a defect after processing stage one and product $P2$ develops a defect after processing stage three. In the second scenario, product $P1$ develops a defect after processing stage two and product $P2$ successfully completes all processing stages without developing a defect. Before the products complete any processing stages the scenarios are not distinguishable. Once product $P1$ completes the first processing stage, these two scenarios become distinguishable. If product $P1$ develops a defect during processing stage one, the underlying value for $\omega_{P1}$ is realized to be $\omega_{P1}^{1}$. If no defect is developed the value of $\omega_{P1}$ is realized not to be $\omega_{P1}^{1}$. To fully realize the value of $\omega_{P1}$, $P1$ must complete more processing stages.

In this example, the events that determine the gradual realization of uncertainty correspond to stage completion, hence we can define the set of events as $e^{k}_{i}$ for all $i,k$ as representing product $i$ completing stage $k$. Note that since the stages are completed in order the set of all permissible events can be defined as $\mathcal{C}=\cup_{p,q=0}^{3}\{e^{P1}_{1},\ldots,e^{P1}_{q},e^{P2}_{1},\ldots,e^{P2}_{p}\}$, where with slight abuse of notation we understand that $p=0$ (or $q=0$) corresponds to no events occurring for product 1 (or 2). Note also that since the products are independent, the variables $\chi$ and $\zeta$ defined above only need to be indexed on the events for the corresponding product. For any two scenarios $(r,s)$ and any set of events, we can establish whether the scenarios are distinguishable given the occurred events. We can then define $\Psi(r,s)$ representing the set of differentiating events that can occur causing scenarios $r$ and $s$ to become distinguishable. For example, if scenario $r$ corresponds to both products passing all stages; and scenario $s$ represents $P1$ failing at stage 2, and $P2$ failing at stage $3$, then the set of differentiating events consists of P1 failing at stage 2, and P2 failing at stage 3, i.e., we can write $\Psi(r,s)=\{e^{2}_{P1},e^{3}_{P2}\}$. We can then use $\Psi$ to define $Z_{rs}^{t}$ and the necessary disjunction based on variables $\zeta$, as in (1). The disjunction can be converted to a linear constraint using standard integer programming techniques.

3 The Proposed Approach for Generation of NACs for Scenario Subsets

3.1 Graph-based model for NACs

Generation of NACs relies on knowledge of the scenario structure and how the groupings of indistinguishable scenarios change with the realization of uncertainty. The algorithm proposed in this section will focus on minimizing the number of pairs of scenarios on which the NACs are enforced. Similarly, when discussing the reduction of the number of NACs we will compare the results in terms of the number of pairs of scenarios. Note that this is a standard approach in the literature, (see, for example Hooshmand and MirHassani (2016)).

In order to construct the algorithm, we will consider a graph $\mathcal{G}=(\mathcal{S},E)$ with vertices representing scenarios $s\in\mathcal{S}\subset\mathbb{S}$ and the edges representing NACs which relate the pairs of scenarios. Hence, a complete graph would correspond to the basic formulation, where NACs are enforced between any two scenarios. A graph with just some arcs, corresponds to a formulation with only some NACs enforced. The graph can be naturally visualized by placing the vertices (scenarios) in a $|\mathbb{I}|+|\mathbb{J}|$ dimensional grid, according to the corresponding realizations of the uncertain parameters. Before discussing the algorithm, we will establish some properties of the graph and how the grouping of vertices on the graph changes with the changes in the knowledge of the underlying values of the uncertain parameters in the system.

Consider the set of permissible event subsets $C$, defined earlier. Any event set $c\in C$ defines a partition of the scenario set $S$,

[TABLE]

where two scenarios $r,s$ are in the same set $S_{j}^{c}$ if the events in $c$ do not differentiate between $r$ and $s$. In the remainder of this paper we will refer to such subsets as cuts. Figure 1 illustrates this construction for the manufacturing example discussed earlier. Specifically, Figure 1a shows the full scenario set arranged in a grid according to realizations of uncertainty, while Figures 1b and 1c give the partitionings for two event sets.

Next, we will describe the relationship between a set of constraints sufficient to enforce nonanticipativity and connectivity of subsets of nodes. At the initial time period no uncertainty has been realized. This implies that all scenarios are identical, thus there is one group of scenarios (the whole set). Throughout the planning horizon uncertainty is realized. As uncertainty is realized, the scenario set is divided into subsets based on the realizations. Figure 2 shows the progression of uncertainty realization in the manufacturing example. Throughout the paper we will refer to the cardinality of event set $c$ as the degree of uncertainty realization, since $|c|$ represents the number of differentiating events which have occurred in order to reach the state of uncertainty realization.

Lemma 1**.**

If $c\subset c^{\prime}$, then partition $\Pi(c^{\prime})$ is a refinement of $\Pi(c)$, that is for every $j^{\prime}$ there exists $j$ such that $S_{j^{\prime}}^{c^{\prime}}\subset S_{j}^{c}$.

Proof 1**.**

*The claim is a direct consequence of construction of partitioning. Indeed, $c\subset c^{\prime}$ implies that the state of the system corresponding to events in $c^{\prime}$ can be achieved as an evolution of the system from a state corresponding to events in $c$. Consequently, scenarios that are indistinguishable in this future state had to be indistinguishable in the past as well. *

If two scenarios $r,s\in S^{c}_{j}$ for some cut $c\in\mathcal{C}$, then they are indistinguishable for some permissible state of the system, and hence, non-anticipativity should be enforced for these two scenarios. On the other hand, a direct non-anticipativity constraint for $(r,s)$ can be eliminated if another path $\mathcal{P}$ exists within nodes in $S_{j}^{c}$ connecting $r$ and $s$ due to transitivity of the constraints. Further, since while uncertainty is gradually realized new scenario sets are constructed as subsets of existing sets (Lemma 1), the sets formed by the differentiating event must exist as connected sets prior to the occurrence of the differentiating event. This observation leads to the following characterization of relationship between NACs and the corresponding graph.

Lemma 2**.**

A set of constraints is sufficient to enforce nonanticipativity, if and only if for all cuts $c\in\mathcal{C}$ all of the corresponding scenario subsets $S^{c}_{j}$ are connected in the corresponding induced subgraph.

Proof 2**.**

Follows by construction, from the discussion above.

Next, we will demonstrate how using the knowledge of the structure of the subsets of vertices formed by each event set and the relationship between event sets, it is possible to develop an algorithm which generates a set of edges for the graph G such that it satisfies Lemma 2 with the minimal number of corresponding pairs of scenarios.

Definition 1**.**

An edge $(r,s)\in E$ will be called necessary for a set of nodes $S^{c}_{j}\in\Pi(c)$ if $r,s\in S^{c}_{j}$, and $r,s\notin S^{c^{\prime}}_{i}$ for all $S^{c^{\prime}}_{i}\in\Pi(c^{\prime})$ and $c^{\prime}\in C$ such that $|c|<|c^{\prime}|$.

By this definition, an edge $(r,s)$ is necessary for $S^{c}_{j}\in\Pi(c)$ if non-anticipativity for scenarios $r,s$ is required (as $r,s\in S^{c}_{j}$), yet it cannot be established based on edges that are within subsets in $\Pi(c^{\prime})$ for $|c^{\prime}|>|c|$, i.e., when more information on uncertainty realization is available. In other words, unless non-anticipativity for $r,s$ is enforced with a path within $S_{j}^{c}$, then it will not follow from any other set of NACs.

This then leads to an intuitive idea for generating a set of NACs. The proposed algorithm is presented in Algorithm 1. It is constructed as a greedy procedure, which begins by considering the fully realized uncertainty set ($k=\max\{|c|,c\in\mathcal{C}\}$) and then generates subsets of necessary edges separately for each possible event set $c\in C$ such that $|c|=k$. Any edge deemed necessary is stored, then the value of $k$ is decreased. The process is repeated until $k=0$. The resulting NAC set is represented by the edges deemed necessary in all iterations.

An example progression of the algorithm is presented in Section 3.2.

Theorem 3**.**

Algorithm 1 terminates with a subset of edges which corresponds to a set of NACs enforced on the minimum number of pairs of scenarios. The algorithm can be implemented in a way that guarantees that the running time grows as $O(|\mathcal{C}||S|^{3})$ as $|\mathcal{C}|,|S|\rightarrow+\infty$, where $S$ is the number of scenarios and $\mathcal{C}$ is the the set of permissible event sets.

Proof 3**.**

*Correctness and optimality of the resulting NAC-set follow by construction and Lemmas 1 and 2. Indeed, the algorithm explicitly enforces connectivity on subsets $S_{j}^{c}$ for all $c$ and $j$. Further, all arcs added at step 6 of the algorithm are necessary in the sense of definition 1 by construction. Hence it is impossible to remove any of the arcs without violating connectivity. Indeed, suppose we remove an arc connecting scenarios $r$ and $s$. In this case, the NACs for these two scenarios are not enforced since the arc was added as a part of the spanning tree, and these NACs could not have been enforced with any of the previously added arcs by construction.

Finally, while each component $S_{j}^{c}$ potentially can be connected in multiple ways, this choice has no effect on the number of arcs needed to be included at future iterations. This claim follows from Lemma 1 and the way each subset $S_{j}^{c}$ is constructed by combining a number of smaller subsets $S_{i}^{c_{i}}$ with $|c_{i}|>|c|$, all of which are already connected since they have been considered at a previous iteration.*

In order to establish running time estimate, observe that the algorithm iterates over all $c\in C$ and then over all $S_{j}^{c}\in\Pi(c)$. For each $S_{j}^{c}$ it then generates the set of scenarios in $S_{j}^{c}$ and introduces arcs to enforce connectivity. The number of subsets $S_{j}^{c}$ for each $c$ is bounded by $|S|$. In a simplest implementation connectivity of a graph can be enforced with a greedy algorithm constructing a spanning tree in $O(n^{2})$, where $n$ is the number of nodes. In our case, the number of nodes in each subset $S_{j}^{c}$ is bounded by $|S|$. Further, generation of each subset $S_{j}^{c}$ can be organized in $O(|S|)$ by enumeration. This then results in overall $O(|C||S|^{3})$ running time.

Remark 1**.**

Each component $S_{j}^{c}\in\Pi(c)$ for $|c|=k$ is treated “in parallel”, i.e., spanning trees to enforce connectivity are constructed without regard to edges that are added at the same iteration. This can lead to creation of cycles in $\mathcal{G}$. When connectivity is enforced on any subset, each connected component is treated as a single node.

The running time is proportional to $O(|\mathcal{C}|)$, which, in the worst case, is a powerset of the set of uncertainty realization events, i.e., it is exponential in the number of events. On the other hand, if the number of uncertain parameters and events is fixed, and $S\ll\mathbb{S}$, which is often the case, then the running time is polynomial in the number of scenarios.

Remark 2**.**

Note that the description of the algorithm is provided in terms of event sets and does not emphasize whether underlying uncertainty is endogenous or exogenous. Indeed, while in the case of exogenous uncertainty, one can trivially avoid having to model gradual realizations, its can still be processed by SNAC in the exact same way. For an exogenous uncertain parameter (eg., demand in the manufacturing example) the corresponding event set $\mathcal{E}_{\sigma}$ consists of a single member, yet the algorithm can still be allied in the same manner as for purely endogenous cases. As an illustrative example, consider a case similar to the manufacturing example with single product, and instead of the uncertainty due to processing stages of product two, we consider uncertain demand. Suppose that demand parameter can take any of the values $\{D1,D2,D3,D4\}$, which are gradually realized in two phases: at time $t=1$ the demand is revealed to be either $\{D1,D2\}$ or $\{D3,D4\}$ and then at time $t=2$ the actual value is observed. Then, following the same construction as in the example above, we can denote events as $\{t=1,t=2,e_{p1}^{1},e_{p1}^{2},e_{p1}^{3}\}$ and then construct set $\mathcal{C}$ and the diagram analogous to Figure 2 as given on Figure 3. Consequently, the algorithm can be applied in the same way as for the purely endogenous case.

3.2 Illustration of the Algorithm Progression for a Simple Manufacturing Example

The manufacturing problem considered is described in detail in Section 2.2. The cardinality of the full scenario set for the manufacturing problem with two products and three processing stages is 16 (calculated as $4^{2}=16$). For this example, we sample 6 random scenarios, $S=\{(\omega_{P1}^{1},\omega_{P2}^{1}),(\omega_{P1}^{4},\omega_{P2}^{3}),(\omega_{P1}^{2},\omega_{P2}^{1}),(\omega_{P1}^{3},\omega_{P2}^{2}),(\omega_{P1}^{4},\omega_{P2}^{1}),(\omega_{P1}^{3},\omega_{P2}^{3})\}$. Figure 4 gives the initial graph corresponding to these scenarios arranged in a grid.

Table 1 shows the event set $C$ (organized by $|c|$, referred to as order).

The algorithm starts by initializing $k$ to six. The event set $c$ with $|c|=6$ produces subsets where each scenario forms its own subset. As a result, no edges are deemed necessary. Reducing $k$ by one and generating new subsets yields the groups given in Figure 5. The corresponding event set with five events is shown below each set of subsets. In Figure 5(A), two scenarios fall into the same subset and the corresponding necessary edges are represented with dotted lines. Figures 6 and 7 illustrate the remaining progression of the algorithm. It generated NACs for a total of five pairs to connect six scenarios. Without reductions the number of pairs of scenarios with enforced NACs would have been 15 calculated by $0.5[|\mathcal{S}|(|\mathcal{S}|-1)]$.

4 Computational Studies

In this section we present results of a computational study. The aim of the study was to evaluate whether the proposed algorithm is able to significantly reduce the number of scenario pairs requiring enforcement of NACs and further, whether it then translates into savings in terms of MSSP solution time. All computational studies were implemented in Python Version 3.5.0 and performed on a 64-bit machine running Windows 7 with Xeon CPU E3-1241 and 32 GB RAM. The resulting optimization problems were solved in IBM ILOG CPLEX.

MSSP instances solved were generated in the same way as the models presented in Colvin and Maravelias (2008). The problems follow the same logic as the manufacturing example we used above with more products, stages and scenarios. For completeness, we present the problem formulation and parameter values in A and B respectively. For a more detailed description of the cases, see Colvin and Maravelias (2008). Note that here we use the same terminology as in the manufacturing example to describe the observed results. We vary the number of products, the number of possible outcomes (processing stages) and the number of scenarios sampled from the full Cartesian product. In each case, random instances are generated 30 times. In each instance generated, we solved both the original MSSPs without identifying the minimum scenario pairs (referred to as “Full”) and the MSSPs constructed using the proposed algorithm (referred to as “SNAC”). Tables 2 and 3 report the information on the instance sizes and running time respectively. It should be noted that the solutions of both models yielded the same objective function value for all instances and scenario subsets.

First, consider Table 2, which reports the number of NACs and the corresponding number of scenario pairs due to SNAC (average and maximum over 30 random instances) and compares it with the full problem. We can observe that both values quickly grow as problem size increases, yet in all cases the algorithm significantly reduces the number of NACs and scenario pairs. This result is not suprizing, since it is clear that full formulation contains a significant number of redundant constraints. Next we investigate whether reduction in NACs leads to a reduction in solution time.

It should be noted that the SNAC only utilizes the uncertain parameter and outcome information without considering the planning horizon or potential completion times. The number of NACs for certain problems including these case studies can be further reduced if the task completion times and planning horizon is considered as shown in the NAC reduction properties by Colvin and Maravelias (2010b).

Table 3 reports running time. Specifically, we report three separate times: SNAC algorithm time, model generation time (i.e., the time required to generate NACs) and CPLEX time (time used by the solver once the model is loaded). Here we can observe a trade-off, whereby on one hand, in all cases the models created by SNAC take significantly less time to both generate and to obtain a solution; yet on the other hand, SNAC itself may require a significant time to complete. For example, in case of 5 products, 4 outcomes and all 1024 scenarios, the algorithm used on average 5140.6 seconds to run, which then allowed to reduce model generation time from 3454.9 to 229.4, and CPLEX time from 4331.7 to only 43.5. This can be very beneficial for two reasons. First, the total time with SNAC ($5140.6+229.4+43.5=5413.5$) is 30% less than the total time for the full problem ($3454.9+4331.7=7786.6$). Secondly, note that SNAC only uses scenario structure and is independent from the rest of model, so if a decision maker is interested in solving a collection of problems, which differ in parameter values (e.g., costs, revenues, demand, etc), they only need to run SNAC once, and all subsequent problems can be solved with the same optimal NAC structure. In the case of full problem formulation, each new problem would need to be solved from from scratch.

To sum up the case study, we conclude that our proposed algorithm can significantly reduce teh number of NACs that need to be enforced in an MSSP with gradual realization of uncertainty. Further, this reduction can directly result in computational time savings. These savings are especially significant if one needs to solve a collection of problems with the same scenario structure.

5 Conclusions

The premise of this paper was to develop an algorithm which would generate the minimum cardinality non-anticipativity constraint set for multistage stochastic programs where the uncertainty is endogenous and the realization of the endogenous parameter is gradual. The algorithm utilized knowledge of the subsets of scenarios formed by the realization of uncertainty in order to construct necessary and minimal sets of non-anticipativity constraints. A running time analysis showed that algorithm scaled as $O(S^{3})$ as $S\rightarrow+\infty$ as long as the number of uncertain parameters is kept constant. Illustrative computational experiments confirm the scalability of the algorithm with respect to the number of scenarios and demonstrate significant reduction in the number of NACs generated and corresponding improvement in MSSP solution time.

Acknowledgments.

We acknowledge support by the US National Science Foundation through the Auburn University Integrative Graduate Research and Education Traineeship (IGERT) program (Award # 1069004).

Appendix A MSSP Formulation used in the study, following Colvin and Maravelias (2008)

[TABLE]

Appendix B Case details

The parameter values used in the case studies are reported in Tables 4, 5, 6, 7, 8.