An approach for analyzing and managing flexibility in engineering systems design based on decision rules and multistage stochastic programming

ABSTRACT This article introduces an approach to assess the value and manage flexibility in engineering systems design based on decision rules and stochastic programming. The approach differs from standard Real Options Analysis (ROA) that relies on dynamic programming in that it parameterizes the decision variables used to design and manage the flexible system in operations. Decision rules are based on heuristic-triggering mechanisms that are used by Decision Makers (DMs) to determine when it is appropriate to exercise the flexibility. They can be treated similarly as, and combined with, physical design variables, and optimal values can be determined using multistage stochastic programming techniques. The proposed approach is applied as demonstration to the analysis of a flexible hybrid waste-to-energy system with two independent flexibility strategies under two independent uncertainty drivers in an urban environment subject to growing waste generation. Results show that the proposed approach recognizes the value of flexibility to a similar extent as the standard ROA. The form of the solution provides intuitive guidelines to DMs for exercising the flexibility in operations. The demonstration shows that the method is suitable to analyze complex systems and problems when multiple uncertainty sources and different flexibility strategies are considered simultaneously.


Introduction
The design of complex engineering systems-such as real estate development projects, airports, bridges, power plants, telecommunication, and Waste-To-Energy (WTE) systems-is challenging, as they typically have a long life span and operate in uncertain environments. This uncertainty inevitably affects their expected life cycle performance. Flexibility in engineering design is an emerging paradigm in engineering design research that is inspired by the economic theory of real options, which is an analytical approach for quantifying the value of flexibility in irreversible investment projects (Dixit and Pindyck, 1994;Trigeorgis, 1996). It aims to provide "the right, but not the obligation" to change a system easily in the face of uncertainty. Flexibility aims to improve the expected life cycle performance by reducing exposure to downside risks (i.e., acting like an insurance policy) and enabling the system to capture upside opportunities (i.e., like a call option on a stock). It is an important property for complex systems and has been shown in many industry case studies to significantly improve expected life cycle performance as compared with standard design and project evaluation approaches (Chen and Yuan, 1999;Nilchiani andHastings 2007, Suh et al., 2007;Alp andTan 2008, Ross et al., 2008;de Neufville and Scholtes, 2011;Jiao, 2012;Koh et al., 2013;Luo, 2015;Subrahmanian et al., 2015).
The motivation for this work is twofold. First, there is a need to provide intuitive guidance to Decision Makers (DMs) on when and how to exercise the flexibility in operations.
CONTACT Michel-Alexandre Cardin macardin@nus.edu.sg This is important because a well-designed system can display a better life cycle performance only if it is implemented following an optimal or best design strategy, in a stochastic sense. As uncertainty evolves over time, it may be challenging to make decisions on when to exercise flexibility so as to achieve better performance, especially when many stakeholders and thousands of design variables are involved. Second, it is difficult to access the value of flexibility in large-scale engineering system projects. These systems typically face multiple uncertainty sources (e.g., market prices and demand, amount of natural resources, regulations, technology, etc.). In addition, several flexibility strategies may be implemented simultaneously to extract additional value from uncertainty. How to analyze and manage different flexibilities concurrently in the face of multiple sources of uncertainty can be a challenging task. Recent efforts have focused on developing practical approaches for assessing the value of flexibility in an engineering setting, to help designers and DMs implement it and manage it in operations. The work on real options and flexibility in engineering design aims to do just that, by relaxing and/or modifying some of the underlying assumptions behind standard Real Options Analysis (ROA) to better suit the needs of industrial and systems engineers. The emphasis is on enabling the ranking and evaluation of different design alternatives (e.g., rigid versus flexible design), based on their expected life cycle performance. In this context, this article explores the question of how one can assess the value of flexibility using an approach different from standard ROA, with a form of solution that is more intuitive to use in operations, while recognizing the inherent value stemming from flexibility. This article proposes an approach based on decision rules and multistage stochastic programming to assess the value of flexibility and determine the best design and exercise strategies in operations. A decision rule is a heuristic or triggering mechanism used by DMs to determine when it is appropriate to exercise flexibility. This approach contrasts with and compares to standard ROA by emulating directly the decision-making process and parameterizing the characteristics as well as the physical design variables, so they can be analyzed thoroughly using optimization techniques. A multistage stochastic programming model is used to find the best decision rules and design variables to implement and manage the flexible system in operations. Such an approach is in line with Simon's view (Simon, 1972) that human decision agents may tend to rely on heuristic rules to achieve a stated satisfactory level of performance when operating under complex and uncertain environments. Under the proposed framework, a typical solution consists of two parts. One part contains the design variables that describe the initial physical state of the system. The other part contains the decision rules that guide decision-making dynamically based on available information at the time of a decision, and as uncertainty unfolds. A hybrid WTE system is used as an example application to compare the performance of the decision rule approach with a standard ROA that relies on multinomial lattice analysis and augmented using Approximate Dynamic Programming (ADP) to address the curse of dimensionality. The complex system involves multiple uncertainty sources, flexibility strategies, and types of decision rules. The approach aims to recognize similar value stemming from flexibility as compared to standard ROA. The form of the solution aims to provide intuitive guidelines on how and when to exercise the flexibility in operations. A special case of the system involving only one uncertainty source and one flexibility strategy is also analyzed to demonstrate the use of the solution in comparison with standard ROA.
The remainder of this article is organized as follows. Section 2 introduces the relevant background to this work. Section 3 introduces the model formulation and solution procedure for the proposed decision-rule-based approach. In Section 4, the standard ROA and the proposed decision rule approaches are used and compared in the analysis of a hybrid WTE system with multiple uncertainty sources and flexibility strategies. In Section 5, a discussion comparing the proposed approach with standard ROA is provided, in addition to a discussion on the limitations and opportunities for future work.

ROA
ROA is a systematic approach that relies on financial options theory to assess the value of flexibility in irreversible investments in real physical assets operating under uncertainty (Mun, 2006). Classic ROA approaches are based on the Black-Sholes options-pricing model and emerge from solving a system of stochastic partial differential equations and using dynamic programming (Black and Scholes, 1973;Dixit and Pindyck, 1994). A simplification of the financial options model was introduced by Cox et al. (1979), who proposed a binomial lattice as a discretetime approximation of the Black-Scholes options pricing formula. The binomial lattice approach has been widely used in capital budgeting and strategic decision-making, as it explicitly accounts for the value of flexibility (Trigeorgis, 1996). A considerable amount of literature exists that is concerned with valuing a capital asset and determining the appropriate time to exercise the option (McDonald and Siegel, 1985;Nembhard et al., 2005).
The ROA literature provides some guidance on how to value flexibility and when to exercise it in operations. Determining the optimal time to exercise, however, may be challenging in practice. For instance, in binomial lattice analysis, one needs to determine the current stage and state in the lattice by fitting past historical data. From there, future evolution of the uncertainty driver is projected, a backward induction process must be applied up to the decision point based on a pre-defined recursive formula, and DMs must choose an exercise policy based on the highest expected reward. On the other hand, wide adoption of ROA techniques has been limited in practice, partly because it relies on advanced mathematical concepts that may not be intuitive to DMs (Engel and Browning, 2008). Another challenge is that standard ROA relies on assumptions that apply well to finance, but not necessarily to an engineering setting. For example, the path independence assumption, crucial to path recombination and computational savings, may not hold (Wang and de Neufville, 2005). Due to path dependencies that are inherent in complex engineering systems, system value after an up-down movement may not be the same as after a down-up movement, because the decision sequences and resulting physical artifacts may differ. The analysis of a flexible system considering two or more independent, but intertwined, uncertainty drivers and flexibility strategies can be a challenging task. The development of a multinomial lattice is required if two or more uncertainty drivers are to be analyzed (Copeland and Antikarov, 2003). The curse of dimensionality that arises in dynamic programming analysis may make such approaches difficult to use when more complex systems and problems are considered (Nembhard et al., 2005).

Flexibility in engineering systems design
Analytical tools better suited for an engineering context have been introduced in the literature, such as decision analysis, and Monte Carlo simulations. These techniques can be separated into three categories. The first one comprises analytical tools based on backward induction and dynamic programming, including binomial lattice and decision analysis, see, for example, Min and Wang (2000), Brandão et al. (2005), and Buurman et al. (2009). The second category relies on mathematical optimization techniques, see, for example, Zhao and Tseng (2003) and Jiao (2012). The third category includes simulation models, that explicitly model the stochastic scenarios and introduce heuristic-based decisions enabled by the flexibility; see, for example, Cheah andLiu (2006), de Neufville et al. (2006), and Zhang and Babovic (2011). However, all of the above techniques have limitations for practical application. All of the decisions on flexibility exercise are determined as outcomes of the model; there are no generic decision rules to guide managers as the uncertainty is resolved throughout the project life cycle.
Exercising the real option is typically based on a process inspired from Bellman's recursive formula, which may not fully capture the actual decision-making process in reality (Bellman, 1952). Although the simulation approach explicitly employs decision rules, there is currently no systematic approach to determine the most value-enhancing decision rule parameters.

Decision rules
A decision rule, or implementable policy, is a function that maps the observations of uncertainty data to decisions (Shapiro et al., 2009). The concept of decision rules has been used in various multistage stochastic programming models (Kuhn et al., 2008;Georghiou et al., 2015). Garstka and Wets (1974) surveyed the different types of decision rules that are used in stochastic programming. Four classes of decision rules were identified and studied via case examples: zero-order, linear, safety-first, and conditional-go decision rules-see Section 3.2 for further details. The authors pointed out that there are difficulties in restricting the class of acceptable rules to those specific functional forms to derive optimal rules, unless some highly structured problems are being considered. Nevertheless, any approximation of the optimal rule by a certain class of decision rules would yield quantitative bounds that have significant theoretical as well as practical implications.
Even for some intractable problems, a good approximation to an optimal solution can be obtained by searching for specific classes of decision rules. This is a utilitarian view in the paradigm investigating how to enable flexibility in design and management of engineering systems. Due to the complexity of the systems and uncertain environment, the design problems are usually of large size and thus cannot be easily solved. By formulating the problem using decision rules, it becomes possible to obtain a good approximation to the optimal solution typically found using standard ROA techniques. Also, generic decision rules have an advantage in that they are practical and intuitive for DMs to use in operations. They can be modeled explicitly to emulate the decision-making process at a firm or organization.
As such, this article recognizes the importance of decision rules as a way to assess the value of flexibility in design and management of complex systems. A systematic approach is proposed to find stochastically optimal decision rules that guide dynamic decision-making based on available information, as uncertainties are resolved.

Multistage stochastic program with decision rules
In this section, a multistage stochastic programming model is proposed to explicitly consider decision rules in the analysis of flexibility in the design of complex engineering systems.

Model formulation
Let ξ = (ξ 1 , . . . , ξ T ) be a scenario of uncertainty, where ξ t is the uncertainty to be observed in period t. Note that ξ t is a vector and thus that it can represent multiple sources of uncertainty. Denote as the set of all possible uncertainty scenarios. For simplicity, it is assumed that is a finite set {ξ 1 , . . . , ξ K } with corresponding probabilities p k ≥ 0, K k=1 p k = 1. There is no operation at the beginning of planning horizon, time 0, which represents the present time. The DM must make a decision x t at the beginning of period t from a set of feasible decisions X t when t ≥ 1. Note that x t is a vector, which means that it can simultaneously consider several flexibility strategies. Let X ⊆ X 1 × · · · × X T denote the set of all feasible decision sequences x, where x = (x 1 , . . . , x T ). A decision rule, or implementable policy, δ, is a function that maps each scenario of uncertainty ξ in to a sequence of decisions x in X (i.e., δ : → X ). Let D denote the set of all mappings from to X. The form of δ varies in different problems, and a vector of parameters θ is used to characterize it. Thus, the decision rule can be represented as δ θ .
In a multistage setting, the uncertain data ξ 1 , . . . , ξ T are revealed gradually over time. Let notation ξ [t] := (ξ 1 , . . . , ξ t−1 ) denote the history of the uncertainty realization up to the beginning of period t. Nonanticipativity in the DM's actions requires that the choice of decision x t at the beginning of period t depends only on the information on uncertain data ξ [t] available up to the beginning of period t. In other words, if the revealed uncertainty data up to the beginning of period t in two scenarios are the same, then the decision made for period t should be exactly the same, no matter how different they will be in the following periods. Exploiting this requirement, δ θ (ξ) can be repre- The DM's goal is to select a feasible decision rule to maximize the total expected profit (or reward function). The profits are determined by a sequence of profit func- where the profit in period t depends on the decision x t and the revealed uncertain data ξ t in period t. Let [t] , ξ t denote the total profit, where π is the discount rate. The problem of choosing an optimal decision rule is then where is a subset of D, the variables with superscript k correspond to the variables in scenario k. If this model can be solved to an optimum value, then an optimal decision rule δ * can be obtained, under which the DM will do best by choosing δ θ * (ξ) to attain the maximum expected total profit as uncertainty ξ unfolds. Even with the finiteness of , however, the problem may be very large. In the case of = D, there is typically no analytical form for the optimal policy δ * . In the design of flexible engineering systems, one is usually interested in certain types of flexibilities embedded in the system, which will lead to some particular patterns of decision rules. Therefore, to obtain the closed-form solution to the decision rule problem, can be restricted to a class of mappings, which is a subset of D. By restricting to be a proper subset of D, it is possible to solve the problem, which is generally computationally intractable.
The choice of the form of the decision rule is based on the characteristics of the system design problem. The decision rule approach is flexible enough to model engineering design problems, considering two important features. First, a decision rule approach can readily handle multiple uncertainty sources. This is because of the structure of the model where ξ t represents a vector of uncertainties. Thus, increasing complexity from a single uncertainty source to multiple uncertainty sources will not significantly increase the complexity of the analysis. Second, it is straightforward to consider multiple flexibility strategies simultaneously. This is because x t is a vector of decisions, with each element representing a decision regarding a given flexibility strategy. This formulation also lends itself naturally to solution approaches such as Lagrangean decomposition and can efficiently handle a considerable number of scenarios.

Types of decision rules
Three generic types of decision rules that can be readily applied to the management of flexibility in engineering systems are considered, as discussed below. Decision rules specific to a given system may be generated using systematic concept generation processes and typically rely on the designer's creativity and expertise with the system. Decision rules can be created by considering generic real option strategies (e.g., capacity expansion, abandonment, switching, phasing capacity deployment, deferring investment, etc.) combined with managerial and/or mechanical design enablers (Cardin, 2014). Thus, many decision rules exist and can be explored, depending on the system of interest, the uncertainty drivers, and system mission or purpose.

... Conditional-go decision rule
With conditional-go decision rules, the decision at each time stage is made based on an estimate of future conditions and past information. The rule is usually expressed as "if the uncertainty realizations in the past satisfy a certain set of criteria, then exercise the flexibility, else maintain the status quo. " There are various ways to formulate the conditional-go decision rules. Binary variables are usually introduced to facilitate this process. A common formulation in the proposed stochastic programming setting is where e t is a vector of binary variables, each element of e t representing the decision to exercise or not a flexibility: if it is equal to one, then the flexibility is exercised; if equal to zero, it is not exercised. 1(·) is an indicator function. g(ξ [t] ) is a function of ξ [t] and can be an estimate of future conditions based on ξ [t] . The φ t are the criteria to be satisfied in order to exercise the flexibility at the beginning of period t.
Due to the use of binary variables, a problem relying on conditional-go decision rules can often be formulated as a Mixed-Integer Linear Program (MILP). The class of conditional-go decision rules is a compromise between precision and computational tractability. Although it cannot necessarily guarantee global optimality, it is often helpful in finding good enough, value-enhancing practical policies. The form of the solution is useful in the management of flexibility in practice, since it emulates a system operator's decision-making processand significantly simplifies the procedure.
To make a connection with standard ROA, the formulation of Bellman's equation in the binomial lattice can be thought of conceptually as a specific implementation of the conditionalgo decision rule. The decision rule can be described as "if the expected value of the next stage is higher when exercising the option, then exercise it, if not, do nothing. " This decision rule is applied recursively in the backward induction process, starting from the last until the initial stage, or the time of exercise. This decision rule is implicitly enforced in the structure of the binomial lattice. There is less freedom to implement different types of decision rules and also decision rules that directly emulate the system operator's decision-making process, such as those considered as examples in this study.

... Linear decision rule
Linear decision rules are a class of decision rules δ t (·) that are linearly dependent on the observed uncertainty data ξ [t] . They can be expressed as δ t (ξ [t] ) = L t ξ [t] , for some matrix L t with proper dimensions. Substituting these linear decision rules into problems (1)-(4) yields the following approximate problem: Using linear decision rules, the size of the original stochastic programming problem grows only moderately with the number of time stages. Therefore, the original problem can be converted into a finite linear program that is amenable to a numerical solution. The form of the solution is also easy to use in operations, since it is a linear combination of the observed uncertainty realizations.

... Constant decision rule
If the decision to be made at time stage t, δ t (·), is independent of the observations of the uncertainty data up to stage t, ξ [t] , then the decision rule belongs to a class called constant decision rules. This is a special case of a multistage problem: the problem of determining a decision rule function degenerates into determining a single sequence of decisions to be made at each stage. Thus, it actually enables the DM to reduce a multistage problem into a two-stage problem. Only one linear program needs to be solved to determine the options to be chosen at each stage. In this way, the computational cost can be reduced dramatically: The shortcoming of this class of decision rules is that it does not follow a fully rational process, as it forces the DM to make decisions for all time stages before realizing any of the uncertainties. In some cases, however, it can provide an insurance policy against the worst possible outcomes.

Solution procedure
An important challenge in multistage stochastic programming is the large dimension that makes most problems NP-hard to solve. Even if the problem can be simplified into a deterministic structure by assuming the finiteness of , the computational effort still increases significantly with the number of scenarios. In order to obtain a good approximation for the uncertainty, however, a large number of scenarios are usually needed. In the cases where integer decision variables are involved, the problem quickly becomes computationally intractable. To overcome this difficulty, decomposition is an effective approach.
It can be observed that the structure of problems (1)-(4) is convenient to be split into K small sub-problems, one for each scenario. The optimal value of problems (1)-(4) is then equal to the weighted sum with weights p k of the optimal value of the sub-problems. It is the nonanticipativity requirement, however, that links the decision sequences associated with different scenarios. The nonanticipativity constraints (3) is crucial to simplify the exercise of flexibility in operations. It can be expressed by ensuring the equality of the parameters θ in each scenario; thus.
where θ k is a replication of the decision rule parameters in scenario ξ k . Another way to write the nonanticipativity is to require that whereθ = K k=1 p k θ k is the average of θ k over its K replications. The size of the original problem grows exponentially in the number of scenarios. If the coupling constraints are relaxed, the problem can be decomposed into K sub-problems, such that each scenario can be solved independently. Lagrangian relaxation is applied to the coupling constraints by allowing the DM to follow different decision rules in different scenarios; however, violations of the decision rule coupling constraints are penalized with Lagrange multipliers in the objective function.
By assigning Lagrange multipliers λ = (λ 1 , . . . , λ K ) to these coupling constraints, the following Lagrangian is obtained: The dualization of problems (1)-(4) with respect to the coupling constraints can be expressed as the following problem: (2) and (4) (13) By general duality theory, the optimal value of problem (13) is greater than or equal to the optimal value of problems (1)-(4). Therefore, an upper bound of problems (1)-(4) is obtained by solving problem (13). The form of D(λ) makes it suitable to be decomposed into K scenario subproblems: The sub-problems are generally of small size in comparison with problem (13), so they can be solved efficiently by applying parallel computation techniques. D(λ) can be obtained by summing up the optimums of the sub-problems with weights p k .
Theoretically, if all of the constraints are convex and all of the variables are continuous, the optimum of Z LD will be equal to the optimum of the original problems (1)-(4). A duality gap may exist, however, due to the existence of integer decision variables. This means that the optimum of Z LD will be strictly larger than the optimum of problems (1)-(4). Therefore, the aim is to obtain a good upper bound on problems (1)-(4) by solving Z LD .
A tight upper bound Z LD can be obtained by using a subgradient method to vary λ. Given an initial value λ 0 , a sequence of multiplier values λ i is generated by the following rule: where θ k i andθ i are optimal solutions to D(λ) with the multipliers set to λ k i and t i is a scalar step size in the ith iteration. For the choice of step size t i , the most commonly used strategy is adopted (Fisher 2004): where D * is the value of the best known feasible solution to problems (1)-(4), D(λ i ) is the optimal solution to the Lagrangian relaxation with multipliers set to λ i , and the scalar ρ i is chosen to be between zero and two, and its value is halved whenever D(λ i ) fails to decrease in a fixed number of iterations. A heuristic method can be used to generate feasible solutions to Z LD , which are lower bounds. This method should be designed based on the characteristics of the problem. For example, in the application of the hybrid WTE problem described in Section 4, the heuristic method used is: Set the maximum capacities x M and y M as the maximum values of x k M and y k M among the K sub-problems, respectively; set the initial capacities ε, ε , and other decision rule variables α, β, γ 1 , γ 2 , γ 3 as the average value of ε k , ε k , α k , β k , γ 1 k , γ 2 k , γ 3 k , respectively. The gap between the upper bound and lower bound determines the quality of the solution. As solving the dual to optimality is not guaranteed, the search is terminated when the gap is lower than a certain value or after a predetermined number of iterations.

Case study: Hybrid WTE system with anaerobic digestion and gasifier technologies
This section presents a comparative analysis of the standard ROA and the proposed multistage stochastic programming approaches with decision rules to the design and deployment of a hybrid WTE system. The goal of this case study is to enable a comparison of both approaches and evaluate their similarities and differences. This case study demonstrates that the proposed approach lends itself naturally to the analysis of more complex systems design problems with multiple sources of uncertainty and decision variables, while alleviating some of the shortcomings of standard ROA.

Problem analysis
Municipal solid waste management is becoming a major issue in the sustainable development of megacities. In the face of increasing waste and limited landfill sites, WTE technologies, which generate energy in the form of electricity and/or heat from waste, are enjoying high favor due to their capacity to recover energy while efficiently disposing of waste. So far, the literature on WTE systems has mainly focused on system optimization and evaluation. Little work has considered the strategic management process of designing and deploying WTE system capacity in an urban environment, especially from the perspective of flexibility. This case study targets the issues of how to manage the decisionmaking process during the implementation of such flexible systems by applying the proposed approach. An upcoming hybrid WTE system with Anaerobic Digestion (AD) and gasifier technology that treat organic waste in Singapore is considered. In the hybrid WTE system, there are two types of technologies to treat organic waste: AD and gasification. Both technologies generate gas from organic waste, and then generate electricity using a gas engine. The AD unit is better at treating high-moisture-content organic waste, such as food waste, whereas the gasifier is better at treating lowmoisture-content organic waste such as paper, wood, and horticultural waste. In this hybrid system, all food waste is treated by AD, and other organic waste is mainly treated by the gasifier. If the capacity of the gasifier is not high enough to handle all other organic waste, the untreated feedstock is transferred to the AD.
This study builds upon and extends the work by Hu and Cardin (2015). Those authors considered embedding flexibility into the system design as a mechanism to deal pro-actively with uncertainty in the amount of organic waste generated and collected. After analyzing the uncertainties and interdependencies between the components of the systems in AD plants, two important components related to capacity were identified as valuable opportunities to embed flexibility (i.e., tipping floor and major tankage). Using these mechanical design enablers, it is assumed that the AD unit is designed with the ability to expand AD capacity, x t , in a flexible manner, as needed. More specifically, one assumes that the AD unit is modularly designed, which means that it can expand the capacity by units of modules, with o u denoting the capacity of an AD module. Exploiting a similar capacity expansion strategy, the capacity of the gasifier unit y t is considered as a continuous decision variable; that is, it can be adjusted approximately continuously. The goal of DMs is to select the appropriate capacity expansion strategy to maximize the expected total profit-that is, Expected Net Present Value (ENPV)-over a life span T = 9 years.
It is assumed that the centralized WTE system is built in an area of Singapore called Tuas, using trucks to collect organic waste from households and food courts all around the city. Denote ξ t as the amount of food waste and η t as the total amount of other organic waste collected in year t. Historical data in the past 10 years show that the amount of food waste and other organic waste collected has a deterministic positive growth rate, but the amount of variation in each period is random. Based on this observation, it is suitable to model the fluctuating waste amount using a standard geometric Brownian motion formulation; that is, dξ t = μξ t dt + σ ξ t dz t . Variable ξ t represents the uncertainty in either food waste or other organic waste (η t ); ξ 0 = 191 per ton per day (tpd) and η 0 = 2823 tpd are the waste amount in year 0 for food and other waste amounts, respectively; μ F = 14.1% and μ O = 6.0% are the mean annual growth rate for food and other waste amounts respectively; σ F = 16.4% and σ O = 4.1% are the volatility of the food and other waste amounts, respectively; and dz t is the basic Wiener process giving a random shock to μ.
The organic waste is turned into biogas, compost, and residue after a series of processes at the plant. The biogas is further used to generate electricity, the compost can be sold for agriculture, and the residue and any unprocessed feedstock is sent to a landfill by paying a disposal fee. The revenue from the WTE system mainly comes from three sources: the sale of electricity to the grid, the sale of compost, and the tipping fee collected from the households and food courts. The cost of the WTE system consists of the installation cost, maintenance cost, land rental cost, labor cost, transportation cost, as well as disposal cost. The profit function r t in year t can be expressed as where T t is the tipping fee; R At and R Gt are the revenue generated by the AD and gasifier units, respectively; C At and C Gt are the cost of AD and gasifier, respectively; and P t is the penalty function incurred by capacity shortage in AD. The definitions of the functions are as below, and technical parameter values from existing WTE technologies are based on (Klein 2002): where x M is the designed maximum number of AD modules, f t = max(0, x t o u − ξ t ) represents AD capacity shortage, p = 70% is the purity ratio of the food waste, and τ = 10% is the residue rate of AD. Parameter z 1 = S$22 469 tpd is the revenue from electricity generation of each ton per day waste treated in the AD, z 2 = S$1336/tpd is the revenue from compost sale of each ton per day waste treated in the AD, z 3 = S$28 105/tpd is the unit revenue from tipping fee, z 4 = S$700/tpd is the unit transportation cost, z 5 = S$75 000/tpd is unit capacity installation cost of the AD, z 6 = S$816/tpd is the unit land rental fee, z 7 = S$204/tpd is the unit reserved land fee, z 8 = S$675/tpd is the unit labor cost of the AD, z 9 = S$225/tpd is the unit maintenance cost of the AD, and z 10 = S$28,105/tpd is the unit disposal cost. Parameter p = 70% is the purity ratio of other organic waste, f t is the capacity shortage for the gasifier, y M is the maximum capacity reserved for the gasifier, z 11 = S$62 678/tpd is the revenue from electricity sale of each ton per day waste treated in the gasifier, z 12 = S$5840/tpd is the total of labor, admin, maintenance cost per ton per day waste treated in the gasifier, z 13 = S$2920/tpd is the cost of refuse-derived fuel process per ton per day waste treated in the gasifier, z 14 = S$96 970/tpd is the capital cost of per ton per day capacity of the gasifier, and τ , = 20% is the residue ratio of the gasifier.

Real options analysis
Based on the binomial lattice approach by Cox et al. (1979), Kamrad and Ritchken (1991) developed a multinomial lattice procedure to value options with multiple sources of uncertainty. In a lattice with two sources of uncertainty, there are five possible movements emerging from each node, each with a probability of occurrence. Similar to the binomial lattice approach, the multinomial lattice is solved using a backward induction process. At the end of each period t, DMs must decide on capacities for the next period x t+1 and y t+1 that will maximize the total expected profit based on the current state of the system, x t and y t , and information on uncertainty realization, ξ t and η t . The main challenge of lattice-based methods, however, is that the computational effort increases dramatically when the dimension of the state variables increases, which makes them difficult to apply when there are more than two state variables. Calculations using standard dynamic programming show that the hybrid WTE problem cannot be solved with two dimensions of uncertainty sources and two flexibility strategies.
To get results from the multinomial lattice model, an approach based approximate dynamic programming is used to solve the problem. ADP is a modeling and algorithmic framework for solving stochastic optimization problems with multidimensional random variables and multidimensional decision variables. The key idea of ADP is to use a value function to approximate the reward function R t+1 (·) defined in the equation: where n = 5 is the number of jumps from each node for a multinomial lattice with two sources of uncertainty, π = 8% is the discount rate, p j represents the probability of each emerging path from a node, and S t, j represents the state corresponding to the emerged path j.
Using an approximation to the value function, the decisionmaking process can be simulated forward in time by generating samples of uncertainty paths. The approximation can then be updated iteratively using regression methods based on observations in the simulation. Linear functions are widely used in the literature, as they are easy to implement, and powerful to approximate complex value functions. In this article, R t+1 (·) is approximated using a weighted linear combination of a set of basis functions φ t+1 (·), associated with weights ρ t+1 ; that is, . Details of the least-square policy iteration (LSPI) algorithm used to solve the hybrid WTE system problem are provided in Appendix A.

Multistage stochastic programming with decision rules
To apply the decision rule approach to analyze the hybrid WTE system, two types of decision rules are used. For the AD unit, a set of conditional-go decision rules are considered. The type of decision rule can be expressed as: At the beginning of year t, if the observed amount of waste collected in the last year is more than a certain level (i.e., (x t−1 − α t )o u ) then expand the capacity by β t o u . Here, x t is an integer representing the number of AD modules installed in period t, α t is a parameter in unit of number of modules representing the severity level of the capacity shortage in year t (i.e., a higher (lower) α t means the DM is more (less) keen to expand AD capacity to prevent capacity shortage), and β t represents the scale of each expansion in year t (i.e., the number of AD modules of capacity o u = 30 tpd deployed). To fit to the framework in Section 3.1, constraint (4) about the decision rule related to AD expansion can be represented by (41) is a MILP. The number of constraints and decision variables increases exponentially as the number of scenarios considered increases. To solve it efficiently, Lagrangian decomposition is used by dualizing the coupling constraint (38) with a set of Lagrange multipliers. The problem is then decomposed into K scenario sub-problems. The algorithm described in Section 3.3 is applied and implemented in C++ by calling CPLEX.

... Numerical results
This section presents the numerical results of the standard ROA using ADP and proposed decision rule approaches. The LSPI algorithm described in Appendix A was programmed using C++. The number of sample paths generated in each iteration was I = 40. The results converged well after H = 1000 iterations.
The multistage stochastic problem (41) was solved under 512 scenarios of food waste amount and other organic waste amount generated via the Monte Carlo approach, respectively (the number of scenarios used is the same as that of a binomial lattice with nine time periods). The algorithm terminated when the gap between the upper and lower bound was smaller than 1% or after 100 iterations. The program was run on a PC with Intel Core i5-2500 @ 3.30 GHz quad core processor and 8 GB memory. Solving the problem with T = 9 year periods shows that the algorithm is able to get feasible solutions with duality gap of 2.1% in about 333 seconds.
To calculate the expected value of flexibility, a baseline inflexible design problem is formulated as where x f , y f are the capacity of the AD unit and gasification unit, respectively, throughout the life span. One should note that building a large rigid system can benefit from economies of scale; that is, the average unit capacity installation cost is lower. In this study, the unit capacity installation cost of the inflexible design is assumed to be 90% of that of the flexible designs. Table 1 shows the solutions for the three models. For the inflexible design, ENPV = S$751 100 000. For standard ROA using the multinomial lattice and ADP, ENPV = S$765 200 000. As for the decision rule model, ENPV = S$765 100 000. As can be seen, the decision rule approach obtains a performance for the flexible system similar to the standard ROA using multinomial lattice and ADP-only a 0.01% difference. Both the decision rule and ROA/ADP approaches highlight the additional value from flexibility in approximately the same amount, as compared to the inflexible system (i.e., benchmark expected value of flexibility is S$765 200 000 -S$751 100 000 = S$14 100 000 and about S$14 000 000 under the decision rule approach). Both approaches can solve more complex problems with multiple uncertainty sources and multiple flexibility strategies. The computation is much faster for the decision rule approach (333 seconds) compared with the standard ROA (4254 seconds).

... Using the solution
This section describes the forms of the solutions under each approach and how to use them to exercise the flexibilities in operations. The comparison aims to support the view that the solution from the decision rule approach is more intuitive to use in operations than from the standard ROA.
In the standard ROA using a multinomial lattice and ADP, a linear function of the state variables is used to approximate the value functions R t (·), and the results are the exact forms of R t (·) at each time period. Based on the forms of R t (·), decisions on the respective capacities for AD and gasifier can be made by solving the optimization problem at each time period t: subject to capacity currently installed x t−1 and y t−1 and uncertainty realization for variables ξ t−1 and η t−1 . There is no concise set of decision variable values that can adequately capture the solution for all time periods. Under the decision rule approach, the AD should be designed with initial capacity of ε = 11 modules (330 tpd) and extra land should be reserved to facilitate a maximum capacity of x M = 1470 tpd in the future. The gasifier should be designed with an initial capacity of ε = 3219 tpd with a maximum capacity of y M = 6382 tpd (values can be rounded for practical considerations). At the end of every year t − 1 ( t ≥ 2), the DM should decide on the capacity of the AD unit for year t based on the amount of food waste collected in year t − 1: if it is higher than the value of the current capacity minus α = 2 modules of capacity (i.e., x t−1 −2 × 30 tpd), then the capacity should be expanded by β = 3 modules (i.e., 90 tpd). At the same time, the DM should adjust the capacity of the gasifier starting from year 4: set the capacity of the gasifier as described by the linear function of the amount of other organic waste over the last three consecutive years; that is, y t = 0.47η t−1 + 0.28η t−2 + 0.35η t−3 .
The decision rule approach provides readily applicable guidance on the decision-making process in operations at each time period. The capacity to deploy for both plants can be determined straightforwardly by following the decision rules, subject to realizations of the two uncertainty drivers. In contrast, the form of the solution under ADP requires solving the optimization program (43) at each time period to determine the capacity for both plants. This is typical of the form of solutions obtained from multinomial and ADP approaches. It is obvious that this form of solutions may be more challenging to use for DMs as compared with the proposed decision rule approach.

Comparison with standard ROA based on a binomial lattice
In the above analysis, a standard ROA based on a multinomial lattice and dynamic programming is unable to solve the complex hybrid WTE problem, due to the curse of dimensionality. It is valuable, nonetheless, to gain further knowledge by comparing the proposed decision rule approach and standard ROA based on dynamic programming. Therefore, in this section, a simplified version of the WTE system is analyzed. In the simplified WTE system, the gasifier is suppressed; that is, it degenerates into an AD plant only. This is a special case of the hybrid WTE system when the capacity of the gasifier is set to zero, and no other organic waste is collected except for food waste. In such a system, there is only one uncertainty source (amount of food waste collected) and one flexibility strategy (expansion of AD capacity). A standard ROA approach based on a binomial lattice and adapted from Nembhard et al. (2005) can be used to analyze this system, where uncertainty is modeled using a binomial lattice. The binomial lattice problem is solved via backward induction to choose optimal options to maximize the total expected profit based on the current state and the realization of the uncertain data in each time period. The results show that both approaches also recognize a similar amount of additional value from flexibility. Standard ROA values the flexible system at ENPV = S$40 800 000, whereas the proposed approach values it at ENPV = S$38 700 000-a 5.1% difference. An out-of-sample study was conducted to numerically validate the performance of the two approaches by generating 10 000 sample scenarios and simulating the decision-making process. The obtained results show that standard ROA values the flexible system at ENPV = S$34 800 000 as compared with the proposed approach that values it at ENPV = S$34 300 000 -a 1.4% difference. In terms of computational performance, the standard ROA is much faster in the optimization part (0.015 seconds) than the proposed approach (73 seconds). It is of similar speed in the out-of-sample study, however, as compared to the proposed approach (0.02 seconds).
The solutions of the standard ROA requires conducting a backward induction process at each node, which can be stored in a lookup table showing the optimal capacity deployment for the next year based on the current state. Theoretically, at the end of each year, DMs can choose the optimal capacity for next year by referring to the lookup table for that node based on the currently installed capacity. One practical disadvantage of using the lookup tables to make decisions, however, is that it does not provide clear guidance when the realization of uncertainty falls outside the values of the lookup nodes in the lattice. One way to approximate the problem is to locate the realization of uncertainty to the nearest corresponding node in the lattice (in fact, the out-of-sample simulation was conducted using this method). Of course, the best way to use standard ROA is to conduct the backward induction process in each year based on realized uncertainty information to find the optimal decisions. This, however, increases computational effort, which is exacerbated in a multinomial setting. Therefore, the correctness of the decision depends on the ability to find the appropriate position in the lattice and correctly apply the recursive process. Comparing to the decision-making process under the decision rule approach described in Section 4.4, one may observe that it is greatly simplified under the proposed approach and intuitive to use as compared with standard ROA.

Discussion and conclusion
This article proposes a stochastic programming approach based on decision rules to analyze flexibility in the design and management of complex engineering systems, with demonstration applications in WTE system technology. The case study demonstrates that the approach is readily applicable to analyze engineering systems with multiple uncertainty sources and flexibility strategies, whereas the standard ROA proves more challenging, due to the curse of dimensionality. The numerical results show that the decision rule approach recognizes a similar amount of additional value stemming from flexibility in terms of ENPV, provides the same rank ordering (i.e., inflexible versus flexible solutions), and requires less computational time as compared with standard ROA based on a multinomial lattice and ADP. An important difference is that the form of the solution under the proposed approach is quite intuitive to use by DMs in operations as compared with standard ROA/ADP techniques. The analysis of a special case of the WTE system shows that similar conclusions can be drawn when comparing the proposed approach to standard ROA based on binomial lattice analysis.

Comparisons
In terms of life cycle performance, both approaches recognize additional value with a similar magnitude stemming from flexibility. This is true in the hybrid WTE system, where the standard ROA based on a multinomial lattice and ADP values the flexible system at ENPV = S$765 200 000, whereas the proposed approach using decision rules values it at ENPV = S$765 100 000-a 0.01% difference. This is also true in the simplified AD plant, where the proposed approach values the flexible system in an out-of-sample study with only a 1.4% difference as compared with the standard ROA approach based on a binomial lattice analysis.
In terms of computational performance, the decision rule approach (333 seconds) is much faster in reaching a solution as compared with the standard ROA using a multinomial lattice and ADP (4254 seconds). Even though in the simplified AD plant analysis the standard ROA based on binomial lattice is much faster than the proposed method in the optimization part (0.015 seconds as opposed to 73 seconds), it is of similar speed in the out-of-sample study as compared with the proposed approach (0.02 seconds).

Limitations and future work
Many opportunities exist for future work to address the limitations of this work. This article only investigates one type of flexibility strategy (i.e., capacity expansion); however, it is clear that many more can be considered and studied (e.g., temporary or permanent shutdown, switching technologies, investment deferral, etc.). Those can be analyzed using the proposed approach, while also considering multiple facilities for deployment in time and space. Investigating other formulations for the decision rule is also possible (e.g., enforcing a number of consecutive years of increase or decrease before making a decision). Although small percentage differences exist between the performance of the flexible systems as recognized under standard ROA and the proposed decision rule approaches, it is possible that decision rules formulated differently could generate better value and performance. After analyzing the forms of decision rules and types of flexibilities, a more generic framework can be developed to guide the choice of decision rules when analyzing different systems and to conduct flexibility analysis in an engineering setting.

Notes on contributors
Michel-Alexandre Cardin is Assistant Professor of Industrial and Systems Engineering (ISE) at the National University of Singapore (NUS) and Lead of the Strategic Engineering Laboratory. He is a research affiliate at the Massachusetts Institute of Technology (MIT) in the United States, principal investigator at the Singapore-ETH Centre Future Resilient Systems project, and has led several projects at the Singapore-MIT Alliance for Research and Technology. His research interests include development, empirical evaluation, and applications of novel methodologies to design and architect complex engineering systems for uncertainty and flexibility-also known as real options. Applications focus on infrastructure systems in domains such as aerospace, defense, energy, real estate, oil and gas, transportation, and water management. He is an Associate Editor of the INCOSE Journal of Systems Engineering, a member of the Editorial Review Board for IEEE Transactions on Engineering Management, and founding chairman of the organizing committee for the conference on Complex Systems Design and Management (CSD&M) Asia. He received a Bachelor's (Hons) degree in Physics from McGill University, a Master of Applied Science degree in Aerospace Science and Engineering from the University of Toronto, a Master of Science degree in Technology and Policy as well as a Ph.D. in Engineering Systems from MIT. He is also a graduate from the Space Science Program at the International Space University.
Qihui Xie is a Ph.D. candidate in the Department of Industrial and Systems Engineering at the National University of Singapore. He received a Bachelor's degree in Industrial Engineering from Shanghai Jiao Tong University in 2008. His current research interests include flexibility and real options analysis, optimization under uncertainty, and engineering systems design.
Tsan Sheng Ng is an Associate Professor in the Department of Industrial and Systems Engineering, National University of Singapore. He received his Ph.D. (ISE) at NUS in 2005. His research interests are in optimization modeling and applications in energy, logistics, and sustainability problems under risk and ambiguity. His work includes robust optimization modeling approaches in energy-economic resilience problems, energy recovery facility location systems, and natural gas production planning problems. He has published in international peer reviewed journals such as Operations Research, Naval Research Logistics, IIE Transactions, European Journal of Operational Research, IEEE Transactions in Reliability, IEEE Transactions in Control Systems Technology, Energy Economics, etc. He has collaborated on industrial research projects with global oil and energy providers and also co-led nationally funded projects in energy and sustainability solutions for megacities, resilience of critical infrastructure, energy security, and robust optimization. He has also conducted training courses in systems modeling for policy-makers and government agencies and supervised numerous students and teams in implementing process improvements projects within multi-national companies and also small-and medium-sized enterprises based in Singapore. He is an Editorial Board member of the International Journal of Automation and Logistics.