Dynamic Operability Analysis of Process Supply Chains for Forest

May 28, 2014 - volatile market conditions are becoming the norm in the chemical process ... strategy to improve the struggling business model entails ...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/IECR

Dynamic Operability Analysis of Process Supply Chains for Forest Industry Transformation Richard Mastragostino and Christopher L. E. Swartz* Department of Chemical Engineering, McMaster University, 1280 Main Street West, Hamilton, Ontario, Canada L8S 4L7 ABSTRACT: An important attribute of a supply chain in a competitive and volatile market environment is the ability to respond rapidly to demand variation. A particularly relevant application is the forest products industry, where a promising strategy to improve the struggling business model entails the shift from commodity products toward high-value specialty products. A key implication is that new process and supply chain designs have sufficient capability to respond quickly to market changes, such that product availability is high. In this study, we develop a computational framework for dynamic operability analysis of process supply chains. A dynamic model of a multiproduct, multiechelon system supply chain system is developed, and incorporated within an optimization framework. A two-stage stochastic programming approach is applied for the treatment of demand uncertainty. A bicriterion optimization problem is formulated for generating the Pareto frontier between an economic and responsiveness criterion. Two case studies are presented to demonstrate the applicability of this framework.

1. INTRODUCTION A supply chain (SC) is a network of facilities that performs the functions of raw material procurement, raw material transformation into intermediate and finished products, and distribution of products to customers.1 Supply chains are traditionally characterized by a forward flow of material, and a backward flow of information in the form of demand and orders. Today’s global companies are becoming ever more aware that they no longer compete as stand-alone entities, but rather as supply chains.2 Process supply chains face an increasing array of challenges, and key to economic sustainability is adequate robustness, since turbulent and volatile market conditions are becoming the norm in the chemical process industry. Flexibility and responsiveness are regarded as two significant attributes of a supply chain for mitigating against uncertainty, particularly when product availability is a fundamental requirement to remain competitive.2,3 The motivating context for this study is the forest products industry (FPI), which has traditionally been primarily commodity based and is facing numerous challenges as a result of global low-cost competition and declining demand. To remain financially viable in the face of low-cost competitors, operating costs, capital costs, and R&D activities have been cut, leading to obsolete mills and operating policies.4 A possible strategy to improve the struggling business model entails revenue diversification of existing pulp and paper mills, through the production of low-volume, high-value specialty products alongside conventional pulp and paper products.4 Various process technologies can be incorporated into existing pulp and paper mills to capitalize on the components in woody biomass under-utilized (e.g., hemicellulose and lignin), effectively transforming them into integrated forest biorefineries (IFBR).4,5 A significant implication of the IFBR paradigm is that new process and supply chain designs should have sufficient resilience to market volatility, feedstock variation, and other sources of uncertainty. Cost reduction has been the major focus to improve profitability in the supply of commodity type © 2014 American Chemical Society

pulp and paper products normally characterized by stable demand; however, this strategy is the incorrect one for the supply of so-called “innovative” (specialty) products where demand is typically highly volatile, and responsiveness should be the major focus.6 Research that addresses the role and properties of the supply chain within this new IFBR paradigm is limited, with some notable recent advances. Mansoornejad et al.7 propose a methodology that links product portfolio selection and supply chain analysis to investigate how manufacturing flexibility can be exploited in a margins-based supply chain operating policy to secure a company’s profitability. Mansoornejad et al.8,9 define two flexibility metrics for supply chain analysis and illustrate their use for analysis of supply chain flexibility under various conditions through application to biorefinery case studies. Supply chain responsiveness has not been directly addressed within the IFBR context. Analogous to supply chains, process plants operate in continually changing environments and are required to satisfactorily attenuate disturbances or transition between operating points rapidly to exploit external opportunities. This may be achieved through an adequate level of operability, defined as the ability of a process to perform satisfactorily under conditions away from the nominal operating or design conditions.10 Two major facets of operability are steady-state flexibility, which reflects a parameter range over which feasible steady-state operation can be attained; and dynamic operability, which indicates the dynamic responsiveness of a system. Both of these measures are in general strongly affected by the process design, which has motivated the development of rigorous and systematic approaches for addressing operability considerations during process design.10,11 Dynamic optimization has been shown to be an effective framework for the inclusion of Received: Revised: Accepted: Published: 9825

February 11, 2014 April 13, 2014 May 9, 2014 May 28, 2014 dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

of both profitability and responsiveness criteria through a bicriterion optimization approach. You and Grossmann23 extend the prior study by considering demand uncertainty and a probabilistic model for stock-out that is dependent on safety stock levels at distribution centers, treated as optimization variables. The expected lead time is the quantitative measure of responsiveness evaluated, linked to the stock-out probability, and a similar definition of worst case lead time as in previous work. The expression for the production delay in multiproduct plants differs depending on a priori specification of large or small changes in a product demand, to account for a different production sequencing policy under high or low demand. A hierarchical algorithm is proposed on the basis of the decoupling of strategic and operational decisions, to efficiently solve the large-scale MINLP problem. You and Grossmann24 further extend this work to address the inventory management of a single-product, multiechelon inventory system, where the maximum (worst case) guaranteed service time (GST) of the last echelon of the supply chain is proposed as a measure of responsiveness. The GST and processing time of the echelon, and the GST of the directly preceding echelon contributes to the net lead time, denoting the time span over which safety stock coverage against demand variation is required. The net lead time is related to the safety stock level through a specified safety stock factor. In this paper we have extended these prior works to evaluate the lead time of a multiproduct, multiechelon inventory supply chain system through a dynamic supply chain optimization formulation that directly links production decisions, transportation quantities, inventory levels and demand. Recourse actions for different demand scenarios are considered through a two-stage stochastic programming formulation. An economic criterion is evaluated, dependent on procurement, production, inventory, and transportation costs. The dynamic operability analysis framework is formulated as a bicriterion, two-stage stochastic, MILP problem, and the Pareto frontier between economic performance and lead time is shown for two case studies within the biorefinery application context. The paper is organized as follows. Section 2 describes the supply chain system considered, with the mathematical formulation of the dynamic supply chain model presented in section 3. Section 4 develops the dynamic operability analysis framework, outlining the methodology for evaluating responsiveness. Section 5 presents case studies to demonstrate the application of the approach and its use to investigate the effect of process supply chain design on responsiveness. Finally, section 6 concludes with some remarks and outlines future research avenues.

dynamic performance considerations in process design formulations.12−17 Typically, an economic objective is maximized or minimized, subject to constraints on the plant’s dynamic response to a specified set of disturbances. This paper will draw upon such formulations for dynamic operability analysis of process supply chains. Notions of operability appear in the supply chain literature largely through the concepts of agility and flexibility and include metrics such as on-time deliveries, back orders, manufacturing lead time, and fill rate. Approaches include conceptual analyses, quantitative metrics, and mathematical programming formulations. Christopher2 suggests that the key to survival in today’s challenging business environment is through agility, where significant routes to agility are improved market sensitivity (i.e., responding to real demand), leveraging information, collaboration between buyers and suppliers, and more efficient coordination and management of the supply chain. Naylor et al.18 discuss the lean and agile manufacturing paradigms within the supply chain context, as well as their integrated application to different segments of a supply chain. Backx et al.19 identify agile, lean, responsive, smooth, and predictive as key objectives of an operable supply chain, and point out the need to cast the assessment in terms of a multicriterion framework with ample consideration of the significance (weight) attributed to each objective. Furthermore, the authors consider the volatility surrounding the plant as an opportunity to be exploited, rather than a disturbance to be compensated for, and propose a strategy for intentionally dynamic operating policies at the manufacturing level. Holweg20 outlines three key dimensions, i.e., product, process, and volume, which influence responsiveness in a supply chain. A conceptual model is proposed to allow enterprises to align their supply chain strategy to the identified factors in order to find a suitable balance between responsiveness and supply chain efficiency. Furthermore, as pointed out in Holweg,20 responsiveness is strongly influenced by manufacturing flexibility, yet the link has not been rigorously explored. A comprehensive review of supply chain performance measures is given in Beamon,21 where performance measures are classified into three categories: resources, output, and flexibility. The first is associated with economics and resource usage, the second includes customer responsiveness and the quality and quantity of product produced, and flexibility includes measures that reflect the ability to accommodate volume and schedule fluctuations from suppliers, manufacturers and customers. Product differentiation postponement is a related concept that involves the practice of delaying the point of differentiation, with the aim of carrying inventory in a generic form, and completing final production/customization once true demand is known. There are several advantages of postponement, including reduction of finished product inventory that is normally more expensive, improving responsiveness, and enabling the ability for “mass customization”.2,6 A few recent studies have addressed responsiveness of process supply chains within mathematical programming frameworks. You and Grossmann22 propose worst case lead time, i.e., length of the longest time path of chemical flow from a supplier to a customer by way of several manufacturing sites assuming zero inventory, as a quantitative measure of supply chain responsiveness. Furthermore, the authors present a framework for strategic decision making, accounting for the selection of suppliers, manufacturing sites and process technology, and operational scheduling, under consideration

2. SYSTEM DESCRIPTION Figure 1 illustrates the supply chain structure considered in the present study. Arrows connecting adjacent echelons illustrate transportation routes for raw material and products. The supply chain comprises suppliers (S), integrated forest biorefineries (IFBR), regional distribution centers (RDC), and markets (M). The model is general in terms of the number of nodes in each echelon, as well as transportation links from one echelon to the next. A raw material order of ORr,h,b,t units is made to supplier h ∈ HRr from biorefinery b ∈ BRr, where HRr is the set of suppliers of raw material r, and BRr is the set of biorefineries which consume raw material r. Similarly, an order of OPp,h,b,t units of product is made to supplier h ∈ HPp from biorefinery b ∈ BPp, 9826

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

production scheduling), based from the state-task-network (STN) formulation proposed in Shah al.26 The supply chain model developed represents a discrete-time formulation, where the horizon is divided into equal length intervals, ΔT, indexed by t = {1, ..., H}, where H is the horizon length. Three types of constraints are included in the supply chain model proposed: mass balance, assignment of production tasks, and capacity constraints. 3.1. Integrated Forest Biorefinery (IFBR). The mass balance of raw material r at biorefinery b ∈ BRr at each time period is given by (1), Figure 1. Diagram of process supply chain system.



OrR, h , b , t − (δ O / ΔT ) + h,b

h ∈ HR r

where HPp is the set of suppliers of product p and BPp is the set of biorefineries that consume or produce product p. This allows a product p to be procured from a supplier or produced at the biorefinery (intermediate product). A delivery delay accounts for the time between when an order is placed to a supplier and when the material arrives on site and is ready to be processed at the biorefinery. Raw material is stored with inventory of IRr,b,t units. A quantity is withdrawn from inventory and converted into Pi,b,t units of main product, through production task i. Biorefineries are multiproduct manufacturing sites, where Ib is the set of production tasks (schemes) assigned to process units at biorefinery b. A process unit k at biorefinery b can be continuous, belonging to the set KCb , or batch, belonging to the set KBb . It is worth clarifying that a process unit can represent a distinct unit operation (performed at a piece of equipment), or more generally a collection of unit operations (section of the biorefinery). This allows sufficient generality in the definition of a “production task”. For all production tasks, mass balances are expressed linearly in terms of the quantity of main product produced. Furthermore, a manufacturing delay is permitted for each production task. Product at each biorefinery is stored with inventory of IPp,b,t units. A quantity of FPp,b,v,t units is withdrawn from inventory and shipped from biorefinery b ∈ BPp to regional distribution center v ∈ VPp, where VPp is the set of regional distribution centers that store and distribute product p. A delay is associated with material shipment between biorefineries and regional distribution centers, which accounts for the transportation and material handling delay. Products distributed from a regional distribution center v are stored with S inventory of IDP p,v,t units. A quantity of Fp,v,t units is withdrawn from inventory and distributed to regional markets to satisfy demand and accumulated back orders. If demand is not met at time period t, a back order of Bp,v,t units has accumulated. The current work does not address different transportation modes, and furthermore we have assumed that sufficient transportation resources (e.g., pipelines, barges, rail cars, vessels) can be acquired to meet any transportation requirement (i.e., infinite capacity). Furthermore, customers that demand the same product, from the same regional distribution center have been grouped into a single regional market for simplicity. Demand for product p at each regional distribution center v is a stochastic parameter.

= IrR, b , t

∑ i ∈ IrCR ,b

ηrR, iPi , b , t + (1 − αr )IrR, b , t − 1

∀ r , b ∈ BR r , t

(1)

δoh,b

where represents the delivery delay between when an order is made to a supplier h and when the material arrives on site and can be processed by a process unit at biorefinery b. Pi,b,t represents the quantity of main product produced by production task i at biorefinery b, beginning at time period t. ηRr,i is the mass balance coefficient of raw material r of production task i, representing the relative proportion of raw material consumed (negative) per unit of main product produced, and ICR r,b represents the set of production tasks that consume raw material r at biorefinery b. The deterioration percentage of raw material inventory during the duration of a time period is modeled using a constant deterioration fraction, αr, as given in Ekşioǧlu et al.25 Raw material purchased from suppliers during a time period is restricted to a maximum quantity as given by (2). The inventory set point of raw material (SRr,b) is equated to the inventory position at t = 0, as given by (3). Additionally, an end point constraint, given by (4), enforces held inventory to be greater than or equal to the inventory set point at the final time period. OrR, h , b , t ≤ λrR, h

∀ r , b ∈ BR r , h ∈ HR r , t

(2)

IrR, b , t = SrR, b

∀ r , b ∈ BR r , t = 0

(3)

IrR, b , t ≥ SrR, b

∀ r , b ∈ BR r , t = H

(4)

Similarly, the mass balance of product p at biorefinery b ∈ BPp at each time period is given by (5),



OpP, h , b , t − (δ O / ΔT ) + h,b

h ∈ HPp

+

∑ i ∈ I pM, b

∑ i ∈ I pC, b

ηpP, iPi , b , t

ηpP, iPi , b , t − (δiM / ΔT ) −

= I pP, b , t ∀ p , b ∈ BPp , t



FpP, b , v , t + I pP, b , t − 1

v ∈ VPp

(5)

δM i ηPp,i

where represents the manufacturing delay of production task i. is the mass balance coefficient of product p of production task i, representing the relative proportion of product produced (positive) or consumed (negative) per unit of main product produced. ηPp,i is 1 if product p is the main product of production task i. The production terms in the mass balance are separated by tasks that consume product p, belonging to the set ICp,b, and tasks that produce product p, belonging to the set IM p,b, similar to as proposed in You and Grossmann.23

3. PROCESS SUPPLY CHAIN MODEL A tactical-level supply chain model with operational elements is presented in this section, with aspects adapted from models developed in You and Grossmann23 and Ekşioǧlu et al.25 However, the model includes constraints for considering the assignment of production tasks to process units (i.e., 9827

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

Pi , b , t ≤ γi , bXi , b , t

The quantity of product purchased from suppliers during a time period is restricted to a maximum quantity as given by (6), and similarly (7) and (8) equate the inventory position to the set point at t = 0 and enforce the inventory level to be greater than or equal to the set point at the final time period. OpP, h , b , t ≤ λpP, h

∀ p , b ∈ BPp , h ∈ HPp , t

I pP, b , t



SpP, b

∀ p , b ∈ BPp , t = 0 ∀ p , b ∈ BPp , t = H



(6) (7)

DP FpP, b , v , t − (δ D / ΔT ) − F pS, v , t + I pDP ,v ,t−1 = Ip,v ,t b,v

∀ p , v ∈ VPp , t

delay from regional distribution center v to satisfy demand and accumulated back orders. Additionally, (13) and (14) are included in the model and have the same purpose as previously defined for raw material and product inventory at biorefineries.

(8)



i ∈ IR k , b

t ′= t

Xi , b , t ≤ 1

i ∈ IR k , b

Xi , b , t ≤ 1

∀ p , v ∈ VPp , t = 0

(13)

DP I pDP , v , t ≥ Sp , v

∀ p , v ∈ VPp , t = H

(14)

R P DP P S OrR, h , b , t , OpP, h , b , t , Pi , b , t , IrR, b , t , I pP, b , t , I pDP , v , t , Sr , b , Sp , b , Sp , v , Fp , b , v , t , F p , v , t ≥ 0

∀ r , p , h ∈ HR r ∪ HPp , b ∈ BR r ∪ BPp , i ∈ Ib , v ∈ VPp , t

Xi , b , t ∈ {0, 1}

∀ b , k ∈ KbB , t (9)

∀ b , i ∈ Ib , t

(15) (16)

4. DYNAMIC OPERABILITY ANALYSIS FRAMEWORK (DOAF) In this section we present an optimization framework for dynamic operability assessment of process supply chains. First, metrics for economic performance and supply chain responsiveness are formulated. Thereafter, a two-stage stochastic programming formulation is presented, in which the expected supply chain operating cost is minimized subject to an upper bound on the system lead time. By variation of the bound on the responsiveness metric, a Pareto-optimal curve representing the trade-off between economic performance and supply chain responsiveness can be generated. The following problem context has been assumed: 1. Inventory set points are design variables in the optimization for hedging against demand uncertainty. 2. At t = 0: process units at the plant are idle, and the quantity of on-route shipments to biorefineries and regional distribution centers are zero. 3. The lead time definition does not address delivery delay between regional distribution centers and individual customers. 4. There is a penalization cost associated with not maintaining inventory at the set point. 5. An end point constraint enforces held inventory to be greater than or equal to the inventory set point. 6. Once demand for product p at regional distribution center v is met, it must continue to be met for the remainder of the optimization horizon. 7. Optimization of the process supply chain is fully centralized. 4.1. Economic Criterion. The operating cost of the supply chain network over the time horizon is defined as the total

3.1.2. Continuous Process. We now define the constraint for representing the assignment of production tasks to continuous process units, where KCb is the set of continuous process units at biorefinery b. For a continuous process, all unit operations are continuous over time; therefore, the manufacturing delay to produce a unit of product depends on the residence time of the process unit. In this work we have assumed this time is small relative to ΔT, therefore δM i is zero for production tasks operated at continuous process units. Equation 10 restricts the assignment of more than one production task i to process unit k ∈ KCb at each time period.



DP I pDP , v , t = Sp , v

3.3. Variable Bounds. All continuous variables in the model are non-negative, as given by (15), and binary variables are defined by (16).

M

t − (δi / ΔT ) + 1

(12)

where δDb,v is the associated shipping and material handling and FSp,v,t represents the quantity of product p shipped

3.1.1. Batch Process. The constraint for assigning production tasks to batch process units will be developed first, where KBb is the set of batch process units at biorefinery b. The basic assignment constraint has been adapted from the assignment constraint presented in Shah et al.,26 where the authors reformulated the constraint originally derived in Kondili et al.,27 to improve computational performance. The constraint, as given by (9), restricts the assignment of a production task i to process unit k ∈ KBb at biorefinery b at time period t, if another task i′ has already been assigned within the M backward interval [t − (δM i /ΔT) + 1, t], where δi is the manufacturing delay or batch duration of production task i. The binary variable Xi,b,t is 1 if production task i begins at time period t, and IRk,b is the set of production tasks that can be performed by process unit k at biorefinery b.



(11)

3.2. Regional Distribution Center. The mass balance of product p at regional distribution center v ∈ VPp at each time period is given by (12), b ∈ BPp

I pP, b , t = SpP, b

∀ b , i ∈ Ib , t

∀ b , k ∈ KbC , t (10)

It is evident that Xi,b,t is not indexed by process unit k. This decouples the production task from the process unit, as proposed in Ierapetritou and Floudas28 for improving computational performance. Therefore, if an equivalent production task i is operated by two process units in the plant, then two tasks i and i′ are defined. Equation 11 confines the production of main product from production task i at biorefinery b beginning at time period t to a maximum quantity (γi,b), if and only if the production task has been assigned to a process unit; otherwise the production is confined to zero. Ib is the set of production tasks that can be performed by process units (continuous and batch) at biorefinery b. For a production task assigned to a batch process unit, the upper bound is dependent on the maximum batch size that can be handled by the process unit, and for a continuous process unit, the upper bound is dependent on the maximum production rate handled by the process unit, multiplied by the duration of a time period (ΔT). 9828

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

in the criterion, inventory will likely be driven to zero during the optimization horizon, resulting in a poor characterization of the inventory response. The economic criterion (Z), a measure of economic performance, is given by (27),

purchasing cost, production costs, transportation cost, and inventory costs as given by (17), Coper = Cpurch + Cprod + C trans + C inven

(17)

Individual costs are computed as follows, Cpurch =

∑ [∑ ∑ ∑ t

+

r

Cprod =

b

(18)

(19)

i ∈ Ib

∑∑ ∑ ∑ t

C inven =

ωpP, hOpP, h , b , t ]

h ∈ HPp b ∈ BPp

∑ [∑ ∑ t

+

r

(20)

ωrIR IrR, b , t +

b ∈ BR r

∑ ∑ p

ωRr,h

ωpT FpP, b , v , t Tb , v

p b ∈ BPp v ∈ VPp

∑ ∑

ωpIPI pP, b , t

p b ∈ BPp

ωpIDI pDP ,v ,t]

v ∈ VPp

(21)

ωPp,h

where and are the costs per unit of raw material r and product p purchased from supplier h and ωIi is the cost per unit of main product produced by production task i. ωTp is the cost per unit, per kilometer of product p shipped from a biorefinery to a regional distribution center, and Tb,v is the distance (km) between biorefinery b and regional distribution center v. ωIR r and ωIP p are the costs per unit of raw material and product inventory at a biorefinery for the duration of a time period, and ωID p is the costs per unit of product inventory at a regional distribution center for the duration of a time period. We also define a penalty term which is included in the economic criterion function as given by (22), Pinven =

∑ φ[∑ ∑ t

+

r

b ∈ BR r

∑ ∑ p

θrR, b , t +

∑ ∑

(27)

It is worthwhile to mention that alternative economic metrics could have been used, such as profit or net present value (NPV). However, NPV would be most applicable if the formulation included supply chain capital decisions, such as whether to invest in a plant and what products to produce and where. This aspect will be considered in future work. 4.2. System Lead Time. A major aspect of this work is the development of a quantitative characterization of dynamic operability. The basis, as previously alluded to, is to evaluate the lead time to respond to a change in demand. For a linear supply chain (i.e., single transportation and processing route), assuming zero inventories, the lead time can be computed directly by summing all transportation and production delays incurred in the route. However, it is a much greater effort to compute the lead time of a more general supply chain network, because material can travel through a number of transportation and processing routes, and raw material and product can be stored. In this work, lead time is computed by implementing a sustained step change in demand at the initial time period, and tracking when demand is fulfilled and sustained at this new level. Delays and inventory are captured through the dynamic model. The inventory set points of raw material and product are time independent decision variables in the optimization formulation for hedging against demand uncertainty, equated to the inventory position at t = 0. Figure 2 illustrates a conceptual continuous-time representation of the approach, where lead time is the delay between the

∑ ∑ ∑ ωiI Pi ,b ,t t

C trans =

ωrR, hOrR, h , b , t

h ∈ HR r b ∈ BR r

∑ ∑ ∑ p

Z := Coper + Pinven

θpP, b , t

p b ∈ BPp

θpDP ,v ,t]

v ∈ VPp

(22)

where φ represents the penalty cost per unit of inventory which R P DP is below the inventory set point and θr,b,t , θp,b,t , and θv,b,t represent the quantity of inventory which drops below the inventory set point, defined through (23)−(25). These variables are non-negative as given by (26). θrR, b , t ≥ SrR, b − IrR, b , t

∀ r , b ∈ BR r , t

(23)

θpP, b , t ≥ SpP, b − I pP, b , t

∀ p , b ∈ BPp , t

(24)

DP DP θpDP , v , t ≥ Sp , v − I p , v , t

∀ p , v ∈ VPp , t

(25)

Figure 2. Basis for evaluating the lead time of process supply chain system.

θrR, b , t , θpP, b , t , θpDP ,v ,t ≥ 0 ∀ r , p , b ∈ BR r ∪ BPp , v ∈ VPp , t

(26) 29

moment the step change occurs, and the time when demand can be served (including accumulated back orders) and sustained at the new demand level for the remainder of the optimization horizon. A binary variable, Ymet p,v,t, is used in the model to indicate this occurrence by taking a value of 1 if demand (and accumulated back orders) for product p at

This construct was posed in Jung et al. for enforcing the maintenance of inventory over the optimization horizon, and in an approximate manner the nonanticipatory aspects of the twostage stochastic formulation. In actual operation, inventory is maintained to remain robust against future changes in demand that are not anticipated. Without considering this penalty term 9829

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

regional distribution center v can be served and sustained at the new demand level at time period t, and 0 otherwise. Demand for product p at regional distribution center v is given by (28), DEp , v , t = (Up , v)St

∀ p , v ∈ VPp , t

met Y pmet , v , t − Y p , v , t − 1 − mp , v , t = 0

∀ p , v ∈ VPp , t = 2, ..., H

Y pmet , v , t − mp , v , t = 0

(28)

where Up,v is the magnitude of the step change and St is a unit step function equal to 1 when time period t is greater than 0, and 0 otherwise. Up,v is a stochastic parameter described by a normal distribution, Up,v ∼ 5 (μp,v,σp,v2), as illustrated in Figure 2. The constraints applied to evaluate the lead time are given by (30)−(34), which have been adapted from a general MILP formulation given in You and Grossmann23 to represent the disjunctive logical condition with a general form ⎡ ¬Y ⎤ ⎡Y ⎤ ⎢ lb ⎥ ⎢ m ⎥ m ub ⎢Q ≤ Q ≤ Q ⎥ ∨ ⎢Q ≤ Q ≤ Q ⎥ ⎢ ⎥ ⎢ ⎥ ⎣ f (Q ) ≤ 0 ⎦ ⎣ f (Q ) ≤ 0 ⎦

Q p2, v , t



DEp , v , t Y pmet ,v ,t

Q p2, v , t ≤ MY pmet ,v ,t

Q 1p , v , t ,

Q p2, v , t

≥0

∀ p , v ∈ VPp , t

(40)

mp , v , t ≤ 1

∀ p , v ∈ VPp , t

(41)

CU

∀ p , v ∈ VPp , t

∀ p , v ∈ VPp , t ∀ p , v ∈ VPp , t

Y pmet ,v ,t = 1

LTp , v =

(31)

In this work, the disjunction formulation allows demand to be partially satisfied at time period t through the variable Q1p,v,t, and restricts the variable Q2p,v,t to be zero when Ymet p,v,t = 0, where M is 1 a sufficiently large number. When Ymet p,v,t = 1 at time period t, Qp,v,t 2 is forced to be 0, and Qp,v,t is greater than or equal to demand at time period t. The back order balance of product p at regional distribution center v ∈ VPp is given by (35),

SLT =

(35)

where Bp,v,t represents the quantity of back orders accumulated for product p at regional distribution center v. Additionally, the end point constraint given by (36) enforces back orders to be fulfilled by the end of the optimization horizon. Bp,v,t is nonnegative as given by (37).

∀ p , v ∈ VPp , t

(37)

∑ (t − 1)mp,v ,t ΔT

∀ p , v ∈ VPp

(44)

1 |V |

⎡ 1 ∑ ⎢⎢ PV | v| v ⎣

⎤ LT ∑ p,v ⎥⎥ ⎦ p ∈ PVv

(45)

4.3. Two-Stage Stochastic Optimization Formulation. A two-stage stochastic framework is proposed for the inclusion of demand uncertainty into the formulation. In two-stage stochastic programming, certain decisions are made “here-andnow” prior to the uncertain event taking place, and the remaining decisions are postponed in a “wait-and-see” manner after uncertainty has been resolved and more information is available; these correspond to first- and second-stage decision variables, respectively.30 Typically, second-stage variables can be interpreted as corrective (recourse) measures applied in response to the uncertainty realization. In our formulation, demand uncertainty, i.e., the magnitude of the step change (Up,v), is resolved after inventory set points are determined in the first stage. Furthermore, operational variables are determined in the second stage after demand uncertainty has resolved. The scenario-based approach is applied for capturing uncertainty in demand. Of significance is that this formulation has the capability to account for recourse actions for different scenarios. Operational decisions made in the second stage, which have an effect on the system lead time, are in response to how uncertainty has resolved, contributing to a more realistic characterization of responsiveness.

∀ p , v ∈ VPp , t

Bp , v , t ≥ 0

(43)

The system lead time (SLT) is given by (45), which reflects the average for all products, over all regional distribution centers, where V is the set of regional distribution centers, and PVv is the set of products distributed and stored at regional distribution center v. An equivalent weighting is assigned to each product, at each regional distribution center in the average calculation; however, the formula can be readily revised to consider a weighted average, or worst case description.

(34)

(36)

∀ p , v ∈ VPp , t = CU , ..., H

t

(32)

∀ p , v ∈ VPp , t = H

(42)

Equation 44 computes the lead time, LTp,v, of product p at distribution center v.

(30)

Bp , v , t = 0

∀ p , v ∈ VPp

t=1

(33)

∀ p , v ∈ VPp , t

DEp , v , t − F pS, v , t + Bp , v , t − 1 = Bp , v , t

(39)

mp , v , t ≥ 0

∑ mp,v ,t = 1

(29)

∀ p , v ∈ VPp , t

Q 1p , v , t ≤ DEp , v , t (1 − Y pmet ,v ,t)

∀ p , v ∈ VPp , t = 1

A tightening constraint is proposed by exploiting system knowledge, as given by (42). The parameter CU ensures demand is met by at least the time period given by CU, where CU ≤ H. If, for instance, we expect the worst case lead time of the system to be smaller than the optimization horizon length (H), CU can be defined appropriately to improve computational performance. Consequently, Ymet p,v,t is fixed as given by (43), from the tightening constraint.

where Q is a continuous decision variable, Y is false when Qlb ≤ Q ≤ Qm and true when Qm ≤ Q ≤ Qub, and f(Q) represents the implicit system inequalities of Q. F pS, v , t = Q 1p , v , t + Q p2, v , t

(38)

Ymet p,v,t

Equations 38 and 39 restrict to not change from 1 to 0, because demand fulfillment is sustained over the optimization horizon, once met. The variable mp,v,t takes a value of 1 to indicate when Ymet p,v,t has switched from 0 to 1. Furthermore, we do not need to treat mp,v,t as a binary variable, it can be relaxed to the interval [0, 1] as given by (40) and (41), which reduces the required number of binary variables in the model and improves computational performance and efficiency. 9830

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

variables d are first-stage decision variables, whereas x, u, xA, and uA are second-stage variables associated with each scenario s, as indicated by the subscript in the DOAF formulation. Note that there is not an explicit cost (first-stage cost) associated with the inventory set points (d); the cost of inventory is accounted for through second-stage variables. Equation 46d restricts d to within an upper bound (dub), where the upper bound could reflect a physical boundary in the system (e.g., warehouse space) or a boundary can be artificially specified to tighten the MILP formulation and improve computational performance. 4.4. Solution Strategy. To compute the range of ε to estimate the Pareto frontier, we first reformulate the optimization problem so that the term [SLT] appears in the objective,

The two optimization criteria derived in sections 4.1 and 4.2 are conflicting, since an increase in economic performance (i.e., holding less inventory) typically results in a decrease in system lead time (dynamic operability). To estimate the Pareto frontier, the ε-constraint method is used, where one criterion is specified as an inequality with a specified bound treated as a parameter (ε), and the optimization problem is solved for different values of ε. In this work the expected system lead time is formulated as an ε-constraint in the optimization problem. The dynamic operability analysis framework (DOAF) is given by min

d , xs , us , xsA , usA

[Z] :=∑ ρs Zs

(46a)

s

s.t. [SLT] =

∑ ρs SLTs

(46b)

s

[SLT] ≤ ε

(46c)

d ≤ d ub

(46d)

h(d , xs , us , xsA , usA , Zs , ps ) ≤ 0

(46e)

us ∈ {0, 1}q

(46f)

usA ∈ {0, 1} f

(46g)

ξ [Z] + [SLT]

(47)

where ξ is a very small value (i.e., order of 10−10) so that it becomes negligible in the objective function. Alternatively, the economic term can be removed all together. We then solve the problem, where (47) is minimized, to compute SLTlb. This represents the lower bound of ε, i.e., the best case system lead time (worst case economic criterion). The upper bound of ε represents the situation where all inventory set points are zero, i.e., dub = 0, because this corresponds to the worst case system lead time. We again solve the problem, where (47) is minimized and dub is 0, to compute SLTub. The DOAF is solved at equal spaced intervals from SLTlb to SLTub. A larger number of intervals results in a more accurate estimation of the Pareto frontier. 4.5. Scenario Generation. Demand as given by (28) is a stochastic parameter, since the step change, Up,v, is stochastic and described by a normal distribution. A Monte Carlo sampling method is applied for generating scenarios,31 combined with a statistical method to evaluate the level of accuracy of the solution (objective function value) through a confidence interval.32−34 The objective function value of the DOAF approaches a normal distribution, and represents a Monte Carlo estimator of the expected cost (economic criterion). Therefore, the Monte Carlo sampling variance estimator, for scenarios considered is given by (48),32,34

s = 1, ..., ns

where [·] represents the expectation operator, the index s represents the scenario, ρs is the probability of scenario s occurring, and ns is the number of scenarios applied for capturing uncertainty. For notational convenience we define the function h to encompass (1)−(15), (17)−(28), and (30)− (45) derived in section 3. The variables and parameters in the above formulation are defined as d = [(SrR, b ∀ r , b ∈ BR r ), (SpP, b ∀ p , b ∈ BPp), T

(SpDP , v ∀ p , v ∈ VPp)]

x = [(OrR, h , b , t ∀ r , h ∈ HR r , b ∈ BR r , t ), (OpP, h , b , t ∀ p , h ∈ HPp , b ∈ BPp , t ), (Pi , b , t ∀ b , i ∈ Ib , t ), (IrR, b , t ∀ r , b ∈ BR r , t ), (I pP, b , t ∀ p , b ∈ BPp , t ), R (I pDP , v , t ∀ p , v ∈ VPp , t ), (θr , b , t ∀ r , b ∈ BR r , t ),

ns

(θpP, b , t ∀ p , b ∈ BPp , t ), (θpDP , v , t ∀ p , v ∈ VPp , t ),

S(ns) =

T

(FpP, b , v , t ∀ p , b ∈ BPp , v ∈ VPp , t ), (F pS, v , t ∀ p , v ∈ VPp , t )]

u = [Xi , b , t ∀ b , i ∈ Ib , t ]T

∑s = 1 ([Z ] − Zs)2 ns − 1

(48)

where Zs is the economic criterion for scenario s. Once the variance is estimated, the confidence interval of the expected economic criterion is given by (49),

x A = [(mp , v , t ∀ p , v ∈ VPp , t ), (Bp , v , t ∀ p , v ∈ VPp , t ),

⎡ zα / sS(ns) ⎤ zα / sS(ns) , [Z ] + ⎢[Z ] − ⎥ ns ⎦ ns ⎣

(LTpPV, v ∀ p , v ∈ VPp), (SLT), T

(Q 1p , v , t ∀ p , v ∈ VPp , t ), (Q p2, v , t ∀ p , v ∈ VPp , t )]

(49)

where the confidence interval is symmetric and zα/2 is the standard normal deviate such that for a standard normally distributed variable, z ∼ 5 (0, 1), Pr(z < zα/2) = 1−zα/2. Alternatively, we can assume the expectation follows a tdistribution to compute a more conservative confidence interval. The confidence interval is then compared with our desired accuracy, and if unsatisfactory, more scenarios are applied to capture uncertainty.

T u A = [Y pmet , v , t ∀ p , v ∈ VPp , t ]

p = [DEp , v , t ∀ p , v ∈ VPp , t ]T

where x and u represent vectors of continuous and binary decision variables, xA and uA represent vectors of continuous and binary auxiliary decision variables for evaluating the system lead time, and p is a vector of demand parameters. The 9831

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

Table 2. Parameter Values Regarding Raw Materials/ Products in Case Study raw material/product a IP ID ωIR r /ωp /ωp : inventory cost ($/unit) αr: inventory deterioration fraction ηRr,i/ηPp,i: mass balanceb production task A production task B production task C production task D production task E ωRr,h/ωPp,h: purchase cost ($/unit) SA SB λRr,h/λPp,h: purchase capacity (unit) SA SB δOh,b: delivery delay (day) S A to IFBR A S B to IFBR A

Figure 3. Diagram of supply chain for dynamic operability analysis case study 1.

RM

IP

PA

PB

4 0.01

7.5 n/a

15 n/a

17 n/a

−1.43 −1.32 −1.1 n/a n/a

n/a n/a 1 −1.3 −1.2

1 n/a n/a 1 n/a

n/a 1 n/a n/a 1

n/a 30

152 n/a

n/a n/a

n/a n/a

n/a 3000

300 n/a

n/a n/a

n/a n/a

n/a 1

4 n/a

n/a n/a

n/a n/a

a

Penalty for not maintaining inventory at set point (φ) is 29.75 ($/unit). bNegative represents consumption; positive represents production; 1 represents main product of the production task.

Table 3. Parameter Values Regarding Regional Distribution Centers in Case Study 1

Figure 4. Manufacturing routes in case study 1.

regional distribution center

Table 1. Parameter Values Regarding Production Tasks in Case Study 1

A shipment delay (day) IFBR A Tb,v: shipment distance (km)a IFBR A μp,v: mean of step change, Up,v (unit) PA PB σp,v2: variance of step change, Up,v (unit2) PA PB upper bound of customer demand (unit) PA PB

production task γi,b: capacity (unit) IFBR A l ωi: production cost ($/unit) δM i : manufacturing delay (day)

A

B

C

D

E

750 35 2

900 50 3

1000 8 1

375 25 1

600 40 2

B

C

δDb,v:

5. CASE STUDIES In this section we present two case studies to demonstrate the applicability of the DOAF. In the first case study, the effect of product customization postponement on dynamic operability is investigated. Postponement allows inventory to be carried in a generic form until demand of finished/final product is known and contributes significantly to an agile supply chain.2 In the second case study the influence of manufacturing flexibility, as defined in this work, on dynamic operability is investigated, where two plant configurations are compared. These case studies are intended to rigorously explore the effect of postponement and manufacturing flexibility on responsiveness. ΔT is 1 day, and H is 21 days for each case study considered. We have assumed in these case studies that sufficient raw material can be procured, and sufficient production capacity exists at the plant to sustain the maximum demand level considered, over the optimization horizon; i.e., flexibility of the supply chain, which is a necessary condition for dynamic operability, is satisfied. The DOAF is modeled with GAMS 23.7.3 and solved using CPLEX 12.3 with the parallel mode enabled (8 threads) to a 3.5% optimality gap. Simulations were performed on a 2.93

a

2

1

3

250

150

500

100 n/a

n/a 75

80 65

100 n/a

n/a 100

225 225

130 n/a

n/a 105

125 110

Shipping cost (ωTp ) for all products: 0.01 ($/unit/km).

GHz IntelCore i7 DELL Vostro machine with 8 GB of RAM, running Windows 7 Professional 64-bit. 5.1. Case Study 1. The process supply chain system considered in case study 1 is illustrated in Figure 3. The system comprises two suppliers (S A and S B), one which supplies a biomass feedstock (RM), and the other which supplies an intermediate product (IP), to the biorefinery (IFBR A). The biorefinery has the capability to produce two products (PA and PB) that are shipped to regional markets, as well as produce an intermediate product (IP). There is a capability to procure IP or produce IP. Three regional distribution centers (RDC A, RDC B, and RDC C) serve the regional market demand of PA and/or PB. Dynamic operability of the process supply chain is explored for two manufacturing routes at the biorefinery, subcases I and II, respectively, as illustrated in Figure 4. Figure 4b represents a more agile route (subcase 1-II) as compared to the route illustrated by Figure 4a (subcase 1-I), due to the 9832

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

Table 4. Model Summaries and Optimization Results for Case Study case 1-I 1-II 1-III 1-IV a

dub a 800 800 800 800

no. of scenarios 50 50 50 50

continuous variables

discrete variable

52 512 56 813 52 512 56 813

6300 7350 6300 7350

obj fnb value

CPU (s) 6.865 1.000 2.861 8.934

× × × ×

3

10 103 104 102

9.667 8.882 9.871 9.188

× × × ×

5

10 105 105 105

95% confidenceb interval ±1.478 ±1.433 ±1.456 ±1.445

× × × ×

104 104 104 104

Equivalent upper bound for each inventory set point in vector d. bAverage over all values of ε considered.

Figure 5. Inventory set points for case study 1: subcase 1-I (left); subcase 1-II (right).

production of an intermediate product through an additional production task (C). This allows postponing the production of finished products PA and PB, and holding inventory in a more generic form (IP), until demand is known. Each production task is assigned to a distinct batch process unit at the plant. We have also considered an extension to subcase I and II, denoted as subcases 1-III and 1-IV, respectively, where the effect of an increase in the delivery delay of raw material RM from S A to IFBR A is explored. The parameter values used in the case study are summarized in Tables 1−3. As indicated in Table 1, the production cost is comparable, and the manufacturing delay is equivalent to

produce a unit of PA and PB through manufacturing route 1 or 2. As summarized in Table 2, the inventory storage cost of PA and PB is greater than IP or RM. Furthermore, the shipping delays between IFBR A and regional distribution centers are summarized in Table 3. Model summaries and optimization results, for each subcase, are given in Table 4. As observed from the table, 50 scenarios have been applied for capturing uncertainty for each optimization problem solved (i.e., for each value of ε considered). The objective function value and 95% confidence interval (average over all problems solved to generate the Pareto frontier) are given in the table. The confidence interval 9833

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

Figure 6. Inventory set points for case study 1: subcase 1-III (left); subcase 1-IV (right).

centers (Figure 5E,F) for a particular system lead time are in general larger with manufacturing route 1 (subcase 1-I) relative to manufacturing route 2 (subcase 1-II). An interesting result is that the inventory set point of RM with manufacturing route 1 (Figure 5A) is zero for a system lead time greater than zero. We suspect that the reason for this result is the relatively short delivery delay for procuring RM, relative to the manufacturing delay to produce PA and PB from RM, whereas the manufacturing delay to produce IP from RM (route 2) is equivalent to the delivery delay to procure RM (1 day). Therefore, for subcase I, it is likely necessary to store PA and PB to meet the system lead time condition defined, and RM is procured and converted to PA and PB to replenish inventory. We now consider the subcases where the delivery delay of RM from S A to IFBR A has increased from 1 to 2 days for both manufacturing routes (subcase 1-III and subcase 1-IV). Figure 6 illustrates the inventory set points computed by the DOAF for subcase 1-III and 1-IV. An increase in the delivery

of the objective function value is within our defined acceptable tolerance. The computational times reported represent the total CPU time (s) to generate the Pareto curve (i.e., the summation of computation times to solve the DOAF for all values of ε considered). Figure 5 summarizes the inventory set points determined by the DOAF for subcases 1-I and 1-II, for an increase in system lead time. In general, the set points have an overall decreasing trend, as we would expect. Because it is more cost efficient to store IP (subcase 1-II), we can observe in Figure 5(D) that the IP inventory set point remains relatively constant, and it only begins to decrease when the system lead time is greater than approximately three. Furthermore, we expect IP to be stored at a greater quantity than PA and PB at the biorefinery, because IP is a more generic product and can be processed further into PA or PB, once demand is realized. Note that the manufacturing delay to produce PA or PB from IP is less as compared to producing PA and PB from RM. From this we can observe that the inventory set points of products at the regional distribution 9834

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

Table 5. Parameter Values Regarding Production Tasks in Case Study 2 production task γi,b: capacity (unit) IFBR A ωIi : production cost ($/unit) δM i : manufacturing delay (day)

A

B

C

500 102 3

525 64 2

300 40 1

Table 6. Parameter Values Regarding Raw Materials/ Products in Case Study 2 raw material/product a IP ID ωIR r /ωp /ωp : inventory cost ($/unit) αr: inventory deterioration fraction ηRr,i/ηPp,i: mass balanceb production task A production task B production task C ωRr,h/ωPp,h: purchase cost ($/unit) SA λRr,h/λPp,h: purchase capacity (unit) SA δOh,b: delivery delay (day) S A to IFBR A

Figure 7. Pareto relationship between economic performance and system lead time for case study 1.

delay of RM from S A to IFBR A by 1 day has a pronounced effect on the inventory set point of RM (Figure 6A,B). It is likely that a larger RM inventory set point is necessary to fulfill the system lead time condition defined, because the delivery delay has increased. This is particularly apparent for manufacturing route 1, where previously (subcase 1-I) the RM inventory set point was zero for a system lead time greater than zero. Figure 7 illustrates the Pareto frontier between economic performance and system lead time for subcases 1-I, 1-II, 1-III, and 1-IV. As we would expect, the overall trend is similar; for a decrease in system lead time (increase in dynamic operability), economic performance worsens. Furthermore, the worst case system lead time (i.e., inventory set points are zero) is the same for subcases 1-I and 1-II, because the total manufacturing delay of routes 1 and 2 is equivalent. However, the Pareto curve for subcase 1-I (manufacturing route 1) is above the Pareto curve for subcase 1-II (manufacturing route 2), which signifies that at an equivalent system lead time economic performance is greater for subcase 1-II. Of importance is that at a high system lead time, the gap between economic performance reduces,

RM

LA

SA

MA

2.05 0.01

20 n/a

35 n/a

50 n/a

−6.6 −5 −5

1 n/a n/a

n/a 1 n/a

n/a n/a 1

25

n/a

n/a

n/a

5000

n/a

n/a

n/a

3

n/a

n/a

n/a

a

Penalty for not maintaining inventory at set point (φ): 100 ($/unit). b Negative represents consumption; positive represents production; 1 represents main product of the production task.

possibly indicating that the additional capital cost for implementation of manufacturing route 2 at the biorefinery is not cost justified if a lower level of dynamic operability is sufficient in the supply chain environment. This analysis is significant because it illustrates that a more agile manufacturing route is not always the best option, and it depends on the lead time that customers are expecting for products. Figure 7 illustrates the effect of an increase in the delivery delay of RM by 1 day on the Pareto curve with manufacturing

Figure 8. Diagram of supply chain for dynamic operability analysis case study 2. 9835

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

of each product reflects a biochemical manufacturing route from RM, where major unit operations include pressurized low polarity water (PLPW) extraction, enzymatic hydrolysis, fermentation (i.e., glucose to acids), and product recovery. The unit operations to produce each product are equivalent, except for the microorganism used in the fermenter and for product recovery. The PLPW extraction operation separates the major components of the woody biomass (i.e., cellulose, hemicellulose, and lignin). The stream containing lignin and cellulose is sent to a enzymatic hydrolysis unit operation step, where cellulose is broken down into glucose for fermentation. In the fermentation unit operation step, different microorganisms are used to produce LA, SA, and MA. Lignin is normally used as a fuel in the boilers. The hemicellulose stream, which contains high levels of xylose and other C5 sugars is sent to a fermentation process to produce xylitol; however, this process will not be considered in the current analysis. The fermentation operation is assumed to contribute to the only significant manufacturing time with respect to the duration of each time period, ΔT. Furthermore, a process unit, as defined in this study, represents the collection of unit operations for producing LA, SA, and MA. Production tasks A, B, and C to produce LA, SA, and MA, respectively, are assigned to process units in the biorefinery. This case study explores the influence of manufacturing flexibility through two distinct process plant configurations, where in a more flexible configuration (subcase 2-I) each production task is assigned to a distinct batch process unit, and in a less flexible configuration (subcase 2-II), production task A is assigned to a distinct batch process unit; however, production tasks B and C are assigned to the same batch process unit, i.e., SA and MA cannot be produced at the same time at IFBR A. Similar to case study 1, parameter values are summarized in Tables 5−7, and optimization results in Table 8. Of note is that the manufacturing delay for production tasks A, B, and C are 3, 2, and 1 day(s) respectively, as provided in Table 5, and are kept the same for each subcase considered. As illustrated in Figure 9, the inventory set points of products at the biorefinery and regional distribution centers depict the same decreasing trend as in case study 1. However, inventory set points at the biorefinery and regional distribution centers are slightly higher for subcase 2-II (less flexible configuration) than for subcase 2-I, as we would expect. From the figure, it would seem that the added manufacturing flexibility in subcase 2-I has a relatively small effect on the inventory set points. An interesting result from the figure is the increase in the inventory set points of raw material (Figure 9A,B) for a system lead time between approximately one and three; likely a result from the notion that it is more efficient to store RM rather than product for an increase in system lead time. The set points will then decrease when a portion of raw material can be procured while still satisfying the system lead time condition defined. Through inspection of the Pareto frontier for each subcase (Figure 10), the economic performance is improved for subcase 2-I, as compared to subcase 2-II for an equivalent system lead time. Furthermore, we can observe from the figure that the

Table 7. Parameter Values Regarding Regional Distribution Centers in Case Study 2 regional distribution center A δDb,v: shipping delay (day) IFBR A Tb,v: shipment distance (km)a IFBR A μp,v: mean of step change, Up,v (unit) LA SA MA σp,v2: variance of step change, Up,v (unit2) LA SA MA upper bound of customer demand (unit) LA SA MA a

B

C

2

1

3

250

150

500

40 25 n/a

35 30 20

n/a 20 20

56.25 25 n/a

49 36 16

n/a 25 25

62.5 40 n/a

56 48 32

n/a 35 35

Shipping cost (ωTp ) for all products: 0.05 ($/unit/km).

route 1 and 2 (subcases 1-III and 1-IV). As expected, both curves have shifted to the right, indicating that the economic performance has reduced for an equivalent system lead time. For subcases 1-III and 1-IV the worst case system lead time has increased by 1 day. This extension illustrates how the proposed framework can be applied for investigating the sensitivity of the system lead time (dynamic operability) to parameters such as delivery delay, among others. We have shown in case study 1 that for a small change in the delivery delay the effect on dynamic operability is appreciable, which may lead to supply chain/manufacturing design modifications to reduce the sensitivity. 5.2. Case Study 2. Case study 2 has been adapted from a case study originally presented in Mansoornejad,35 which describes a potential implementation of the intergrated forest biorefinery (IFBR) paradigm at an existing pulp and paper mill. The product portfolio option proposed in ref 35 represents the production of lactic acid (LA), succinic acid (SA), and malic acid (MA), reflecting a set of low-volume, high-value bioproducts. Benefits for producing high-value, low-volume specialty bioproducts, as compared to commodity type bioproducts (e.g., cellulosic ethanol) include low competition in the market, small quantities of biomass feedstock required, high margins, low price volatility, and common process elements between similar products.35 Of importance is that this option represents an adjacent or decoupled biorefinery; i.e., new process units for producing bioproducts are not tightly integrated with existing units in the mill that produce conventional pulp and paper products. The process supply chain system considered for case study 2 is illustrated in Figure 8, where a woody biomass feedstock (RM) is procured from a single supplier (S A). The production

Table 8. Model Summaries and Optimization Results for Case Study 2 case 2-I 2-II a

dub a 750 750

no. of scenarios 40 40

continuous variables

discrete variables

CPU (s)

8400 8400

1.016 × 10 9.491 × 104

67 896 67 896

Equivalent upper bound for each inventory set point in vector d.

b

obj fn 3

b

value

1.524 × 10 1.592 × 106 6

95% confidence

b

interval

±2.758 × 104 ±2.773 × 104

Average over all values of ε considered. 9836

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

Figure 9. Inventory set points for case study 2: subcase 2-I (left); subcase 2-II (right).

worst case system lead time has increased from subcase 2-I to subcase 2-II, as a result of the reduced manufacturing flexibility in the biorefinery. This case study illustrates the importance of this analysis framework during supply chain design, since it captures quantitatively the relationship between manufacturing flexibility and dynamic operability. This allows a decision-maker to rigorously analyze the added benefits obtained from improving manufacturing flexibility at the plant, to justify the increased capital cost of a biorefinery with more flexibility.

6. CONCLUSIONS AND FUTURE WORK A quantitative approach for dynamic operability analysis of process supply chains is developed. Responsiveness, a key attribute of a dynamically operable (agile) supply chain, is evaluated through a novel dynamic supply chain optimization formulation, where the response dynamics are captured through the model. Inventory set points, which are equated to the initial inventory position, are integrated with system lead time through a two-stage stochastic programming approach.

Figure 10. Pareto relationship between economic performance and system lead time for case study 2.

9837

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

KBb KCb BRr BPp VPp

Inventory set points are determined prior to the resolution of demand uncertainty in the first stage, and operational decisions are determined in the second stage in response to how uncertainty resolves. The dynamic operability analysis framework is formulated as a bi-criterion optimization problem to estimate the Pareto frontier between economic performance and system lead time. Two case studies are presented to demonstrate the applicability of the proposed framework. In case study 1 the influence of product customization postponement on dynamical operability is investigated, and in case study 2 the influence of manufacturing flexibility on dynamical operability is investigated. Although the effects of these concepts on responsiveness have been discussed in the literature, the relationship between them has not been rigorously explored, as in this work. An important extension for this work is to integrate the analysis formulations with a comprehensive framework for designing dynamically operable supply chains. Future work should investigate whether design decisions can be integrated within the existing framework, or whether an iterative design and analysis approach is more practical. Other avenues for investigation include decomposition approaches and alternative scenario reduction techniques, to more efficiently solve largescale stochastic mixed-integer programming problems.



PVv

Binary Variables

Xi,b,t one if production task i at biorefinery b begins at time period t; and 0 otherwise Ymet p,v,t one if demand (and accumulated back orders) for product p at regional distribution center v can be served and sustained at the new demand level at time period t; and 0 otherwise Continuous Variables

ORr,h,b,t purchase order of raw material r to supplier h from biorefinery b at time period t OPp,h,b,t purchase order of product p to supplier h from biorefinery b at time period t Pi,b,t quantity of main product produced by production task i at biorefinery b beginning at time period t IRr,b,t inventory of raw material r at biorefinery b at time period t IPp,b,t inventory of product p at biorefinery b at time period t IDP inventory of product p at regional distribution center v p,v,t at time period t θRr,b,t quantity of raw material r inventory which drops below the inventory set point at biorefinery b at time period t θPp,b,t quantity of product p inventory which drops below the inventory set point at biorefinery b at time period t θDP quantity of product p inventory which drops below the p,v,t inventory set point at regional distribution center v at time period t SRr,b inventory set point of raw material r at biorefinery b SPp,b inventory set point of product p at biorefinery b SDP inventory set point of product p at regional distribution p,v center v FPp,b,v,t quantity of product p shipped from biorefinery b to regional distribution center v at time period t FSp,v,t quantity of product p distributed from regional distribution center v to satisfy demand and accumulated back orders at time period t Bp,v,t quantity of back orders accumulated for product p at regional distribution center v at time period t mp,v,t indicates the time period when Ymet p,v,t switches from 0 to 1 Z economic criterion LTp,v lead time of product p at distribution center v SLT system lead time

AUTHOR INFORMATION

Corresponding Author

*C. L. E. Swartz: e-mail, [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS Financial support for this research through the Natural Sciences and Engineering Research Council of Canada (NSERC) Strategic Network on Value Chain Optimization (VCO), Ontario Graduate Scholarship Program (OGS), and the McMaster Advanced Control Consortium (MACC) is gratefully acknowledged.



NOMENCLATURE

Sets/Indicies

r p h b i k v s t HRr HPp ICR r,b ICp,b IM p,b IRk,b Ib

set of batch process units at biorefinery b set of continuous process units at biorefinery b set of biorefineries that consume raw material r set of biorefineries that consume or produce product p set of regional distribution centers which store and distribute product p set of products which are stored at regional distribution center v

raw material product suppliers integrated forest biorefinery production task process unit regional distribution center scenario time period set of suppliers of raw material r set of suppliers of product p set of production tasks which consume raw material r at biorefinery b set of production tasks which consume product p at biorefinery b set of production tasks which produce product p at biorefinery b set of production tasks which can be assigned to process unit k at biorefinery b set of production tasks performed by units at biorefinery b

Parameters

ηRr,i ηPp,i λRr,h λPp,h δOh,b δM i 9838

mass balance coefficient of raw material r of production task i mass balance coefficient of product p of production task i maximum quantity of raw material r that can be purchased from supplier h during a time period maximum quantity of product p that can be purchased from supplier h during a time period delivery delay between when an order is made to supplier h and when the order arrives at biorefinery b and is ready to be processed manufacturing delay of production task i dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research δDb,v αr γi,b Up,v DEp,v,t μp,v σp,v2 ΔDp,v ωRr,h ωPp,h ωIi ωTp ωIR r ωIP p ωID p φ Tb,v ρs ΔT H



Article

Computer-Aided Process Design; Westerberg, A. W., Chien, H. H., Eds.; AIChE: New York, 1984; pp 931−1030. (11) Morari, M.; Perkins, J. D. Design for Operations. In Foundations of Computer-Aided Process Design; Biegler, L. T., Doherty, M. F., Eds.; AIChE Symposium Series, No. 304; 1995; Vol. 91, pp 105−114. (12) Lenhoff, A. M.; Morari, M. Design of resilient processing plants−I. Process design under consideration of dynamic aspects. Chem. Eng. Sci. 1982, 37, 245−258. (13) Mohideen, M. J.; Perkins, J. D.; Pistikopoulos, E. N. Optimal design of dynamic systems under uncertainty. AIChE J. 1996, 42, 2251−2272. (14) Schweiger, C. A.; Floudas, C. A. Interactions of design and control: Optimization with dynamic models. In Optimal Control: Theory, Algorithms and Applications; Hager, W. W., Pardalos, P. M., Eds.; Kluwer Academic Publishers: Amsterdam, 1998; pp 388−435. (15) Bansal, V.; Perkins, J. D.; Pistikopoulos, E. N.; Ross, R.; van Schijndel, J. M. G. Simultaneous design and control optimization under uncertainty. Comput. Chem. Eng. 2000, 24, 261−266. (16) Baker, R.; Swartz, C. L. E. Simultaneous Solution Strategies for Inclusion of Input Saturation in the Optimal Design of Dynamically Operable Plants. Optimization and Engineering 2004, 5, 5−24. (17) Sakizlis, V.; Perkins, J. D.; Pistikopoulos, E. N. Recent advances in optimization-based simultaneous process and control design. Comput. Chem. Eng. 2004, 28, 2069−2086. (18) Naylor, J. B.; Naim, M. M.; Berry, D. Leagility: Integrating the lean and agile manufacturing paradigms in the total supply chain. Int. J. Production Economics 1999, 62, 107−118. (19) Backx, T.; Bosgra, O.; Marquardt, W. Towards intentional dynamics in supply chain conscious process operations. AIChE Symp. Ser. 1998, 320, 5−20. (20) Holweg, M. The three dimensions of responsiveness. International Journal of Operations & Production Management 2005, 25, 603− 622. (21) Beamon, B. M. Measuring supply chain performance. International Journal of Operations & Production Management 1999, 19, 275− 292. (22) You, F.; Grossmann, I. E. Optimal Design and Operational Planning of Responsive Process Supply Chains. In Process Systems Engineering: Supply Chain Optimization; Papageorgiou, L. G., Georgiadis, M. C., Eds.; Wiley-VCH Verlag GmbH & Co. KGaA: Weinheim, 2007; Vol. 3, Chapter 4, pp 107−134. (23) You, F.; Grossmann, I. E. Design of responsive supply chains under demand uncertainty. Comput. Chem. Eng. 2008, 32, 3090−3111. (24) You, F.; Grossmann, I. E. Balancing responsiveness and economics in process supply chain design with multi-echelon stochastic inventory. AIChE J. 2011, 57, 178−192. (25) Ekşioǧlu, S. D.; Acharya, A.; Leightley, L. E.; Arora, S. Analyzing the design and management of biomass-to-biorefinery supply chain. Computers & Industrial Engineering 2009, 57, 1342−1352. (26) Shah, N.; Pantelides, C. C.; Sargent, R. W. H. A general algorithm for short-term scheduling of batch operations−II. Computational issues. Comput. Chem. Eng. 1993, 17, 229−244. (27) Kondili, E.; Pantelides, C. C.; Sargent, R. W. A general algorithm for short-term scheduling of batch operations−I. MILP formulation. Comput. Chem. Eng. 1993, 17, 211−227. (28) Ierapetritou, M. G.; Floudas, C. A. Effective Continuous-Time Formulation for Short-Term Scheduling. 1. Multipurpose Batch Processes. Ind. Eng. Chem. Res. 1998, 37, 4341−4359. (29) Jung, J. Y.; Blau, G.; Pekny, J. F.; Reklaitis, G. V.; Eversdyk, D. A simulation based optimization approach to supply chain management under demand uncertainty. Comput. Chem. Eng. 2004, 28, 2087−2106. (30) Birge, J. R.; Louveaux, F. Introduction to stochastic programming; Springer: New York, 1997. (31) Hammersley, J. M.; Handscomb, D. C. Monte Carlo Methods; Methuen & Co Ltd.: London, 1964. (32) Liu, M. L.; Sahinidis, N. V. Optimization in Process Planning under Uncertainty. Ind. Eng. Chem. Res. 1996, 35, 4154−4165.

shipping delay between biorefinery b and regional distribution center v deterioration percentage of raw material r during the duration of a time period maximum production capacity of main product of production task i at biorefinery b magnitude of demand step change for product p at regional distribution center v demand for product p at regional distribution center v at time period t mean of demand (step change) for product p at regional distribution center v variance of demand (step change) for product p at regional distribution center v range from mean demand of product p at regional distribution center v cost per unit of raw material r purchased from supplier h cost per unit of product p purchased from supplier h cost per unit of main product p produced by production task i cost per unit, per km of product p shipped from biorefinery b to regional distribution center v cost per unit of inventory of raw material r over the duration of a time period cost per unit of inventory of product p at the biorefinery over the duration of a time period cost per unit of inventory of product p at the regional distribution center over the duration of a time period penalty cost per unit of product or raw material inventory which drops below the inventory set point distance (km) between biorefinery b and regional distribution center v probability of the occurrence of scenario s duration of each time period optimization horizon length

REFERENCES

(1) Papageorgiou, L. G. Supply Chain Optimisation for the process industries: Advances and opportunities. Comput. Chem. Eng. 2009, 33, 1931−1938. (2) Christopher, M. The Agile Supply Chain: Competing in Volatile Markets. Industrial Marketing Management 2000, 29, 37−44. (3) Christopher, M.; Towill, D. An integrated model for the design of agile supply chains. International Journal of Physical Distribution & Logistics Management 2001, 31, 235−246. (4) Wising, U.; Stuart, P. Identifying the Canadian forest biorefinery. Pulp Pap. Can. 2006, 107, 25−30. (5) Chambost, V.; Stuart, P. Selecting the most appropriate products for the forest biorefinery. Ind. Biotechnol. 2007, 3, 112−119. (6) Fisher, M. L. What is the right supply chain for your product? Harvard Business Review 1997, 75, 105−116. (7) Mansoornejad, B.; Chambost, V.; Stuart, P. Integrating product portfolio design and supply chain design for the forest biorefinery. Comput. Chem. Eng. 2010, 34, 1497−1506. (8) Mansoornejad, B.; Pistikopoulos, E. N.; Stuart, P. Metrics for evaluating the forest biorefinery supply chain performance. Comput. Chem. Eng. 2013, 54, 125−139. (9) Mansoornejad, B.; Pistikopoulos, E. N.; Stuart, P. R. Scenariobased strategic supply chain design and analysis for the forest biorefinery using an operational supply chain model. Int. J. Production Economics 2013, 144, 618−634. (10) Grossmann, I. E.; Morari, M. Operability, resiliency and flexibility: Process design objectives for a changing world. In Proceedings of the International Conference on Foundations of 9839

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840

Industrial & Engineering Chemistry Research

Article

(33) Mak, W.-K.; Morton, D. P.; Wood, R. K. Monte Carlo bounding techniques for determining solution quality in stochastic programs. Operations Research Letters 1999, 24, 47−56. (34) You, F.; Wassick, J. M.; Grossmann, I. E. Risk management for a global supply chain planning under uncertainty: Models and algorithms. AIChE J. 2009, 55, 931−946. (35) Mansoornejad, B. Design for flexibility in the forest biorefinery supply chain. Ph.D. Thesis, École Polytechnique de Montréal, 2012.

9840

dx.doi.org/10.1021/ie500608w | Ind. Eng. Chem. Res. 2014, 53, 9825−9840