Ind. Eng. Chem. Res. 2010, 49, 11491–11505
11491
Optimal Scheduling of Pipeline Systems with a Resource-Task Network Continuous-Time Formulation Pedro M. Castro Unidade Modelac¸a˜o e Optimizac¸a˜o de Sistemas Energe´ticos, Laborato´rio Nacional de Energia e Geologia, 1649-038 Lisboa, Portugal
This work addresses the short-term scheduling problem involved in the pipeline transportation of refined petroleum products to distribution centers serving local markets. A new continuous-time formulation is proposed that can handle complex treelike systems involving multiple refineries and depots. The novelty concerns the modeling of the complex storage policy associated with the flow of multiple products inside the pipeline. This is accomplished through the use of real and accumulated volume resources for every product, coupled with filling, moving, and emptying tasks, under the scope of the resource-task network process representation. Computational studies are performed to compare the performance of the new approach to three models from the literature that rely on the concept of batches rather than events to model what is in reality a continuous system. The results show that the proposed approach is wider in scope with the typical drawback of being less efficient for structural subclasses of the general problem. 1. Introduction Pipelines represent the most reliable and cost-effective way of transporting large quantities of refined petroleum products from major supply sources such as refineries to local distribution centers.1,2 Pipelines operate continuously unlike other transportation modes (e.g., road, railroad, and vessel) and can carry a wide range of products including various grades of gasoline, diesel, liquefied petroleum gas (LPG), and jet fuel. Different products are pumped back to back in the pipeline usually without any device separating them, resulting in some product contamination due to the occurrence of mixing at the interface. The extent of mixing is a function of the flow rate and of the products involved, so the pipeline should operate in turbulent flow with limited stops and with a subset of product sequences being forbidden for quality reasons.3 The primary aim of the pipeline scheduler is to meet the demands established by the markets while keeping the inventory levels at the distribution depots within given limits. It is critical that correct decisions are made to avoid major delays resulting from the large distances involved. Ideally, one should minimize the operating cost by considering the daily fluctuation in electricity price in the pumping costs, as well as interface and inventory costs.1 Planning and scheduling of multiproduct pipelines has received some attention in the past decade. From a structural point of view, pipeline systems range from a straight pipeline with no branches, to a treelike structure where products can take alternative routes before arriving at a particular destination,4 or even to networklike systems where pipelines connect nodes, with these being refineries or depots. Earlier optimization approaches have focused on straight pipeline systems consisting of a single refinery source with one3,5 or multiple depots located along the line.1,2 The scope is now being widened through the consideration of multiple refinery sources feeding a straight pipeline6 or a single source serving a treelike distribution sytem.4 Yet, to the best of our knowledge, no rigorous approach has been proposed for the most general case. Rejowski and Pinto1 have proposed a discrete-time, mixed integer linear programming (MILP) model where the pipeline is divided into as many segments as the number of depots. Each segment is then further divided into a few constant-volume packs
that can hold one product at a time. While one can define diverse volumes for packs belonging to different segments, to better adapt to varying total capacities, the model can be classified as being discrete volume. The discrete-time model was later shown7 to be a particular case of a continuous-time mixed integer nonlinear programming model that was able to achieve better results while covering hydraulic aspects of pipeline operation. Cafaro and Cerda´2 were the first to propose an MILP model using both a continuous-time and continuous-volume representation. The pipeline contents are viewed as a set of slugs (batches), each holding a particular product, which are characterized by length and completion timing variables as well as volume variables. The latter indicate the location coordinate and size, with respect to the end of a new pumping run. Overall, the optimization model can determine the set of slug injections and their volumes together with the sequencing and scheduling of such slugs. The authors later8 proposed a rolling-horizon9 algorithm to efficiently handle large time horizons and multiple due dates, which acts as a rescheduling procedure whenever new demand data from the customers is received. Relvas et al.5 extend the model by Cafaro and Cerda´2 to handle daily due dates rather than just at the end of the time horizon, while integrating operational aspects at the downstream distribution center. In particular, each new lot that arrives at the tank farm must respect a settling period for quality improvement and go through a set of batch control and approval tests.3 Within the context of a real-world case study from a Portuguese oil distribution company, the authors recognized that the new features resulted in a lower performance when compared to the base model. Tests were then performed to measure the trade-off between the decrease in solution quality and increase in performance when going from free to mixed to fixed product sequences, with the latter corresponding to the strategy used by the plant schedulers. New features were later10 included to account for pipeline stoppages and settling periods that vary between products, while a reactive scheduling strategy to deal with several situations that may occur within the system was proposed. This is concerned with generating a new solution with minimum changes from the previous one. Pipeline systems are not the only ones receiving attention from the scheduling community. The past two decades have
10.1021/ie1010993 2010 American Chemical Society Published on Web 09/17/2010
11492
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 1. Treelike multiproduct pipeline distribution system.
seen major developments with emphasis on the so-called general approaches that can virtually handle any problem despite probably not being the most computationally efficient alternative for specific problem classes. The reader is directed to the excellent paper by Me´ndez et al.11 for an overview of the most relevant models for batch plants, some of which are also suitable or can be extended to continuous plants. The key element is that general approaches rely on unified frameworks for process representation: either the state-task network12 (STN) or the resource-task network13 (RTN). RTN-based formulations have been shown to be very versatile, successfully handling a variety of real-life problems.14-16 As an example, the successful case reported by Wassick16 is just one of 12 applications addressed at Dow Chemical Co., demonstrating that industry is slowly acknowledging their potential. One important aspect is that the process representation can be made mostly independent of the scheduling model’s underlying time representation.14,15 In other words, one process model will be suitable for both discrete-time and continuoustime approaches, which is particularly relevant given the fact that there is still a lot of uncertainty concerning the most suitable time representation for a particular problem. This paper presents a new RTN-based continuous-time and continuous-volume MILP formulation for the short-term scheduling of multiproduct pipeline systems that can handle multiple discrete due dates as well as continuous demand rates. The novelties are the following: (i) the description of the system as a collection of building blocks consisting of refinery, depots, and pipeline segments, which allows consideration of virtually any structure with multiple refineries, depots, and pipeline branches; (ii) the pipeline segment model, involving a complex shared storage policy featuring a first-in, first-out strategy, which requires the use of a few pipeline state and volume resources together with a variety of tasks. The remainder of the paper is structured as follows. Following the problem definition in section 2, we give in section 3 a detailed description of the system and pipeline segment models as resource-task networks. The continuous-time mathematical model is then given in section 4, with the indication of the changes that need to be made to transform it into a discretetime but still continuous-volume model. The approach is illustrated in section 5 with three example problems taken from the literature, leaving the detailed performance studies for section 6. The latter includes a comparison to previous continuous-time and continuous-volume models. Finally, we give the conclusions in section 7.
2. Problem Definition A set of refineries RF must distribute a set P of petroleum products (e.g., gasoline of different grades, diesel, LPG, jet fuel) to a few depots DP through a pipeline system. The pipeline system can comprise a single pipeline or a more complex treelike structure with multiple branches. Assuming that a pipeline segment s ∈ S ends with a refinery, depot, or junction, we arrive at the general system representation given in Figure 1. A few product dedicated storage tanks are typically available at the refinery and depots, with at most one being connected to the pipeline at each time. Rather than considering all available tanks individually (tank farm management), we consider a single aggregate tank per product at each depot and refinery. Given are the following: (i) the maximum pumping flow rate at the refineries Frfmax (m3/h); (ii) the initial product inventory at 3 the refineries Vir,0 p,rf (m ); (iii) the capacity of the different pipeline 3 segments Vs (m ), which are always full; (iv) the minimum, Vmin p,dp, 0 maximum, Vmax p,dp, and initial, Vp,dp, storage capacities for product p at depot dp (m3); (v) the initial product inventory for each r,0 , together with exact product location, pipeline segment, Vp,s 3 which can be specified through an accumulated volume Va,0 p,s (m ); (vi) the product demands at the depots. The latter will be modeled either as a constant continuous removal of the total demand, dtp,dp (m3), throughout the duration of the time horizon, or as an instantaneous removal of material that can occur at any hour of the day, dp,dp,dy,hr (m3). Additional data include forbidden product sequences to avoid excessive product contamination and information about the structure of the pipeline system. The objective is to optimize the short-term operation of the multiproduct system for a given time horizon H. We will be focusing on the minimization of the total amount pumped from the refineries since the most challenging part of the problem is the correct management of the distribution depots. Note that there are upper and lower bounds on the tank capacities and it may take a considerable amount of pumping time to start filling them again. Lack of products in the tanks will affect local consumers, whereas excess may paralyze the transfer in the pipeline and even interrupt production in the refinery.17 2.1. Remarks. Minimizing the total pumped volume is directly related to pumping and inventory costs, main contributors to the typical1,2,4 objective of minimizing total operating cost. It is straightforward to consider such objective since (i) transition costs can be easily incorporated given the fact that the pipeline segment model features changeover tasks to explicitly identify the product sequence, (ii) there are variables
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
11493
Figure 2. RTN representation of a treelike multiproduct pipeline system.
that give the inventory at key points in time, which can be used to account for average inventory, and (iii) variable electricity pricing can be incorporated by following the guidelines given by Castro et al.14 for continuous processes. 3. Resource-Task Network Process Representation The system description in Figure 1 involves real entities such as pumps, tanks, pipeline segments, and products. Existing scheduling models for multiproduct pipeline systems1,2,4 feature parameters and variables that are related to these entities and thus can be natural alternatives for someone working in the field. The other possibility is to rely on a more general approach for process representation, the resource-task network,13 which has the advantage of being associated with powerful mathematical models that can handle both batch and continuous processes18 featuring a discrete13,16 or continuous representation of time, involving single time grid18 or unit-specific19 approaches. Provided that the underlying RTN-based formulation can handle all problem features efficiently, the only challenge that remains is identifying appropriate sets of resources and tasks. In other words, we need to systematically convert the real process entities into the RTN virtual entities. 3.1. System Overview. We consider the first equipment within the system boundary to be the upstream refinery pump (Pump_RF1). Pumps are the most critical resources since they are responsible for all the movement until the products reach the depots. In particular, the upstream pump will be responsible for taking the products out of the refinery tanks into the beginning of the pipeline. If contamination is neglected, products remain unchanged during the pipeline, so the problem here is keeping track of product location throughout the duration of the time horizon. For that, we need location resources:14,20 (i) inside refinery (e.g., P1_iRF1); (ii) leaving refinery (e.g., P1_lRF1); (iii) leaving pipeline segment (e.g., P1_lS1); (iv) inside depots (e.g., P1_DP1). Note that resources are represented as circles/ovals in Figure 2. The rectangles represent the tasks that consume and produce resources. Pumping tasks use the corresponding pump resource to produce leaving refinery or leaving segment resources, depending on whether we are dealing with the upstream or an intermediate refinery. They are continuous tasks that can be executed at a rate up to the maximum pumping flow rate. We define one pumping task per refinery and per product (e.g., Pump_P1_RF1). Notice in Figure 2 that refinery RF2 may
Figure 3. Alternative operating modes and paths for a product inside a pipeline.
produce P1_lS2 just like the empty tasks inside segment S2. Since for operational reasons6 this cannot happen simultaneously, we need to define valve resources (e.g., Valve_S2) and link them to both pumping and empty tasks. For a maximum availability of 1, one can model multiple entries to a pipeline segment with a single active entry at a time since each task will consume one unit of the valve resource. Furthermore, there is no need to explicitly define opening/closing valve tasks since the valve status will be known implicitly from the tasks being executed. At the end of a certain pipeline segment there may be a depot, so we need to define send to depot tasks (e.g., Send Depot_P1_DP1). Each processing rate will correspond to that of the active upstream refinery pump. Note that it is possible to process a particular product to a depot and to the following pipeline segment simultaneously. The leaving segment resources are thus also consumed by the pipeline segment filling tasks, which will now be described. 3.2. Pipeline Segment. The large majority of the modeling effort was concerned with the pipeline segment building block. It can be viewed as a shared storage tank being able to store multiple products simultaneously rather than just one a time.14 Furthermore, it follows a first-in first-out (FIFO) queuing procedure where a product starts leaving the segment right after going through the segment volume. This can be viewed as a complex storage policy, which to the best of our knowledge has not been previously handled by general scheduling formulations. The pumping schedule can be viewed as a sequence of product batches entering the pipeline. Depending on the relation between the batch and pipeline volumes, there are three operating paths; see Figure 3. For large batch volumes, the product enters the pipeline while pushing some other product out. The filling phase continues until all the pipeline is filled with the product. At that time, we enter the fill and empty phase
11494
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 4. Tasks and resources within a pipeline segment (illustration for product P1).
where there is no net accumulation of product inside the pipeline. Once the full batch has been pumped, there will eventually be another product pushing it out of the pipeline (empty phase). If the batch size is lower than the pipeline volume, the product will enter a moving phase right after filling is completed and before discharge is started. Finally, if the two volumes match, the end of the filling phase coincides with the beginning of the empty phase. We can now explain the rather complex superstructure given in Figure 4 for pipeline segment s. For simplicity, the segment index has been removed from the tasks and resources despite being required. The easiest operating mode occurs with the execution of the fill and empty task, in the center of the diagram, since there is no net product accumulation. Thus, the task is continuously consuming the outlet resource from the previous segment, e.g., product P1 leaving segment s - 1 (P1_lSs-1), while producing its own output resource, e.g., P1_lSs. The filling task (on the upper left) consumes the exact same volume resource (linked to tasks by solid lines) but produces four other sets of resources of the same type: (i) the inside pipeline location resource P1_iP that holds the product volume inside the pipeline; (ii) the accumulated volume resource, which holds the volume that entered the pipeline after start pumping P1 right until it starts to be discharged, P1_aV; (iii) the buildup volume resource P1_bV, required to model the case when there is another product responsible for the movement; (iv) the global resource, which although not necessary leads to a reduction in the integrality gap by forcing the active empty task to be processed at the same rate as the active filling task. This is
Figure 5. Illustrative example to identify active sets of tasks and resources.
because the inside pipeline and accumulated volume resources will be allowed to exist in excess contrary to the other two sets (see eq 24). To better understand the complex interactions going on, consider the example in Figure 5, where there is an additional product (P4) in the system besides those caught in the snapshot. The active tasks are the following: Fill_P1, Move_P2, Empty_P3, and Do Nothing_P4. Task Fill_P1 is producing volume resources G, P1_iP, P1_aV, P2_bV, P3_bV, and P4_bV. Resources P1_iP and P1_aV are accumulating, while the others are being consumed at the same rate. Resource G by task Empty_P3, task Move_P2 is converting P2_bV into P2_aV, P3_bV is being consumed by Empty_P3, while Do Nothing_P4 consumes P4_bV. Task Move_P2 also hides temporarily (discrete consumption/ production at its start/end) a small amount of P2_iP to prevent it from being executed when there is no P2 in the pipeline and hence distinguish it from Do Nothing_P2. In terms of the actual value of the resources in Figure 5, we have P1_iP ) P1_aV ) 2000, P2_iP ) 3000, P2_aV ) 5000, P3_iP ) 1000, and P3_aV ) 0 m3. It is important to highlight that, immediately before the start of the discharge of P3, we had P3_aV ) 6000 m3, which by being the capacity of the
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 6. Different batches of the same product are not allowed simultaneously in the same pipeline segment.
pipeline, triggered the execution of an instantaneous switch empty task, which made P3_aV ) 0. Besides the continuous tasks, Figure 4 shows six instantaneous tasks that make the switch between different states of the pipeline. One important aspect with respect to other problems dealing with changeovers21 is that we can no longer assume that the equipment is at a single state at a certain point in time (e.g., clean, empty, ready to process P1, etc.). Nevertheless, we can assume that the pipeline is at a single state concerning a particular product of the portfolio. Of the five possible states (unitary resources), four are related to the operating modes listed in Figure 3: (i) F_P1 stands for fill P1; (ii) M_P1 stands for move P1; (iii) E_P1 means while emptying; (iv) FE_P1 means during simultaneous filling and emptying of P1; N_P1 is active when P1 is not inside the pipeline segment. Note that unitary resources interact discretely with processing tasks (dashed line) being consumed at their start and produced at their end. The instantaneous tasks account for sequence-dependent changeovers explicitly, enabling the straightforward incorporation of transition costs in the objective function. For instance, the pipeline condition in Figure 5 resulted from pumping products in the sequence P3-P2-P1, meaning that earlier in time we had the execution of switch fill tasks P3_P2 and P2_P1 and we will have in the future switch empty P3_P2 and P2_P1. Yet, there are two possible types (A and B) of switch fill and switch empty tasks that differ in the nature of virtual resources involved. The distinction is that tasks of type B involve one fill and empty resource (e.g., FE_P1). Thus, since the batch of P2 is lower than the pipeline volume, we know that Switch Fill A_P2_P1 was executed, transforming F_P2 into M_P2 as well as N_P1 into F_P1. However, there is insufficient information to determine the type of switch fill P3_P2. With a large P3 batch, the task was of type B with P3 switching from a fill and empty state (FE_P3) to an empty state (E_P3). The same reasoning can be applied to the switch empty tasks. 3.3. Remarks. The proposed RTN representation of the pipeline system implicitly assumes that there can be at most a single batch per product in a certain pipeline segment; see Figure 6. This is a very important restriction that may prevent feasibility in some problems. Nevertheless, the representation is general if looked from a batch or slug (continuous stream of a single homogeneous product within the pipeline2) perspective. One will then need to associate slugs with products by extending the model presented in the following section with additional sets of variables and constraints. This will be the subject of future work. The other possibility is to define more segments, each with a volume smaller than the lowest possible batch volume. It is an option that will lead to a significant increase in the number of tasks and resources, which may compromise tractability. 4. Continuous-Time Formulation The continuous-time formulation for the short-term scheduling of multiproduct pipeline systems is given below. It gives rise to problems of the mixed integer linear programming (MILP) type and is essentially a simplification of the recent model by
11495
Figure 7. Single time grid used by the continuous-time model.
Figure 8. The model makes no distinction in the resource balances among three extreme cases.
Castro et al.14 for continuous plants, which can be considered an upgraded version of the general model for batch and continuous plants by Castro et al.18 A novel aspect concerns the incorporation of continuous demand rates, which can be used in conjunction with multiple due dates that are modeled as instantaneous removals of material. 4.1. Time Representation. The continuous-time formulation uses a single time grid to keep track of events taking place. It consists of |T| - 1 time slots (see Figure 7), with the interval boundaries being named event points or global time points.11 The absolute time of event point t, variable Tt (h), is going to be determined by the optimization. Instantaneous tasks occur at the event points, and the resource balances assume that continuous tasks start at event point t and end at t + 1. While there is no ambiguity if the amount processed by the continuous task, ξi,t (m3), is equal to its maximum processing rate times the slot duration, there is still some freedom to schedule the task within the slot otherwise; see Figure 8. Wherever a more rigorous approach is needed, other sets of variables and constraints can easily be added to the formulation.14,22 Bear also in mind that processing tasks last a single slot, hence the use of a single time index, which is not limiting in any way since multiple instances of the task can be performed to meet large volumes. 4.2. Parameters. The RTN representations given in Figures 2 and 4 are brought into the model by the structural parameters. Of the five sets discussed in Castro et al.14 only three are required in this problem. Parameter µr,i gives the total resource consumption (negative sign) at the start of the task, while µ j r,i gives the amount produced at its end (positive sign). Resources and tasks involved are the ones linked by dashed lines, and such parameters are linked to the binary variables Ni,t. As examples we have µPump_RF1,Pump_P1_RF1 ) -1; µ¯ Pump_RF1,Pump_P1_RF1 ) 1; µE_P1,Empty_P1 ) -1; µ¯ E_P1,Empty_P1 ) 1; µP1_aV,Switch Empty A_P2_P1 ) -Vs (1) On the other hand, for a particular resource r there may be a continuous interaction with task i during its execution, which is linked to variables ξi,t and modeled using parameters λr,i. Continuous interactions are represented by solid lines in Figures 2 and 4, and such parameters can take three possible values: -1, 0, and 1. For instance
11496
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
λP1_lS1,Send Depot_P1_DP1 ) -1; λP1_DP1,Send Depot_P1_DP1 ) 1; λP1_lSs-1,Fill and Empty_P1 ) -1; λP1_lSs,Fill and Empty_P1 ) 1
drr ) dtp,dp /H ∀p ∈ P, dp ∈ DP, r ∈ RFP p,dp
Resource availability at the start of the time horizon, Rr0, is known for all resources. Initially, all refinery pumps as well as 0 0 ) 1; RValve_S ) connection valves are available (e.g., RPump_RF 1 2 1). The real importance of these parameters is however related to the definition of the initial pipeline segment states, one for each product. Based on the capacity of the pipeline segment (Vs), the initial product inventory, and accumulated volumes r,0 a,0 , Vp,s , recall Figure 5), the following conditions apply. (Vp,s If there is no product inside the pipeline segment, the initial state is clearly do nothing. r,0 0 If Vp,s ) 0, then RN_P ) 1 ∀s ∈ S, p ∈ P p,s
(3)
If the product inventory and accumulated volume take the same positive value, then we have a filling state.
Alternatively, there may be multiple demand points associated with an instantaneous removal of material. Assuming TD as the full set of demand points, with subset TDdy,hr including the one linked to day dy and hour hr, the discrete amount of resource r removed from the system at demand point td is given by Πout r,td; see eq 13. The actual time of occurrence of point td is in turn assumed to be given by parameter tfxtd. out Πr,td ) dp,dp,dy,hr ∀p ∈ P, dp ∈ DP, dy ∈ DY, hr ∈ HR, td ∈ TDdy,hr, r ∈ RFP p,dp (13)
Finally, the maximum rates of processing tasks (subset Ic) are set to the maximum pumping flow rate, eq 14, unless they are pumping tasks (Ipump), for which the processing rate should be lower than the refinery’s maximum pumping rate.
∑
Fmax ) i
r,0 a,0 0 If Vp,s ) Vp,s * 0, then RF_P ) 1 ∀s ∈ S, p ∈ P p,s
If the accumulated volume is greater than the real product volume, then another product p′ has entered the pipeline after filling product p and so the slug of p is moving.
Frfmax + max Frfmax ∀i ∈ Ic
(14)
rf∈RF i∉Ipump
rf∈RF i∈Irfpump
(4)
a,0 r,0 0 If Vp,s > Vp,s , then RM_P ) 1 ∀s ∈ S, p ∈ P p,s
(12)
(2)
4.3. Constraints. We consider as objective function the minimization of the total amount of products pumped into the pipeline (m3).
(5)
min
∑
∑
(15)
ξi,t
i∈Ipump t∈T∧t*|T|
The empty state becomes active following the full depletion of the accumulated volume to zero, provided that the real volume is positive and lower than the pipeline segment’s volume.
r,0 IP Rr0 ) Vp,s ∀s ∈ S, p ∈ P, r ∈ Rp,s
(8)
Resources include equipment, the different product locations, and the possible pipeline segment states. The multiperiod excess resource balances keep track of resource availability over time and include two sets of constraints. Equation 16 deals with the discrete interactions with tasks. Volume resources (RCT) are continuously produced/consumed by processing tasks (Ic), so there can be an inventory change between t and right before the end of slot t (see eq 17). For this reason, we use in the right-hand side of eq 16 the excess amount just before the end end (m3), instead of the value at the of the previous slot, Rr,t-1 beginning of the slot (Rr,t-1) like for the unitary resources. Notice that resource inventory at event point t can change due to the start of tasks at t and the end of continuous tasks started at t 1. On the other hand, the start and end of instantaneous tasks (It) affect a single event point. The last term of eq 16 deals with the discrete removal from the system of the prespecified product demands at the depots. It involves final product resources (RFP) and is active if the event point is linked to one of the demand points.
a,0 AV Rr0 ) Vp,s ∀s ∈ S, p ∈ P, r ∈ Rp,s
(9)
end Rr,t ) Rr0|t)1 + Rr,t-1 |r∈RCT + Rr,t-1|r∉RCT +
a,0 r,0 0 If Vp,s ) 0 ∧ 0 < Vp,s < Vs, then RE_P ) 1 ∀s ∈ S, p ∈ P p,s (6)
Otherwise, there is no net accumulation of product inside the pipeline and the active state is fill and empty. r,0 0 a,0 If Vp,s ) 0 ∧ Vp,s ) Vs , then RFE_P ) 1 ∀s ∈ S, p ∈ P p,s (7)
The actual real and accumulated volumes will affect other resources. Considering RIP p,s as the subset that includes the inside pipeline resource linked to product p in segment s (e.g., P1_iP AV as the one that includes the corresponding in Figure 4) and Rp,s accumulated volume resource (e.g., P1_aV), we have
The remaining excess resource variables needing initialization are those linked to the inside refinery and final product resources. FP Being RIR p,rf and Rp,dp the resources linked to product p in refinery rf and depot dp (e.g., P1_iRF1 and P1_DP1 in Figure 2), eqs 10 and 11 result. ir,0 IR Rr0 ) Vp,rf ∀p ∈ P, rf ∈ RF, r ∈ Rp,rf
(10)
Rr0 ) V0p,dp ∀p ∈ P, dp ∈ DP, r ∈ RFP p,dp
(11)
Final product resources can be removed from the system in two different ways. First, there may be a continuous removal rate, drr (m3/h):
µ¯ r,iNi,t-1) +
∑ (µ
r,i
+ µ¯ r,i)Ni,t -
i∈It t*|T|
(∑
∑ (µ
r,iNi,t|t*|T|
i∈Ic
out out Πr,td Yt,td
td∈TD
)|
+
r∈RFP∧t*1
∀r ∈ R, t ∈ T
(16)
Equation 17 accounts for the continuous interaction with tasks. To account for the amount of products continuously removed from the system, we need to multiply the demand rate by the slot duration (last term on the right-hand side). end Rr,t ) Rr,t +
∑λ
r,iξi,t
- drr(Tt+1 - Tt)| r∈RFP
i∈Ic
∀r ∈ RCT, t ∈ T, t * |T|
(17)
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
The first set of timing constraints ensures that the duration of slot t is greater than the processing time of all processing tasks executed at t. The latter is determined by dividing the volume processed by the maximum processing rate. Equation 18 is to be written just for the pump resources (see Figure 2), which are the only equipment resources or states involved in the timing constraints (RTC) since all system movement is due to the pumping tasks. Recall that pump resources have a discrete interaction with pumping tasks, which is why parameter µ j r,i is used instead of λr,i.
∑
Tt+1 - Tt g
µ¯ r,iξi,t
i∈Ic
Fmax i
∀r ∈ RTC, t ∈ T, t * |T|
11497
Figure 9. Single, uniform time grid used by the discrete-time model.
are produced. In this particular case, we do this for the leaving refinery (RLR) and leaving segment resources (RLS), buildup volume (RBV), and global volume resources (RGV), which have been identified in section 3. end Rr,t ) 0 ∀r ∈ RLR ∪ RLS ∪ RBV ∪ RGV, t ∈ T
(18)
(24)
In the case where there are discrete demand points, we need to make the correspondence between event points and demand points. Each demand point td must be associated with a single event point t, which is certainly not the first, the origin of the time horizon; see eq 19. In other words, all demand dates will be event points.
In cases where there is no lower bound on the amount to process on a particular slot, we can make the excess resource variables of the associated equipment units equal to zero since it is possible to have ξi,t ) 0 for Ni,t ) 1. Equation 25 does this for all possible pipeline segment states (RPS) and valves (RVV).
∑
Yout t,td
) 1 ∀td ∈ TD
If a certain event point is in fact a demand point, then they must share the same time coordinate. We can make Tt ) tfxtd out whenever Yt,td ) 1 with eqs 20 and 21. Together with eqs 18 and 19, they ensure that the initial sequence of demand dates is verified when translating it to the event time scale. As an example, demand dates td and td + 1 can respectively be linked to events t and t + 1 or t and t + 2 but not to t and t - 1.
∑
tfxtdYout t,td ∀t ∈ T, t * 1
∑
tfxtdYout t,td + H(1 -
td∈TD
We can also set bounds on the excess resource variables to enforce capacity constraints. On the one hand, the available amount of the inside pipeline (RIP) and accumulated volume resources (RAV) cannot exceed the corresponding segment volume, eq 26. On the other hand, we have the minimum and maximum capacities at the depots to be enforced for final product resources. Since the discrete interactions are all removals of material (e.g., at event point t), maximum inventory occurs at the end of t - 1 while minimum inventory occurs at t. Thus, different sets of excess variables are involved in eqs 27 and 28.
(20) end Rr,t e
td∈TD
Tt e
∑
td∈TD
(21) The two sets of extent variables are related through capacity constraints. Clearly there can only be volume processed if the task is executed (Ni,t > 0). Furthermore, the volume must be lower than the maximum processing rate multiplied by the time horizon, the upper bound on the duration of every slot.
(22)
The switch fill tasks in Figure 4 are only present in the superstructure if the product sequence is not forbidden. However, these are instantaneous tasks and it is possible to execute as many as required on a particular event point. In other words, the sequence p-p′ can be obtained by switching from p to p′′ and then p′′-p′. To avoid this, we ensure the execution of at most one switch fill task in segment s at time t. Note that it is not required to do the same for the switch empty tasks since they also consume the pipeline volume, which is not readily available.
∑N
i,t
e 1 ∀s ∈ S, t ∈ T, t * |T|
(23)
i∈IsSF
The function of some of the resources of the RTN is to ensure that one task immediately follows another. To accomplish this, we need to ensure that such resources are never available in excess; i.e., they will be consumed at the exact same rate they
∑ ∑
p∈P
Yout t,td) ∀t ∈ T, t * 1
c ξi,t e Ni,tFmax i H ∀i ∈ I , t ∈ T, t * |T|
(25)
(19)
t∈T∧t*1
Tt g
Rr,t ) 0 ∀r ∈ RPS ∪ RVV, t ∈ T, t * |T|
end Rr,t e
Vs ∀r ∈ RIP ∪ RAV, t ∈ T
(26)
FP Vmax p,dp ∀r ∈ R , t ∈ T, t * |T|
(27)
s∈S IP ∪RAV r∈Rp,s p,s
∑ ∑
p∈P dp∈DP FP r∈Rp,dp
Rr,t g
∑ ∑
FP Vmin p,dp ∀r ∈ R , t ∈ T
(28)
p∈P dp∈DP FP r∈Rp,dp
Finally, there are the bounds on the timing variables, with eq 29 reflecting Figure 7. T1 ) 0;
T|T| ) H;
0 e Tt e H ∀t ∈ T
(29)
4.4. Remarks. The RTN representation given in Figures 2 and 4 is not a function of the time representation used by the mathematical model. In fact, it can also be used as the basis for a discrete-time model. Maravelias and Grossmann23 have proved that a state-task network continuous-time model for batch plants reduces to a discrete-model under the assumption of a fixed time grid and constant processing times that are multiples of the uniform interval length. For continuous plants, where tasks span just one time slot, Castro and co-workers14,15 have shown that RTN discrete and continuous-time formulations are very similar. The required minor changes are the following. With a uniform time grid of size δ (see Figure 9) we know a priori the location of the discrete events resulting from the out out can be replaced by Πr,t in eq multiple due dates. Hence Πr,td 16 and the extra set of binary variables can be dropped. Then,
11498
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 10. Initial inventory in pipeline, refinery supply, and depot demand data for Ex0.
the last term of eq 17 becomes the product of two parameters; see eq 31. Finally, the relation between the binary and continuous extent variables can be tightened; see eq 32. The full set of constraints of the discrete-time model is then eqs 15, 23-28, and 30-32. end Rr,t ) Rr0 | t)1 + Rr,t-1 | r∈RCT + Rr,t-1 | r∉RCT +
∑ (µ
r,iNi,t | t*|T|
+ µ¯ r,iNi,t-1) +
i∈Ic
∑ (µ
r,i
+ µ¯ r,i)Ni,t -
i∈It t*|T|
out | r∈RFP∧t*1 ∀r ∈ R, t ∈ T Πr,t end Rr,t ) Rr,t +
∑λ
r,iξi,t
(30)
- drrδ| r∈RFP ∀r ∈ RCT, t ∈ T, t * |T|
i∈Ic
(31) c ξi,t e Ni,tFmax i δ ∀i ∈ I , t ∈ T, t * |T|
(32)
One reviewer asked the question if it is possible with the RTN approach to ensure that any pipeline segment stoppage only occurs if a single product is stored inside, so that the idle time does not increase the size of the interfaces. In principle, meaning that this was not actually tried out, this could be achieved with the same superstructure and the following actions (recall Figure 4): (i) set minimum processing rates for the fill, empty, move, and fill and empty tasks using a constraint similar to eq 22 (this will ensure consumption of the alternative pipeline segment states only if there is pumping going on); (ii) change eq 25 by removing from the domain the pipeline segment states corresponding to the case of a single product inside the pipeline, i.e., states fill and empty (for the product inside the segment) and do nothing (for the others). In other words, during stoppage, only the latter states would be allowed to exist in excess. 5. Illustrative Examples We now consider three example problems from the literature to illustrate the capabilities of the new continuous-time formulation for the short-term scheduling of multiproduct pipeline systems. They range from linear to treelike structures, feature one or two refineries, feature three to four products possibly with a few forbidden sequences, and feature a few depots with either discrete (at the end of the time horizon) or continuous demand patterns. Focus is set on identifying the active tasks at the refineries and depots as well as within the pipeline segments, together with the inventory profiles. The computational statistics and comparison to previous approaches can be found in section 6. 5.1. Example 0 (Ex0). Before moving on to the more interesting examples, we present a trivial problem consisting of a single refinery, pipeline segment, and depot, as well as two products (Figure 10). Its purpose is to show the required model parameters, as illustrated in Figure 11. A total of 22 tasks and 26 resources are required that are properly identified next to the names used in Figures 2 and 4. To make the diagram
simpler, some tasks and resources have been duplicated, with each copy filled in gray. In order to meet the depot demands over a time horizon H ) 24 h, three event points suffice to find the optimal solution, which involves pumping 9000 m3 of P1. The initial inventories at the refinery and depot are specified as initial availability (Rr0) for resources 1, 2, 25, and 26 (see eqs 10 and 11). For the latter (located on the right of Figure 11), we also need the continuous demand rate drr, calculated through eq 12. In the beginning of the time horizon the following resources are also available: (i) resource 5, the refinery pump; (ii) resource 10, since there is no P1 inside the pipeline segment (eq 3); (iii) resource 13, since the pipeline is completely filled with P2 and ready to discharge (eq 7); (iv) resource 20, the inside pipeline resource for P2, in an amount equal to the pipeline volume, 5000 m3 (eq 8). As for the nonzero structural parameters, these are mostly -1 or 1 depending on the direction of the arrow linking tasks with j r,i, or λr,i) is resources. In turn, the actual parameter involved (µr,i/µ a function of the line type, respectively dashed or solid. The exceptions occur for the moving tasks (resources 4 and 10) that deal with a specified minimum, 100 m3, of the corresponding inside pipeline resource (resources 16 and 20), and for the switch empty A/B tasks (resources 17-18 and 21-22) that consume at their start 5000 m3 of the product’s accumulated volume resource (respectively resources 23 and 19). 5.2. Example 1 (Ex1). The first example is taken from Cafaro and Cerda´6 and involves finding the optimal schedule for a linear pipeline system featuring two refineries, three products, and three depots, over an H ) 120 h time horizon. The pipeline state at the start of the horizon is given in Figure 12 together with the remaining problem data. Notice that the pipeline segment volumes and refinery pumping flow rates are normalized and the total demands at the depots can be found inside the vessels. The latter are modeled as an instantaneous removal of material at the end of the time horizon (see eq 13). One important operating condition is that there can be no simultaneous movement in segment S2 and pumping in refinery RF2. In other words, the two refineries can only be operating simultaneously if all the outlet flow from segment S1 is going into depot DP1. We consider two cases. In Ex1a, there is an unlimited supply at the refineries. In contrast, in Ex1b, the refineries are product dedicated (RF1 handles P1 and P2 while RF3 deals with P3) and the supplies are limited. In both, the optimal solution corresponds to the pumping of 140 units. The optimal schedule for Ex1a is given in Figure 13, which shows the pumping tasks at the refineries as well as the send to depot tasks (recall Figure 2). Inside the rectangles we give the amount handled by the task. To construct the charts, it was assumed that the pumping tasks were always processed at the maximum processing rate, followed by an idle period whenever the duration of the slot is greater. This corresponds to the middle option in Figure 8, and the reason is to better identify the remaining pumping capacity. In practice, however, the best option would be the top alternative in Figure 8, which reduces settling effects by ensuring that the pipeline never remains idle; i.e., we want to keep moving the products on turbulent flow.24 With respect to the bifurcations (into two segments or into a depot and a segment), we assume a constant volume ratio throughout the duration of the pumping task; i.e., when the two tasks are active in the same slot, they take place simultaneously. From the optimal solution one can see that there are two independent subsystems, which is hardly surprising. Refinery RF1 meets the demand of DP1, while RF2 meets the demands of DP2 and DP3. As a result, there is no movement inside segment S2. Overall, there is considerable free capacity in the
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
11499
Figure 11. Full set of resources, tasks, and associated RTN structural parameters for Ex0.
Figure 12. Initial inventory in pipeline, refinery supply, and depot demand data for Ex1.
system since RF1 is only pumping for 25 h and RF2 is pumping for 91.667 h. In Ex1b, the system becomes fully integrated with the majority of the workload shifting to RF1 as a consequence of
RF2 being restricted to P3. There is now movement in all pipeline segments, with Figure 14 showing the active tasks inside them. Initially, segment S1 is entirely filled with P1, see Figure 12, so when RF1 starts pumping P2, there will be the
11500
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 13. Optimal solution for Ex1a.
execution of a filling task for P2 and of an empty task for P1; recall Figure 4. Since the amount handled is 20 units, the segment volume, all P1 will have been replaced by P2 in S1 at the end of the pumping task. Since RF1 keeps pumping P2 (second instance of the task, now corresponding to 10 units), the next active task in S1 is a fill and empty task that is represented as a gradient-filled rectangle. Concerning the other end of the pipeline, DP1 and DP2 take all the material until the end of the third day, so there will be no movement in segment S4 until roughly the 72 h mark. At that time, segment S3 is already full with P2, like S4, so the first task to be executed in the latter is fill and empty for P2. 5.3. Example 2 (Ex2). Example 2 is adapted from Rejowski and Pinto1 and involves a linear pipeline with a single refinery. Initial inventory data in segments S3 and S4 were slightly changed to overcome the model limitation of having a single batch of the same product per segment, as discussed in section 3.3. Four products are involved, and there are four forbidden sequences: P2-P3, P3-P2, P3-P4, and P4-P3. The other difference from Ex1 is that we assume a continuous removal of material from the depots, at a constant rate throughout the duration of the time horizon, H ) 75 h (see eq 11). Besides the data in Figure 15, we need to know the initial, minimum, and maximum capacities at the depots, as well as total demands (Table 1). The optimal solution for Ex2 is given in Figure 16, and there are a couple of aspects worth discussing. First, the majority of the demand at the depots is met by depleting the initial inventory. This can be seen from the lack of active send to depot tasks for DP4 and from the inventory profiles for DP1, where there are strong reductions in the amount in storage for P1-P3. While this is acceptable for the current scheduling horizon, one may run into problems in subsequent periods, where it may be impossible to reach the minimum capacity constraints given the delay between start pumping and the product reaching the final destination. Clearly, one should use the remaining free capacity at the refinery to prevent such large differences between the initial and final inventories. Overall, modeling demands as continuous rather than discrete removals of material brings us much closer to reality. Second, while it seems that P3 follows P4 in the schedule (this is a forbidden sequence), there is a small amount (100 m3) of P1 in between. This was the minimum batch volume defined for the movement task in Figure 4. As a consequence, there are sometimes three products on a particular pipeline segment as can be seen from the number of simultaneous tasks inside segment S1 or, more easily, from its inventory profiles. The latter confirm that there is no net accumulation or depletion of material inside any given segment. It is important to highlight that inventory profiles inside the segments and at the depots
Figure 14. Optimal solution for Ex1b.
can be easily generated from the values of the excess resource variables corresponding to the inside pipeline and inside depot resources, respectively. Bear in mind that a minimum batch volume of 100 m3 is not realistic since interfaces between two products can grow as much as 200 m3, as highlighted by one reviewer. Using the latter value leads to an optimal solution that is essentially the same despite pumping a total volume equal to 22 200 instead of 22 100 m3. 5.4. Example 3 (Ex3). Example 3 is adapted from the recent paper by MirHassani and Jahromi,4 which deals with a complex treelike structure featuring two branches from the main pipeline. The structural and initial inventory data are given in Figure 17, while the demands and capacities at the depots can be found in
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
11501
Figure 15. Initial inventory in pipeline for Ex2. Table 1. Capacities and Product Demands at the Depots for Ex2 (m3)
demand
t (dp,dp )
max maximum (Vp,dp )
min minimum (Vp,dp )
0 initial (Vp,dp )
P1 P2 P3 P4 P1 P2 P3 P4 P1 P2 P3 P4 P1 P2 P3 P4
DP1
DP2
DP3
DP4
DP5
10000 7000 6000 6000 40000 40000 7000 40000 9000 9000 1000 9000 19000 18000 5000 12000
11000 9000 4000 5000 40000 40000 7000 40000 9000 9000 1000 9000 23000 21000 6500 14000
12000 10000 4000 5000 40000 40000 7000 40000 9000 9000 1000 9000 20000 18000 6000 19000
12000 8000 5000 40000 40000 7000 40000 9000 9000 1000 9000 24000 18000 6000 19000
15000 10000 2000 5000 40000 40000 7000 40000 9000 9000 1000 9000 19000 18000 6000 17000
Table 2. When compared to the original problem, the differences are in the initial condition of segments S5 and S7 and in the lower demand of P4 at DP4, which enables it to be met exclusively by its initial inventory. These changes are motivated by the need to reduce the number of events required to generate the solution, which have a great influence on computational effort as will be seen in section 6, besides the already mentioned model restriction of a single batch per product on a particular segment. The forbidden product sequences are now P3-P4 and P4-P3. As in Ex2, the time horizon (H ) 90 h) seems too long for the maximum pumping rate at the refinery (1000 m3/h). In the optimal schedule given in Figure 18, tasks executed in the first four time slots are next to each other, while those of the fifth slot are executed almost 2 days later. Naturally, there is no reason for this behavior since those tasks could follow immediately after the others. It is just evidence of solution degeneracy. During the first time slot, the 3000 m3 of P3 that is pumped is being used to drive P1 into DP6. Interestingly, 2000 m3 would suffice to keep the inventory above the minimum at the end of the time horizon. The reason is that DP5 needs just 2000 m3 of P3, so the other 3000 m3 initially in S7 moves into S8 and stays there until the end. The second time slot is the longest, and all other depots are being serviced. Notice that the extent of the pumping task in the Gantt chart is equal to the sum of the extents of the send to depot tasks, which confirms that the overall mass balance is met. It is interesting to observe that three different products are arriving at their destinations (P1, P3, and P4) by pumping just one (P3). Overall, DP4 receives the largest amount of material followed by DP2 and DP1. 6. Computational Results After showing that the proposed mathematical formulation can cope with a wide variety of problem features, we now focus on computational performance. In addition to the test cases just discussed, three simpler versions of Ex3 are considered to
Figure 16. Optimal solution for Ex2.
illustrate the increase in size and complexity with the increase in the number of pipeline segments and depots. These are named Ex4a, Ex4b, and Ex4c. The algorithm to systematically generate the tasks and resources, and associated RTN parameters, was implemented in GAMS 23.2 together with the mixed integer linear programming formulation. This was solved using CPLEX 12.1 with two parallel threads and default options except for the relative optimality tolerance equal to 10-6. The hardware consisted of a laptop with an Intel Core2 Duo T9300 2.5 GHz processor with 4 GB of RAM running Windows Vista Enterprise. 6.1. Performance Indicators. The computational statistics are given in Table 3. The second column gives the number of process entities in terms of numbers of products, segments, depots, and refineries. The fifth column gives the size of the time grid for a particular optimization run in terms of number
11502
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 17. Initial inventory in pipeline for Ex3. Table 2. Capacities and Product Demands at the Depots for Ex3 (m3)
t demand (dp,dp )
P1 P2 P3 P4 max maximum (Vp,dp ) P1 P2 P3 P4 min minimum (Vp,dp ) P1 P2 P3 P4 0 initial (Vp,dp ) P1 P2 P3 P4
DP1
DP2
DP3
DP4
DP5
DP6
6000 11500 7000 12500 100000 100000 100000 100000 9000 12000 4000 8000 15000 22000 8000 19000
7000 3000 5000 3000 8000 20000 23000 18000 1000 11000 12000 5000 3000 14000 17000 6000
1000 3750 3000 4000 100000 100000 100000 100000 2000 10000 9000 5000 3000 12000 12000 9000
4200 6000 9000 14000 100000 100000 100000 100000 4000 1000 11000 14000 7000 5000 20000 28000
12000 10000 4500 1000 20000 24000 9000 21000 9000 8000 2000 11000 19000 18000 5000 12000
14000 14000 1000 5000 30000 40000 13000 37000 11000 7000 5000 9000 23000 21000 6000 14000
of event points. Columns RMIP and MIP give the values of the relaxed solution and final solution from the MILP, with their difference defining the integrality gap. The first aspect to highlight is the relative large number of RTN entities that are required to model the process system, with a maximum of 329 tasks and 508 resources for Ex3. It is a direct consequence of the sequence dependent changeover tasks used to change between products on the two ends of a pipeline segment, which tend to increase such entities by roughly 1 order of magnitude,21 when compared to other previously addressed industrial cases.14,15 A typical behavior of time grid based continuous-time formulations is the strong dependence of the computational statistics on the number of event points. The most important aspect is that there are a minimum number of events for which
the problem is feasible. This is related not only to the process structure but also to the other data, as can be seen comparing Ex1a to Ex1b. With two separate subsystems six events are enough, while EX1b requires the products to go through more segments triggering two more events for the same total displacement of 140 units. An interesting result is that the minimum number of events for which the problem becomes feasible provided the real optimal solution for all examples but Ex4a. From past experience, we know that the search from the global optimum solution may involve up to a few iterations, with gradual improvements in the objective function past the minimum. The explanation for this awkward behavior is that the objective function is now the minimum amount of volume pumped from the refinery, which will generally involve fewer events than pumping larger volumes. This insight also serves to explain why it is also particularly difficult to find the first feasible solution. The associated drawback is that one will not know if the model is in fact infeasible for the specified number of events or if it is just too hard to find a feasible solution. This is the practical reason behind the simplifications that were made for Ex3. The motive why Ex4 has a higher computational effort than Ex3, having less complexity in the system topology, is also related to the difficulty of finding feasible solutions. While the solution of the relaxed problem increased to 39 250 m3 in 100 CPUs, the first feasible solution was 42 000, which then improved to 39 500 around 250 CPUs and eventually to the optimal solution at 2797 CPUs. Bear in mind that branch and bound is a heuristic process, so this behavior is not uncommon. With the exception of Ex1, the problems feature a large integrality gap. This is a direct consequence of the several
Table 3. Computational Statistics for Proposed Continuous-Time Formulationa (P,S,DP,RF)
resources
tasks
|T|
DV
SV
SE
RMIP
MIP
CPUs
nodes
Ex1a
(3,4,3,2)
133
171
Ex1b
(3,4,3,2)
133
171
Ex2
(4,5,5,1)
214
284
Ex3 Ex4a
(4,8,6,1) (4,4,3,1)
329 169
508 256
Ex4b
(4,5,4,1)
210
320
Ex4c
(4,7,5,1)
288
444
5 6 7 8 6 7 6 4 5 6 5 6 6
688 860 1032 1204 1420 1704 2540 768 1024 1280 1280 1600 2220
1939 2390 2841 3292 3896 4632 6301 2001 2611 3221 3252 4012 5510
1279 1565 1851 2137 2500 2957 3800 1244 1602 1960 1991 2436 3324
140 140 140 140 13 000 13 000 23 450 14 750 14 750 14 750 17 950 17 950 21 450
infeasible 140 infeasible 140 infeasible 22 100 39 250 22 500 17 750 17 750 23 250 23 250 39 250
2.80 4.5 16.1 27.1 87.9 299 901 0.56 11.4 40.7 16.4 347 2797
499 1220 6984 4958 0 1104 1109 1102 6229 14795
a
DV ) discrete variables; SV ) single variables; SE ) single equations.
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Figure 18. Optimal solution for Ex3.
resources that are involved in modeling plug flow inside the pipeline segments (recall Figure 4). While a few attempts were made to derive alternative formulations, none was able to close the gap and the one proposed exhibited the best overall performance. Since the branch-and-bound search procedure will normally be more time-consuming with an increase in the integrality gap, there is still plenty of room for improvements. Bear in mind that there may be more than a suitable process representation, with some being more computationally efficient than others as also noted by Chen and Chang when tackling heat-integrated batch plants.25
11503
6.2. Comparison to Alternative Formulations. It is now worthwhile to provide a comparison to previous continuoustime and continuous-volume formulations that have been used to solve pipeline scheduling problems. We consider the ones by Cafaro and Cerda´2,6 and MirHassani and Jahromi.4 These are suitable for linear systems with single2 or multiple sources6 or treelike systems with a single refinery.4 Despite the structural differences in scope, such formulations can be viewed as being conceptually similar. The number of pumping batches is one of the main elements of the formulations, and there are variables that give their completion time, pumping length, and volume. Then, a set of binary variables links batches to products. In the new formulation a different concept is used. Now, there are pumping tasks linked to products instead of batches, meaning that a so-called batch may involve the execution of a couple of tasks, as previously seen in section 5. As a consequence, there are no explicit variables to characterize the batch, with such information being implicit in the tasks continuous extent variables as well as in the timing variables. These differences have a profound impact on the computational statistics as can be seen in Table 4. The first aspect is the larger number of binary variables employed by our new approach, which ranges between factors of 4 and 8. It can be explained by the use of several sequence dependent changeover tasks to switch between pipeline states, whereas in the other models continuous variables4 coupled with traditional constraints2 suffice to prevent forbidden sequences. In contrast, the total number of constraints is roughly similar. In terms of computational effort, it is clear that previous approaches are faster considering the advances in hardware and the MILP solver together with the use of two parallel threads in the current work. In addition, we made changes to Ex2 and Ex3 to reduce the complexity, so the performance gap is definitely higher than the one apparent from direct comparison of the CPUs. On the other hand, the use of an alternative objective function is not expected to play a major role. Overall, the advantage here is the high generality of the new approach, which is not structure dependent unlike the more efficient, problem-specific formulations. 6.3. Continuous vs Discrete Time. To conclude the analysis, we show in Table 5 the performance of the mathematical model for a discrete-time representation, where the time horizon is divided into equal length slots of size δ. The results show that this is very restricting for most example problems. While the number of event points that ensures feasibility increased every single time, doubling the number may even not enough (e.g., Ex1b). This makes the iterative search procedure for the global optimal solution considerably harder. In general, the performance of the discrete-time formulation was always worse than its continuous-time counterpart. This contradicts the results in Castro et al.14 that also involve a short-term scheduling problem of a continuous system. In there, the discrete-time formulation was able to handle 1 order of magnitude more slots and performed much better. Perhaps the difference in behavior is due to the fact that feasible solutions are now considerably harder to obtain and no reduction in integrality gap can be observed when switching from continuous to discrete time. It serves to show that there is no such thing as an overall best time representation. 7. Conclusions This paper has focused on the modeling of a multiproduct pipeline system with multiple sources and sinks connected by
11504
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010
Table 4. Computational Results from Different Sources
hardware solver discrete variables single variables equations CPUs discrete variables single variables equations CPUs discrete variables single variables equations CPUs
Ex1b
Ex2
Ex3
a
Cafaro and Cerda´6
Cafaro and Cerda´2
MirHassani and Jahromi4
this work
Intel 2.8 GHz CPLEX 11.0 297 1076 2180 12.5a -
Pentium IV, 2 GHz CPLEX 9.1 214 3000 2349 34.98a,b -
Pentium 4, 2.4 GHz CPLEX 10.0 444 2095 3810 48a,b,c
Intel T9300, 2.5 GHz CPLEX 12.1 1204 3292 2137 27.1 1704 4632 2957 299 2540 6301 3800 901
Minimization of total operating cost. b Slightly different initial inventory in pipeline. c Slightly different product demands.
tion of having a single batch of a particular product on the same pipeline segment and testing other objective functions.
Table 5. Computational Statistics for Proposed Discrete-Time Formulation δ (h) |T| Ex1a 24 20 Ex1b 10 8 Ex2 7.5 6.25 Ex3 11.25 Ex4a 15 11.25 Ex4b 15 11.25 Ex4c 11.25
6 7 13 16 11 13 9 7 9 7 9 9
DV
SV
SE
RMIP
MIP
CPUs
nodes
855 1026 2052 2565 2840 3408 4064 1536 2048 1920 2560 3552
2379 2828 5522 6869 7565 9035 9874 3824 5042 4765 6283 8633
1544 1826 3518 4364 4775 5687 5874 2312 3026 2875 3763 5137
140 140 140 140 13 000 13 000 23 450 14 750 17 950 21 450
infeasible 140 infeasible infeasible infeasible infeasible 40 125 infeasible 17 750 infeasible 23 250 39 250
11.3 13.1 313 83.7 146 503 7200a,b 0.41 460 0.46 1308 7200a,c
524 7050 3645 4506 7257
a Maximum resource limit. possible solution ) 37 250.
b
Best possible solution ) 37 750.
c
Best
a treelike structure. The challenge has been to find an appropriate set of virtual entities that are the tasks and resources of the resource-task network process representation and linking them in a systematic way to the real system entities that are the products, refineries, and depots. The overall system has been constructed as a collection of elementary building blocks with the pipeline segment being the most important. It comprises fill, empty, and moving tasks together with changeover tasks that change the state of the pipeline while preventing forbidden product sequences. More importantly, it features three different sets of volume resources per product to successfully model plug flow. The advantage of relying on a unified framework for process representation is that, once the system model is found, one can take advantage of existing general mathematical formulations to solve the problem at hand. A very efficient single time grid continuous-time formulation has been used that allows handling both discrete and continuous product demands. The performance of the new formulation was tested on a set of three example problems taken from the literature and compared to three conceptually different approaches that are suitable for specific subclasses of the problem considered. The results have clearly shown that the higher generality comes at a cost in terms of efficiency, which rises with the number of events that are needed to represent the solution. Thus, it is fair to say that problems of industrial relevance featuring long time horizons should not be tackled with the full-space approach. Instead a decomposition approach involving a rolling-horizon algorithm8,9,26 should be sought. Future work will also involve overcoming the model limita-
Acknowledgment The author gratefully acknowledges financial support from European project COMET (FP7-ENERGY-2009-1 Grant Number 241400). Nomenclature Sets DP, dp ) depots DY, dy ) days of the week Hr, h ) hours of the day I, i ) tasks Ic ) continuous processing tasks Ipump ) pumping tasks at refineries IsSF ) instantaneous switch fill tasks in segment s It ) instantaneous tasks P, p ) products R, r ) resources RF, rf ) refineries RAV ) accumulated volume resources AV Rp,s ) accumulated volume resource linked to product p in segment s RBV ) buildup volume resources RCT ) resources continuously produced/consumed RFP ) final product resources FP Rp,dp ) final product resource linked to product p at depot dp GV R ) global volume resources RIP ) inside pipeline resources IP Rp,s ) inside pipeline resource linked to product p in segment s RIR ) inside refinery resources IR Rp,rf ) inside refinery resource linked to product p at refinery rf LR R ) leaving refinery resources RLS ) leaving segment resources RPS ) resources corresponding to pipeline segment states RTC ) resources involved in timing constraints RVV ) connecting valve resources S, s ) pipeline segment T, t ) time slots (event points) TD, td ) demand points TDdy,hr ) demand point associated with hour hr of day dy Parameters t dp,dp ) total demand product p from depot dp during the time horizon (m3)
Ind. Eng. Chem. Res., Vol. 49, No. 22, 2010 dp,dp,dy,hr ) instantaneous removal of product p from depot dp at hour hr of day of the week dy (m3) drr ) continuous removal rate of resource r throughout the duration of the time horizon (m3/h) H ) time horizon Rr0 ) availability of resource r at start of time horizon tfxtd ) absolute time of demand point td (h) 0 Vp,dp ) amount of product p in storage at depot dp at start of time horizon (m3) max Vp,dp ) maximum capacity of product p at depot dp (m3) min Vp,dp ) minimum allowable capacity of product p at depot dp (m3) a,0 Vp,s ) accumulated volume (after product p entered the pipeline) inside segment s at the start of time horizon (m3) Vir,0 p,rf ) amount of product p inside refinery rf at start of time horizon (m3) r,0 Vp,s ) amount of product p inside pipeline segment s at start of time horizon (m3) Vs ) capacity of pipeline segment s (m3) δ ) duration of time slots in uniform time grid (discrete-time model) λr,i ) continuous interaction of resource r during execution of task i µr,i ) discrete interaction of resource r with task i at its start µ j r,i ) discrete interaction of resource r with task i at its end Πout r,t ) amount removed from the system of resource r at time point t (discrete-time model) out Πr,td ) amount of resource r removed from the system at demand point td Fimax ) maximum processing rate of task i (m3/h) Frfmax ) maximum pumping flow rate at refinery rf (m3/h) Variables Ni,t ) binary variable identifying the execution of task i during slot t or instantaneously at event t Rr,t ) excess amount of resource r at event point t end Rr,t ) excess amount of resource r immediately before the end of interval t (m3) Tt ) absolute time of event point t (h) out Yt,td ) binary variable making the correspondence between event point t and demand point td ξi,t ) amount handled by task i during slot t or instantaneously at event t (m3)
Literature Cited (1) Rejowski, R.; Pinto, J. M. Scheduling a Multiproduct Pipeline System. Comput. Chem. Eng. 2003, 27, 1229. (2) Cafaro, D. C.; Cerda´, C. Optimal Scheduling of Multiproduct Pipeline Systems using a Non-Discrete MILP Formulation. Comput. Chem. Eng. 2004, 28, 2053. (3) Relvas, S. Optimal Pipeline Scheduling and Inventory Management of a Multiproduct Oil Distribution Centre. Ph.D. Thesis, Instituto Superior Te´cnico, Universidade Te´cnica de Lisboa, Portugal, 2008. (4) MirHassani, S. A.; Jahromi, H. F. Scheduling Multi-Product TreeStructure Pipelines. Comput. Chem. Eng. [Online early access]. DOI: 10.1016/j.compchemeng.2010.03.018. Published Online: April 8, 2010. (5) Relvas, S.; Matos, H. A.; Barbosa-Po´voa, A. P. F. D.; Fialho, J.; Pinheiro, A. S. Pipeline Scheduling and Inventory Management of a Multiproduct Distribution Oil System. Ind. Eng. Chem. Res. 2006, 45, 7841.
11505
(6) Cafaro, D. C.; Cerda´, C. Optimal Scheduling of Refined Product Pipelines with Multiple Sources. Ind. Eng. Chem. Res. 2009, 48, 6675. (7) Rejowski, R.; Pinto, J. M. A novel continuous time representation for the scheduling of pipeline systems with pumping yield rate constraints. Comput. Chem. Eng. 2008, 32, 1042. (8) Cafaro, D. C.; Cerda´, C. Dynamic Scheduling of Multiproduct Pipelines with Multiple Delivery Due Dates. Comput. Chem. Eng. 2008, 32, 728. (9) Dimitriadis, A. D.; Shah, N.; Pantelides, C. C. RTN-based Rolling Horizon Algorithms for Medium Term Scheduling of Multipurpose Plants. Comput. Chem. Eng. 1997, 21, S1061. (10) Relvas, S.; Matos, H. A.; Barbosa-Po´voa, A. P. F. D.; Fialho, J. Reactive Scheduling framework for a Multiproduct Pipeline with Inventory Management. Ind. Eng. Chem. Res. 2007, 46, 5659. (11) Me´ndez, C. A.; Cerda´, J.; Grossmann, I. E.; Harjunkoski, I.; Fahl, M. State-of-the-art Review of Optimization Methods for Short-Term Scheduling of Batch Processes. Comput. Chem. Eng. 2006, 30, 913. (12) Kondili, E.; Pantelides, C. C.; Sargent, R. A General Algorithm for Short-Term Scheduling of Batch OperationssI. MILP Formulation. Comput. Chem. Eng. 1993, 17, 211. (13) Pantelides, C. C. Unified Frameworks for the Optimal Process Planning and Scheduling. In Proceedings of the Second Conference on Foundations of Computer Aided Operations; Cache Publications: New York, 1994; p 253. (14) Castro, P. M.; Harjunkoski, I.; Grossmann, I. E. New ContinuousTime Scheduling Formulation for Continuous Plants under Variable Electricity Cost. Ind. Eng. Chem. Res. 2009, 48, 6701. (15) Castro, P.; Westerlund, J.; Forssell, S. Scheduling of a Continuous Plant with Recycling of Byproducts: a Case Study from a Tissue Paper Mill. Comput. Chem. Eng. 2009, 33, 347. (16) Wassick, J. Enterprise-wide Optimization in an Integrated Chemical Complex. Comput. Chem. Eng. 2009, 33, 1950. (17) Sasikumar, M.; Prakash, P. R.; Patil, S. M.; Ramani, S. Pipes: a Heuristic Search Model for Pipeline Schedule Generation. Knowledge-Based Syst. 1997, 10, 169. (18) Castro, P.; Barbosa-Po´voa, A.; Matos, H.; Novais, A. Simple Continuous-time Formulation for Short-Term Scheduling of Batch and Continuous Processes. Ind. Eng. Chem. Res. 2004, 43, 105. (19) Shaik, M.; Floudas, C. A. Unit-specific event-based continuoustime approach for short-term scheduling of batch plants using RTN framework. Comput. Chem. Eng. 2008, 32, 260. (20) Castro, P.; Barbosa-Po´voa, A.; Novais, A. Simultaneous Design and Scheduling of Multipurpose Plants Using Resource Task Network Based Continuous-Time Formulations. Ind. Eng. Chem. Res. 2005, 44, 343. (21) Castro, P. M.; Novais, A. Q. Scheduling Multistage Batch Plants with Sequence-Dependent Changeovers. AIChE J. 2009, 55, 2122. (22) Gime´nez, D.; Henning, G.; Maravelias, C. A Novel Network-based Continuous-Time Representation for Process Scheduling: Part I. Main Concepts and Mathematical Formulation. Comput. Chem. Eng. 2009, 33, 1511. (23) Maravelias, C.; Grossmann, I. E. On the Relation of Continuousand Discrete-Time State-Task Network Formulations. AIChE J. 2006, 52, 843. (24) Hane, C. A.; Ratliff, H. D. Sequencing Inputs to Multi-Commodity Pipelines. Ann. Oper. Res. 1995, 57, 73. (25) Chen, C. L.; Chang, C. Y. A Resource-Task Network Approach for Optimal Short-term/Periodic Scheduling and Heat Integration in Multipurpose Batch Plants. Appl. Therm. Eng. 2009, 29, 1195. (26) Castro, P. M.; Harjunkoski, I.; Grossmann, I. E. Optimal Scheduling of Continuous Plants with Energy Constraints. Comput. Chem. Eng. [Online early access]. DOI: 10.1016/j.compchemeng.2010.05.008. Published Online: May 19, 2010.
ReceiVed for reView May 14, 2010 ReVised manuscript receiVed August 17, 2010 Accepted August 25, 2010 IE1010993