A General Modeling Framework for the Long-Term Scheduling of

Article pubs.acs.org/IECR

A General Modeling Framework for the Long-Term Scheduling of Multiproduct Pipelines with Delivery Constraints Hossein Mostafaei* and Alireza Ghaffari Hadigheh Department of Applied Mathematics, Azarbiajan Shahid Madani University, Tabriz, Iran S Supporting Information *

ABSTRACT: Common carrier pipelines are one of the most economic modes for transportation of petroleum refined products over land, especially when huge amounts of these products have to be pumped toward long-distance terminals. This paper introduces a novel Mixed Integer Linear Programming (MILP)-based approach for the long-term scheduling of a real world multiproduct pipeline connecting a unique refinery to several distribution centers. This approach allows consideration of multiple due dates for demands at period ends, flow rate limitation on pipeline segments, and simultaneous deliveries at distribution centers. The proposed model results in substantial reduction in pump operation and maintenance costs in comparison with the available models. Computational results and data are reported.

1. INTRODUCTION Several million barrels of oil derivatives such as gasoline, kerosene, aviation fuel, diesel, home heating oil, and liquid petroleum gases (LPG) are moved daily around the world in imports and exports. These products are usually carried by roads, railroads, vessels, and pipelines, among which the latter is the most reliable and economical one for the transportation of huge amounts of oil products between production areas and distribution centers. Suffice it to say that common carrier pipelines convey nearly 70% of petroleum products in the United States.1,2 The pipelines must always be completely filled with products; hence, assuming incompressible fluids, a basic pipeline operation can be described as pumping a given quantity of product into the pipeline at the origin and simultaneously receiving the same volume of product at the other end of the pipeline. In the last two decades, the scheduling of dispatching refined products through the multiproduct pipeline has been studied from various perspectives, mostly based on optimization techniques such as knowledge-based search heuristics and mixed integer linear programming (MILP) formulations. The MILP formulations can be classified as discrete and continuous time-based models. In the discrete approach, the pipeline is divided into packages of uniform sizes, each transporting a single product. In this approach, the planning horizon is also divided into time intervals having a fixed duration. These limitations are not considered in the continuous methods. The most important aim of the common-carrier pipeline scheduling is to provide products for different demand points at the right time with the lowest operational costs. This work is a very complex task with many operational constraints. Milidiu and Liporace3 have demonstrated that the problem of oil product transportation in a multiproduct pipeline is an NPcomplete problem when interface constraints are taken into account. Sasikumar et al.4 developed an expert system for the scheduling problem that concerns a pipeline connecting a unique refinery to several distribution centers. The authors © 2014 American Chemical Society

considered the allowable and forbidden sequence of products inside the pipeline which present a complicating factor concerning the feasibility of the schedules generated. Their objective function was to minimize the pumping cost, together with the product contamination inside the pipeline over a monthly time horizon. One of the first approaches was presented by Shah,5 who developed an MILP discrete time formulation for the transportation of crude oil from a port to a refinery through a pipeline where the planning horizon was divided into intervals of an equal duration. Due to nonlinearity issues, he divided the problem into two models. The first one considers the loading of crude oil from the docks to the storage tanks in a port, and the second one determines these operations at the refinery to unload the crude oil from tanks into the distillation units. Their objective was to minimize the crude oil heels of charging tanks over a monthly time horizon. Techo and Holbrook6 presented a simplified model for the scheduling of a multiproduct pipeline system with multiple destinations to minimize the interface losses. They mentioned that the reprocessing cost concerning the interfaces losses is very high. Hane and Ratliff7 presented an MILP formulation based on discrete time representation for transporting commodities in a large multiproduct pipeline from a unique refinery to several destinations to minimize the pumping and maintenance costs. The problem was broken up into some subproblems that can be easily solved by a branch-and-bound algorithm. Reddy et al.8 introduced a continuous time formulation with the purpose of scheduling of activities in a refinery that receives crude oil from vessels through a pipeline. Neiro and Pinto9 depicted a mathematical representation for modeling petroleum supply chains. The system comprises pipelines, crude oil suppliers, oil refineries, and distribution centers. The pipelines convey crude oil from petroleum Received: Revised: Accepted: Published: 7029

November 10, 2013 February 19, 2014 April 2, 2014 April 2, 2014 dx.doi.org/10.1021/ie4038032 | Ind. Eng. Chem. Res. 2014, 53, 7029−7042

Industrial & Engineering Chemistry Research

Article

moving in either direction. The pipeline system connects a harbor to an inland refinery. In a subsequent work,19 the authors considered the same problem with improved efficiency using a combined approach based on the usage of constrained logic programming and MILP. This procedure reduces the computational burden of the previous work by 1 or 2 orders of magnitude. Mirhassani and Ghorbanalizadeh20 presented an MILP discrete time formulation for the scheduling of a realworld tree-structure pipeline network. For their problem, a novel MILP formulation based on a continuous time representation was introduced by Mirhassani and Fani Jahromi21 that achieved better results in solution quality and computational performance. However, the branching line can receive material from only a single lot. Later, Cafaro and Cerdá22 developed an improved MILP continuous time formulation for the operational scheduling of a multiproduct pipeline system with branching structure. The model permits the branching line to receive material from different lots during the pumping runs. Castro23 developed a continuous time formulation for the scheduling of complex treelike unidirectional pipeline systems with multiple sources and depots. The pipeline is divided into segments that connect two consecutive points (refinery, depot, split point) in the pipeline network. However, the model assumes that different lots of same product are not permitted in the same pipeline segment at the same time. Such an assumption may prevent feasibility in some real world problems. The objective was to minimize the total amount of products pumped into the pipeline. In another study, Herrán et al.24 developed a mathematical model for the scheduling of a multiproduct, multipipeline network over a short-term planning horizon featuring a length of 100 h. The set of pipelines is connected together with a branching structure. The approach divides each pipeline into packages of identical volume each containing a single product. They further divided the time horizon into the intervals of equal length. Their objective was to minimize the pumping and start/stop costs, interface losses, and inventory carrying charges. The authors presented two variants of the model called the simplified and complete models. The simplified model can only be utilized for a high-demand scenario and low pumping cost. Otherwise, the complete version should be implemented. When the latter model is used, the number of variables and constraints are dramatically increased and, consequently, the CPU solution time will be high. Due to this reason, the same authors25 developed several metaheuristic methods such as multistart search, variable neighborhood search, taboo search, and simulated annealing for the same problem24 to overcome its complexity. Recently, Cafaro and Cerdá26 presented an MILP continuous time formulation for the short-term scheduling of a multiproduct pipeline with multiple sources. A case study involving a pipeline system that transports three products from two sources (refinery) to three distribution centers on two alternative pipeline operation modes was explored. They referred to these modes as fungible and segregate modes. In both modes, at any pumping operation only one source can inject material into the pipeline. Afterward, they27 extended the previous work26 for simultaneous injections at two or more input terminals. An intermediate input node can inject products to the pipeline and receive products from it at the same time. The results show that the execution of simultaneous batch injections leads to a substantial reduction on the overall time required to satisfy

terminals to refineries and refined products from refineries to distribution centers. Rejowski and Pinto10 presented a discrete MILP model to schedule a real-world multiproduct pipeline transporting four oil products from a unique refinery to five distribution centers. A volume of the pipeline located between two adjacent nodes is considered as a segment and each segment is divided into packages. They devised two models based on disjunctive programming and discrete time representation. In the first model, the volume of all packages was equal, while in the second one they were different. Batch sequencing is optimally selected so that an objective function accounting for inventory, pumping, and interface costs is minimized. The solution time of the model for a real world case study was extremely high because of too many binary variables. Later they proposed a new model for the same problem by adding special constraints and valid integer cuts to the earlier model to improve the model’s performance.11 Afterward, Rejowski and Pinto12 introduced a hybrid optimization approach based on continuous-time MINLP formulation for the same problem considered in their prior works. The novel MINLP formulation yields a better solution compared to the discrete time representation for the same time horizon. Cafaro and Cerdá1 studied a problem similar to the one considered by Rejowski and Pinto10 and presented a highly efficient MILP continuous time formulation that optimally chooses the sequence of products inside the pipeline, length of pumping operations, pumping rate, volumes, and types of products transferred from the pipeline up to depots. The model handles the problem with the smaller number of binary variables and constraints, with respect to that of Rejowski and Pinto.11 As a result, an optimal solution is obtained about 3-fold faster than that of Rejowski and Pinto.10 In a subsequent work, they generalized the previous model to consider due dates on demands over a monthly moving horizon.13 In their model, the time horizon is partitioned into several time periods of equal or unequal sizes, and product demands for each period were deterministic data. At each time period, new data of the problem are employed to provide a rescheduling procedure. Relvas et al.14 developed a continuous-time MILP model for the scheduling of the multiproduct pipeline from a single refinery to a unique distribution center based on the model presented by Cafaro and Cerdá.1 In their model, pumping rate is a fixed value, and the time horizon comprises a length of 720 h. The model is concentrated upon inventory management at the distribution center, especially settling a period effect. Cafaro and Cerdá15 studied the problem similar to that of Relvas et al.,14 presenting smaller MILP formulation and providing better schedules in a very short computational time. Relvas et al.16 extended the previous work to explore the variable flow rate, pipeline stoppage, and a special settling period for any product. The authors also developed a novel methodology for updating scheduling that handles unpredictable events. Due to the choice of the sequence of products and the settling period constraints, the number of binary variables is dramatically increased. As a result, the solution time of the presented model is high, and the optimum cannot be discovered after a CPU time limit of 2 h. To improve the model performance, Relvas et al.17 presented a heuristic approach to determine the product sequence in the pipeline. Magatão et al.18 presented an MILP discrete-time approach for scheduling batch injections in a bidirectional pipeline (reversible pipeline) with a fixed pumping rate for batches 7030

dx.doi.org/10.1021/ie4038032 | Ind. Eng. Chem. Res. 2014, 53, 7029−7042


Article

Figure 1. Pipeline structure for the problem operating system.

Figure 2. Pipeline schedule for the motivating example.

depot requirements. For the problem of Cafaro and Cerdá,27 a smaller continuous time MILP model based on formulation of Cafaro and Cerdá1 was proposed by Mirhassani et al.28 The pipeline is divided into few segments that connect two consecutive input nodes. Unlike the Cafaro model, 27 an intermediate segment can simultaneously receive materials from an adjacent segment and the input source at the segment origin. However, all the different approaches mentioned above have focused on sequence deliveries, i.e., at any time only one distribution center can receive material from the pipeline while injecting a new product to the line from a unique input node. Actually, a real pipeline does not operate in this manner, and pipeline operators often carry out simultaneous deliveries to reduce the number of pumping operations and maintenance costs. More recently, Cafaro et al.29,30 presented a continuoustime MILP model for operational scheduling of a single-source pipeline system that allows the execution of simultaneous deliveries at two or more distribution centers. In this case, the number of flow restarts and stoppages in pipeline segments is significantly reduced. However, their model assumes that, at most, a single product can be diverted to a distribution center during the execution of a pumping run. For this reason, the number of pumping runs increases to a considerable degree. To overcome such a drawback, we introduce a new mathematical model that permits transferring lots of different products to a distribution center during the same pumping run. As a result, the number of pumping operations through the planning horizon is substantially decreased. Moreover, by cutting down the number of stoppages and restarts of the pipeline segments, considerable savings in maintenance costs are achieved. Two real-world examples have been solved, and the results have been compared with previous works to illustrate that the new MILP formulation provides much better pipeline operational schedules. The remainder of the paper is organized as follows: section 2 gives a short explanation about the considered problem. Section 3 describes the model assumptions. The detailed mathematical formulations of the problem are presented in section 4. The numerical results concerning this allocation are reported in section 5. The last section puts forward the conclusions and sums up the paper.

2. PROBLEM STATEMENT This paper focuses on the long-term scheduling of activities in a specific pipeline system. It connects a unique refinery to several distribution centers (removal terminals). Figure 1 depicts the pipeline physical structure overview. The pipeline is composed of few segments that may present a lower diameter. The first segment connects the refinery at the origin to the first removal terminal. Other segments lie between two successive removal terminals. As a result, the number of segments is equal to the number of removal points. The pipeline carries batches of different types of products (e.g., gasoline, kerosene, liquefied petroleum gas (LPG), aviation fuel), which are oil derivatives. The batches are injected uninterruptedly, and there is no physical barrier between adjacent products as they move in the pipeline. The problem’s goal is to timely meet all product demands at distribution centers by performing simultaneous deliveries and preventing unnecessary stoppages and restarts of pipeline segments, but taking into account the flow rate limitations on pipeline segments. Remark 1: The existent batches inside the pipeline at the starting time of the planning horizon are regarded as old batches throughout the paper (Iold = {i1, i2, ..., i|Iold|}). The first element of Iold occupies the end of the pipeline. Batch i2 is flowed immediately after batch i1 inside the pipeline and so on. 2.1. Motivating Example. To exemplify the problem under study, let us consider a motivating instance as shown in Figure 2. At the starting time of the planning horizon, there are five old batches i5, i4,i3, i2, and i1 inside the pipeline containing products P1, P3, P1, P2, and P4, respectively, as shown at the top of Figure 2. The admissible flow rate ranges at pipeline segments, given in (m3/h) are: [540−800] for Segment 1, [400−720] for Segment 2, and [280−540] for the last one. The pipeline operators have planned to supply products P15000 and P315000 to the removal terminal j1 for the next 30 hours, product P25000 to the removal terminal j2 and product P410000 to the removal terminal j3 for the next two days (with the subscripts pointing out the product sizes). Therefore, there are two delivery due dates taking place at the end of a time period (EPt1 = 30; EPt2 = 48). In fact, the planning horizon (Tmax = 48 h) is divided into two periods of unequal length. 7031



Article

To reach the best pipeline schedule, the operator injects further amounts of product P120000 in batch i5 from time 00.00 h to time 25.00 h (see the second line of Figure 2). During the pumping of batch i5 at a flow rate of 800 m3/h, two delivery operations at removal terminal j1 are scheduled. Tanks at j1 first receive 5000 m3 of product P1 from batch i3, and then 15000 m3 of product P3 from i4 . In the second line of Figure 2, it can be observed that, while executing the stripping operation at removal terminal j1, there is no flow in the pipeline segment connecting removal terminals j1 and j3. Through the first pumping run (see the second line of Figure 2) product demands at the first removal terminal j1 are fully satisfied. Therefore, the operator must supply the product to other terminals (j2 and j3). To better use the pipeline transport capacity and reduce the pump operating costs, the operator decides to transfer products to removal terminals j2 and j3 at the same time. To do this, a new batch i6 containing 15000 m3 of product P4 at a flow rate of 720 m3/h from time 25.00 to 45.83 h is inserted. To meet the demands, the pipeline operator diverts 5000 m3 of batch i2 containing product P2 to removal terminal j2 at a flow rate of 240 m3/h and 10000 m3 of batch i1 conveying product P4 to removal terminal j3 at a flow rate of 480 m3/h. As can be seen in Figure 5, product deliveries to removal terminal j2 diminish the streamflow rate to the relevant value of 480 m3/h in the last segment. It can be derived from the third line of Figure 2 that, when batch i2 is diverted into the removal terminal j2 at a flow rate of 240 m3/h, simultaneously, the upper coordinate of the same batch proceeds in the last segment at a flow rate of 480 m3/h.

12. The planning horizon is divided into time periods, and product demands at removal terminals due before the period end are deterministic data. 13. There is a valve at any removal point as well as at the inlet of any segment of the pipeline. The valves of active terminals and segments are kept open from the beginning of the pumping operation up to its end.

4. OPTIMIZATION MODEL The problem considered in this work is based on the work by Cafaro et al.30 (cited hereafter as the CC model) and the same nomenclature will be applied in some parts as listed at the end of the paper. Although the current paper and the CC model address the same problem with the same objective, there are major differences between them on modeling issues: (i) In the CC model, an active off-take terminal can extract a product from an in-transit batch through a pumping operation only if the batch has previously reached the location of the removal point of that terminal, before starting the pumping run. Due to this requirement, each active terminal only receives product from a single batch during a pumping operation. This limitation has been relaxed, and it is possible that each active terminal can be fed by multiple batches during the same pumping run. To this end, the additional constraints are devised in subsections 4.14 and 4.15. In these sections, the innovative constraints of the proposed model will be differentiated, in comparison to the CC model. Such constraints also enable to reduce the number of on/off pumping switchings. Pipeline operators mostly aim to reduce the number of interfaces in the pipeline. The size of the interface becomes extensive when the flow in the pipeline is started or stopped. To restart and stop the flow in the network, pumps should be turned on and off. (ii) In the CC model, the sequence and the volume of batches discharged to the removal terminals should be determined before solving the model. To satisfy these necessities, the model presented by Cafaro and Cerdá13 must be solved to specify how much of each batch should be discharged to each removal terminal. Such a methodology is feasible, but in many cases it is not an optimal solution. In contrast with the CC model the sequence and the volume of batches discharged to the removal terminals are both model decisions. Such conditions are imposed by the constraints presented in subsections 4.6, 4.7, and 4.8. In particular, the proposed model is capable of considering multiple due dates at period ends (constraints in subsection 4.17). In the CC model, the model only manages the upper and the lower bounds on delivery sizes. (iii) Analogous to case (ii), in the CC model, the sequence and the volume of batches injected into the pipeline at the origin must initially be discovered by solving the aggregated plan provided by Cafaro and Cerdá.13The CC model only controls the upper and the lower bounds on batch injections. While in our formulation, the model will decide which size of the given product to inject in each batch. To this end, the appropriate constraints are put forth in subsections 4.3, 4.4, and 4.5. Since batch sizing constraints and mass balance constraints do not depend on the problem configuration, such constraints are roughly similar to those first presented by the CC model.

3. MODEL ASSUMPTIONS To present a mathematical model based on MILP the following assumptions are considered: 1. All products have constant densities and move along the pipeline without any physical barrier in between. 2. The pipeline is always full, and the only way to unload a volume of product (to removal terminals) is to pump an equal volume into the pipeline at the origin. 3. At any pumping operation, only a single product is injected into the pipeline. 4. The sequence and the volume of old batches are known. 5. The inventory level of products in the refinery to be injected in the pipeline through the planning horizon is given. 6. The upper and the lower bounds on the size of batch injections are known. 7. The upper and the lower bounds on the volume of product deliveries to distribution centers are also given. 8. When two batches with different products meet each other, some of them are mixed together. This part is regarded as an interface. The interface volume between each particular pair of products is a known constant. 9. Due to the vast product contamination (excessive interface), some specified products cannot be pumped back-to-back into the pipeline (forbidden product sequences). 10. The flow rate may change with the pipeline segment within the allowable range. 11. The flow rate from the pipeline to the removal terminal can vary in the permitted range. 7032



Article

4.1. Constraints for Length and Start Time of Each Pumping Run. Let STk be the starting time of the pumping run k and Lk represent its duration. The starting of the pumping run k should not be before the end of the previous one (k − 1). Hence, STk = Lk − 1 + STk − 1 ′ ∀ k ∈ K (k ≥ 2)

(1)

Besides, the duration of any pumping operation should not exceed the overall length of the time horizon (Tmax). Thus,

∑ Lk ≤ Tmax

(2)

k∈K

4.2. Product Allocation to Batches. Let the binary variable yi,p show that product p is contained in batch i whenever yi,p = 1. Otherwise yi,p = 0, and no product is allocated to batch i. Every batch can, at most, carry one product, so

∑ yi ,p ≤ 1,

∀i∈I Figure 3. A simple example illustrating constraints in parts 4.1−4.8

(3)

p∈P

Remark 2: The set Inew given by Inew = {i|Iold|+1, i|Iold|+2,..., i|Inew|}, the set of new batches to be injected into the pipeline during the time horizon. The cardinality of set Inew is arbitrarily adopted. Always, batch in is injected right after batch in−1 into the pipeline. If the number of new batches to be injected in the pipeline is lower than the cardinality of set Inew by some elements in, some in predefined new batches are never inserted into the line. They are considered as dummy batches. For dummy batches i ∈ Inew newer injected into the line, the binary variable yi,p takes a value of 0 for every p ∈ P. Some model equations may force the dummy batches to be the last elements of the set Inew. To reduce the size of the feasible region of the problem and accelerate the search for the optimal schedule, dummy batches i ∈ Inew featuring yi,p = 0 for all p ∈ P, must be placed at the end of batch sequences. Therefore,

∑ yi ,p ≤ ∑ yi− 1,p , p∈P

∑ i ∈ I , i ≥|I old|

∑

i ∈ I , i ≥|I old|

i ∈ I , i ≥|I old|

zik ≤ |K | ×

k∈K

(6)

∑ yi ,p ∀

i ∈ I new , k ∈ K (7)

p∈P

Abki

Let denote the amount of batch i injected at the origin of the pipeline by the pumping run k. Abki for any nonempty batch i is limited to its boundaries. Thus, zikAbmin ≤ Abik ≤ zikAbmax , ∀ i ∈ I (i ≥ |I old|), k ∈ K (8)

where (Abmin, Abmax) represent the minimum and maximum batch sizes that potentially can be injected during a pumping run. 4.4. Relation between Volume and Pumping Duration. Using a nonconstant injection rate within the interval [ vbmin; vbmax], the relationship between the volume of lot (Abki ) and the length of the pumping k is given by:

(4)

Remark 3: A batch can be injected into the pipeline in consecutive pumping runs (see Figure 3). Such a condition is formulated in the sequel by constraint (5) and further by eq (25). 4.3. Size of Batch Injected into the Pipeline. Let us define the binary variable zki to represent the pumping run k injects a batch i ∈ I in the pipeline whenever it equals 1. By assumption 3, at most a single batch can be inserted into the pipeline through any pumping run, so

∑

zik − 1 , k ∈ K (k ≥ 2)

∑

And by Remark 4:

∀ i ∈ I new

p∈P

zik ≤

Lk vbmin ≤

Abik ≤ Lk vbmax , ∀ k ∈ K

∑ ∀ i ∈ I(i ≥|I old|)

(9)

4.5. Inventory Control in Refinery. To manage the product inventories at the refinery, the model must identify the volume and the type of any product injecting to the line through any pumping run k. To do that, let APki,p be the volume of product p contained in batch i inserted from the origin of the pipeline during the pumping run k. At most, one product can be injected from refinery the pipeline by a pumping run k, so

zik ≤ 1, ∀ k ∈ K (5)

∑ APik,p = Abik ,

Remark 4: Empty batches are never injected into the line. Since the number of pumping operations is highly dependent on the overall terminal requirements and the size of injected batches, the cardinality of the set K must be initially surmised. It should be as low as possible to reduce the size of the feasible region. When the elements in K to discover the optimal solution are not enough, its cardinality is increased by one, and the model must be resolved. However, in some cases, the number of pumping operations to be performed is lower than | K|. In this case, some predefined elements of the set K are never carried out and must be confined to the end of injections to reduce the solution degeneracy. Hence,

∀ i ∈ I(i ≥ |I old|), k ∈ K (10)

p∈P

On the other hand, when product p is not allocated to batch i, the value of APki,p should be equal to zero. Thus,

∑ APik,p ≤ |K | × yi ,p × Abmax ,

∀ i ∈ I(i ≥ |I old|), p ∈ P

k∈K

(11)

Let us assume that parameter Invp presents the total amount of product p that can be pumped from the refinery to the line during the planning horizon. So, 7033



∑ ∑ APik,p ≤ Invp,

Article

material injected to batch i through the pumping run k (Abki ), minus the overall volume transferred from batch i to the removal terminals through run k (∑jAdki,j). Hence,

∀p∈P (12)

k ∈K i∈I

4.6. Stripping Material from the In-Transit Batch. Let the binary variable wkj denote that removal terminal j is active during the pumping run k, whenever wkj = 1, otherwise wkj = 0. There is at least one batch that is withdrawn by removal terminal j during the pumping run k, when wkj = 1. As a result, wkj = 1 if and only if ∑i∈Ixki,j ≥ 1 Another way of stating this requirement is, wjk ≤

∑ xik,j ≤ |I |wjk ,

W ik = W ik − 1 + Abik −

∑ Adik,j ,

∀ i ∈ I, k ∈ K

j∈J

Initial volumes of old batches that existed in the pipeline at the start of the scheduling horizon are determined by the parameter IWi. W ik = IW,i ∀ i ∈ I old , k = 0

∀ j ∈ J, k ∈ K (13)

i∈I

(20)

(21)

Figure 4 is used to illustrate constraint (20).

xki,j

where is a binary variable denoting that terminal j receives material from batch i during the pumping run k, whenever xki,j = 1. Now, let Adki,j be the size of batch i transferred to removal terminal j by the pumping run k. Adki,j is positive if and only if xki,j = 1 and it is fixed to zero whenever xki,j = 0. Then, Adminxik, j ≤ Adik, j ≤ Admaxxik, j , ∀ i ∈ I , j ∈ J , k ∈ K

(14)

where Admin and Admax respectively stand for a lower and upper bound on the amount of material that can be diverted from batch i to the removal terminal j. 4.7. Activity Duration of Removal Terminals. Let LDkj be the activity length of removal terminal j through the pumping run k. Hence, LDkj × vdjmin ≤

∑ Adik,j ≤ LDjk × vdjmax ,

Figure 4. Size of each batch inside the pipeline at time Sk + Lk.

4.10. Batch Location in the Pipeline. When a pumping run injects a new batch into the pipeline, the location of some batches inside the pipeline may change. Therefore, the model must be capable to update the coordinates of any batch at the end of any new injection. Let Fki be the upper coordinate of batch i at the end of pumping run k i.e., at time Sk + Lk, Such a continuous variable denotes the pipeline volume between the origin and the farthest extreme of batch i at time Sk + Lk. According to the definition, the upper coordinate Fki is equal to the overall volume of batches pursuing batch i at time Sk + Lk, plus the size of batch i at the completion time of pumping run k (see Figure 5).

∀ j ∈ J, k ∈ K

i∈I

(15) max [vdmin j ;vdj

The interval ] stands for the admissible feed rate from the pipeline into the removal terminal j. If no material is transferred into the removal terminal j during the pumping run k, then wkj = 0, and from constraint (16), LDkj = 0. Further, LDkj max must be kept within the allowed range [LDmin j ; LDj ] only if terminal j is activated by run k. k k max k LDmin j wj ≤ LD j ≤ LD j wj , ∀ j ∈ J , k ∈ K

Fik =

(16)

W ik′ , ∀ i ∈ I , k ∈ K

∑ i ′≥ i , i ′∈ I

(22)

4.8. Volume of Product Discharged from Batches to Removal Terminals. Let VPki,p,j represent the volume of product p transferred from batch i to terminal j through the pumping run k. If yi,p = 0 the value of VPki,p,j will be equal to zero. Otherwise yi,p = 1, and the volume of all products that are diverted from batch i to removal terminal j is equal to Adki,j.

∑ VPik,p,j ≤

|K | × yi , p × Admax , ∀ i ∈ I , p ∈ P , j ∈ J

k∈K

∑ VPik,p,j = Adik,j ,

(17)

Figure 5. Upper coordinate of batch i2 inside the pipeline.

(18)

As can be observed in Figure 5, Fki represents both the upper coordinate of batch i and the lower coordinate of batch i − 1 at the end of pumping run k at the same time. 4.11. Mass Balance Constraints. A crucial aspect of the multiproduct pipelines is that they must be completely full with product at any time. So at the end of any pumping operation k, the size of all batches moving along the pipeline must be equal to the overall volume of the pipeline from the origin to the farthest removal terminal (PV). Hence,

∀ i ∈ I, j ∈ J, k ∈ K

p∈P

Besides, the overall volume of product p that can be stored in distribution center j are bounded.

∑ ∑ VPik,p,j ≤ Inp,j , k∈K i∈I

∀ p ∈ P, j ∈ J (19)

Figure 3 illustrates constraints in parts 4.2− 4.8. 4.9. Size of Batch i in the Pipeline at the End of the Pumping Run k. Let Wki denote the size of batch i at the end of the pumping run k (Sk + Lk) The value of Wki is equal to the size of batch i at time Sk (Wk−1 i ), plus the amount of the

∑ Wik = PV, ∀ k ∈ K i∈I

7034

(23)



Article

coordinate of that batch at time Sk, i.e (Fk−1 i ) has not surpassed the volumetric coordinate of the removal terminal j(σj). Now, suppose that the batch i is diverted to the removal terminal j during the pumping run k (xki,j = 1). Then, the volume of the material transferred from batch i to removal k−1 terminal j is limited by σj − Fk−1 (see Figure 7). In i+1 or Wi

According to the assumptions 1 and 2, during any pumping run k, the overall volume discharged from batches to removal terminals (∑i∈I∑j∈J Adki,j) must be equal to the total volume of product injected to the pipeline at the origin (∑i∈I Abki ).

∑ Abik = ∑ ∑ Adik,j , ∀ k ∈ K i∈I

i∈I j∈J

(24)

4.12. Condition for Supply Material to the Batches. Through the pumping run k, a batch i can receive material from the refinery only if the lower coordinate of that batch at time Sk is equal to zero. As a result, zki = 1 yields Fk−1 i+1 = 0. Based on this condition we have constraint (25). Fik+−11 ≤ PV(1 − zik), ∀ i ∈ I(i ≥ |I old|), k ∈ K (k ≥ 1) (25) Figure 7. Simple example illustrating constraint (31).

4.13. Interface Volume between Two Successive Batches. Since there is no physical separation between different products inside the pipeline, there is always a certain volume of intermixing between the two adjacent batches at the interface, i.e. the point where they meet. The interface volume depends on the physical properties of the products and is assumed to be a constant value, MIXp,q (see Figure 6). The

some situations, it is possible that the pumping operation k injects a further amount of product to the existing batch i while some portion of this batch is extracted by removal terminal j. As k a result, xki,j = 1 implies Adki,j ≤ σj − Fk−1 i+1 + Abi . Thus the following constraint is valid. i

∑ Adik′ ,j + Fik+−11 ≤ σj + Abik + PV (1 − xik,j), i ′= 1

Figure 6. Interface material between two products.

∀ i ∈ I, j ∈ J, k ∈ K

4.15. Controlling the Stream Flow Rate at Pipeline Segments. Let the binary variable vkj denote that segment j is active during the pumping run k whenever vkj = 1. Otherwise, vkj = 0, and no material can be transferred from segment j − 1 to segment j. Clearly, segment j + 1 is out of work, when segment j is idle. As a result, vkj = 0 yields vkj+1 = 0, and vkj+1 = 1 implies vkj = 1 for all j ∈ J, k ∈ K Based on these conditions we have constraint (32).

lower bound of continuous variable INVi,p,q, which presents the interface volume between batches i and i + 1 when they transport product p and q respectively, is obtained by the following constraint. INVi , p , q ≥ MIX p , q × (yi , p + yi + 1, q − 1), ∀ i ∈ I , p , q ∈ P (26)

Remark 5: Since the value of ∑i∈I ∑p∈P(p ≠ q) ∑q∈PINVi,p,q will be minimized during the planning horizon, the upper bound of INVi,p,q is restricted by the objective function. Because of high product contamination (vast interface), some product sequences are not allowed. This restriction is described by the following constraint: yi , p + yi + 1, q ≤ 1 + Seq p , q, ∀ i ∈ I , p , q ∈ P

vjk+ 1 ≤ vjk , ∀ j ∈ J , k ∈ K

≥

σjxik, j ,

∀ i ∈ I, j ∈ J, k ∈ K

(27)

v1k ≤

∑ zik , ∀ k ∈ K i∈I

(33)

Recall that in this work a single source unidirectional pipeline system is considered, thus whenever the removal terminal j is receiving material from the pipeline, all of the pipeline segments between the refinery and the receiving terminal j are active. Hence, wjk ≤ vjk , ∀ j ∈ J , k ∈ K

(34)

LSkj

Now, let be the activity duration of segment j through the pumping run k. The value of LSkj should be between the minimum and maximum bounds, as indicated in constraint (35).

(28)

Fik+−11 ≤ σj + PV (1 − xik, j), ∀ i ∈ I , j ∈ J , k ∈ K (k ≥ 1) (29)

vjk LDmin ≤ LSkj ≤ vjk LDmax , ∀ j ∈ J , k ∈ K

Fik − 1 ≥ σjxik, j − PV (1 − wjk′), ∀ i ∈ I , j , j′ ∈ J(j′ > j), k ∈ K (k ≥ 1)

(32)

In addition, when the first segment is operating through the pumping run k (vkj = 1), there is a batch i ∈ I that is inserted into the first segment. This condition is mathematically written as follows,

If product p cannot be adjacent to product q, the value of parameter Seqp,q is zero; otherwise Seqp,q = 1 4.14. Condition for Diverting Material from the Pipeline into Removal Terminals. During the pumping run k, a batch i can be withdrawn by the removal terminal j only if the lower and upper coordinates of that batch at two consecutive pumping runs k and k − 1 satisfy the following constraints. Fik

(31)

(35)

Remark 6: As stated before, there is a valve at the inlet of any segment of the pipeline. Inlet valves of active segments are kept open from the beginning of the pumping operation up to its end.

(30)

These constraints state that a batch i can be transferred to the removal terminal j during the pumping run k only if, the lower 7035



Article

By Remark 6, it can be drawn that through any pumping operation, the activity duration of all active segments between refinery and the farthest active removal terminal is of the same size. Such a condition is considered by the following constraints. Tmax(2 − vjk − vjk′) + LSjk ≥ LSjk′ , ∀ j , j′ ∈ J , k ∈ K (36)

Figure 8. Simple example illustrating constraints (39)−(41).

In addition by Remark 6, the activity duration of an active segment j through a pumping run k must be equal to the length of the same pumping operation. Hence, Tmax(1 − vjk) + LSkj ≥ Lk ≥ LSkj , ∀ j , k ∈ K

farthest active terminal, then AP1 = σj2 = 60 m3. Thus, constraints (39)−(41) must result in AP1 = 60 m3. Constraint (39) reads for j = j1 and k = 1, AP1 ≥ σj1 × w1j1 = 25 × 1 = 25 m3. Similarly for j = j2 and k = 1, AP1 ≥ σj2 × w1j2 = 60× 1 = 60 m3; and for j = j3 and k = 1, AP1 ≥ σj3 × w1j3 = 90 × 0 = 0. Thus AP1 ≥ 60 m3. On the other hand, by constraint (40) and for j = j1 and k = 1, AP1 ≤ PV + (PV − σj1)(v1j2 − v1j1) = PV + (PV − 25) × (1 − 1) ≤ PV = 90 m3. Similarly, for j = j2 and k = 1, AP1 ≤ PV + (PV − σj2)(v1j3 − v1j2) = PV + (PV − 60) × (0 − 1) ≤ 60 m3, and consequently, AP1 ≤ 60 m3. Furthermore, through the constraint (41), AP1 ≤ PV = 90 m3. Therefore, for the first pumping run, we have AP1 = 60. So constraint (42) suggests the pipeline volume activated by the first pumping run, i.e. (AV1) being 60 m3. (AV1 = AP1 − AP0 = 60 − 0 = 60). In the second pumping run (the third line of Figure 8), the removal terminal j3 is the farthest active terminal receiving material from the pipeline, and consequently AP2 = σj3 = 90 m3. Through the constraint (39), AP2 ≥ σj1 × w2j1 = 25 × 0 = 0, AP2 ≥ σj2 × w2j2 = 60× 1 = 60 m3, and AP2 ≥ σj3 × w2j3 = 90× 1 = m3. Thus AP2 ≥ 90 m3. On the other hand, by constraint (40), AP2 ≤ PV + (PV − σj1)(v2j2 − v2j1) = PV + (PV − 25)(1 − 1) = PV = 90 m3. Similarly AP2 ≤ PV + (PV − σj2)(v2j3 − v2j2) = PV + (PV − 60)(1 − 1) = PV = 90 m3. Besides, through the constraint (41), AP2 ≤ 90 Therefore, AP2 = 90 as it should. Thus, the pipeline volume activated by the second pumping run, i.e. (AV2), takes a value of 30 m3. (AV2 = AP2 − AP1 = 90 − 60 = 30 m3). 4.17. Delivery Due Date Constraints. Suppose that the planning horizon Tmax has been partitioned into several time periods t with equal or unequal sizes and EPt presents the utmost of period t. Moreover, suppose that the demand for each removal terminal due at period t is given. The aim is to satisfy all product demands due at period t. To manage this requirement, let us define the binary variable mkt denoting that the pumping run k is completed within the period t whenever mkt = 1. Every pumping operation must be completed at one of the periods, and for a dummy run k, mkt must be zero. Hence,

(37)

Through any pumping run k, the streamflow rate at the operational segment j should be set within the feasible range [vsjmin;vsmax j ]. The streamflow rate limitation for the first segment is controlled by constraint (9). For other segments, the following constraint is imposed. LSkj vsjmin ≤

∑ ∑ Adik,j ′ ≤ LSjk vsmax j , j ′∈ J , j ′≥ j i ∈ I

∀ j ∈ J(j ≥ 2), k ∈ K

(38)

Remark 7: To reduce the size of the problem feasible region, constraints (30) can be replaced by the following equation: Fik − 1 ≥ σjxik, j − PV(1 − vjk+ 1), ∀ i ∈ I , j ∈ J , k ∈ K (k ≥ 1)

4.16. Computing the Active Volume of Pipeline for Run k. Let ACPk be the active volume of pipeline through the pumping run k. This continuous variable presents the overall volume between the pipeline system origin (refinery) and the farthest active terminal during the pumping run k. In other words, if depot j is the farthest active terminal through the pumping run k, the active volume of the pipeline is determined by ACPk = σj. Equivalently, ACPk ≥ σjwjk , ∀ j ∈ J , k ∈ K

(39)

ACPk ≤ PV + (PV − σj)(vjk+ 1 − vjk), ∀ j ∈ J(j < |J |), k ∈ K

(40)

On the other hand, the active volume of the pipeline at any pumping run must never be greater than PV. Thus, ACPk ≤ PV, ∀ k ∈ K

(41)

Due to the operational strategies, the pipeline segments must operate at the minimum flow restarts and stoppages. To this end, the model must minimize the total activated and stopped pipeline volumes during a limited time horizon (Tmax). Let AVk and SVk be the activated and stopped volume of pipeline through the pumping run k, respectively. Hence, AVk ≥ ACPk − ACPk − 1 , ∀ k ∈ K (k ≥ 1),

(42)

SVk ≥ ACPk − 1 − ACPk , ∀ k ∈ K (k ≥ 1)

(43)

∑ mtk = ∑ zik , t∈T

i∈I

∀k∈K (44)

On the other hand, when the pumping run k is completed in the period t, the starting time of the pumping run k + 1 must take a value in range [EPt−1; EPt]. As a result, mkt = 1 yields EPt−1 ≤ Sk+1 ≤ EPt. Another way of stating this requirement is: EPt + Tmax(1 − mtk ) ≥ Sk + 1 ≥ E Pt − 1mtk , ∀ t ∈ T , k ∈ K (45)

demtp,j

Let represent the overall demand of product p to be satisfied at the removal terminal j through the time period t. Product demands at removal terminals should be satisfied on time, and unsatisfied demand undemtp,j results in penalty cost. Hence,

Since we use the variables AVk and SVk in the objective function, their upper bounds are controlled by the model. To justify the validity of constraints (39)−(41), let us consider an expository example as shown in Figure 8. At the first pumping operation (Run 1), the removal terminal j2 is the 7036



Article

t ′∈ T , t ′= 1

(

∑

4.19. Problem Objective Function. The problem objective function, given in eq (53) below, comprises four major terms. The first one is the on/off pump switching costs. The second one stands for the cost of activated/stopped pipeline volumes during the planning horizon. The third term represents the cost of unsatisfied demands at removal terminals. The last one points out the reprocessing cost of the interface volume between two adjacent products.

dempt ′, j)(mtk − mtk + 1) + undempt −, j 1

t k ′∈ K , k ′= 1

≤

∑ ∑ VPik, p′ , j + undempt , j k

(46)

i∈I

where p ∈ P, j ∈ J, t ∈ T, k ∈ K. When = 1 and = 0, the value of ∑kk′=1 ∑i VPki,p,j ′ stands for the overall amount of product p diverted from the pipeline into the removal terminal j within the interval [0,EPt]. Demand for period t must be satisfied during the time interval [EPt−1; EPt], in the light of the fact that unsatisfied demands at the end of each period are usually highly penalized. Therefore, to meet the demand due at period t, at least one pumping operation should be completed in this period. mkt

∑ mtk ≥ 1, ∀

mk+1 t

Here, the coefficient f co represents the fixed cost to fulfill a pumping operation. The parameters ca and cs stand for the unit flow restart and stoppage costs, respectively. The coefficient cutp,j is the unit backorder penalty cost to tardy meet demand due at period t. In the end, the parameter crp,q is the unit reprocessing cost of the interface material including products p and q.

t∈T

k∈K

(47)

4.18. Simultaneous Delivery Constraints. In this work, an active removal terminal must receive material from the pipeline while the refinery starts to inject material into the pipeline. Thus, batch i should be withdrawn by an active removal terminal j providing the following condition: Fk−1 i+1 ≤ σj k−1 ≤ Fk−1 It follows that, if a batch i, which satisfies F i i+1 ≤ σj, is not extracted by the removal terminal j, any batch succeeding batch i cannot be transferred to the same removal terminal. As a k k result, Fk−1 i+1 ≤ σj and xi,j = 0 implies xi′,j = 0 for all i′ ∈ I(i′ > i) . Such a condition is described by the following constraint.

5. RESULTS AND DISCUSSION To validate the proposed mathematical model, two real-life case studies are explored. In both examples, the optimal solution was found at a rational CPU time. The MILP mathematical formulation was solved to optimality on an Intel Core 2 Duo CPU 2.10 GHz Processor, 4 GB of RAM using CPLEX 12.2 with 6 parallel threads as an MILP solver. 5.1. Example 1. This example is similar to the one presented by Cafaro et al.,30 in which a unidirectional pipeline transports four refined products (gasoline, P1; diesel, P2; LPG, P3; and jet-fuel, P4) from an oil refinery to five distribution centers (D1−D5) over a time horizon comprising four 1-week periods. The pipeline length is 955 km, and the pumping rate at the refinery should be kept between 700 and 1200 m3/h. The pipeline network structure for Example 1 is depicted in Figure 9. Similarly to Cafaro et al.,30 the lower bound on the size of

σj − PV(1 − xik+ 1, j + xik, j) ≤ Fik+−11, ∀ i ∈ I , j ∈ J , k ∈ K (48)

Further, suppose that a batch i ∈ I is transferred to the removal terminal j during the pumping run k . If batch i + 1 which directly follows batch i inside the pipeline is not transferred to the same removal terminal, the upper coordinate of this batch should not surpass the volumetric coordinate of the removal terminal j. As a result, xki,j = 1 and xki+1,j = 0 yields Fki+1,j ≤ σj for all i ∈ I, j ∈ J, k ∈ K. On the basis of this condition the following constraint is considered. Fik+ 1 ≤ σj + PV(1 + xik+ 1 − xik), ∀ i ∈ I , j ∈ J , k ∈ K (49)

Figure 9. The pipeline distribution network considered in Examples 1 and 2.

Recall that there is a valve at any removal point on the pipeline and the valves of active terminals are kept open from the beginning to the end of the pumping operation. Such an operating feature can be mathematically written as follows,

product delivery to distribution center is assumed to be 250 m3 to compare model outcomes. Moreover, we use the values f co = 1000 ($/run), ca = 0.10 ($/m3) and cs = 0.00 ($/m3) for the parameters employed in the objective function. The pipeline system is composed of five segments, and the volumes of the segments (refinery−D1, D1−D2, D2−D3, D3− D4, and D4−D5) are 400, 250, 250, 600 and 135 [ × 102 m3], respectively. Moreover, the admissible flow rate ranges at the pipeline segments are as follows: 700−1200 m3/h for the first segment (refinery−D1) and 0−1200 m3/h for other segments. The volume of products available in the refinery and the interface material cost and volume for two different products are all listed in Table 1. Forbidden product sequences are indicated with an “X” in this table. Moreover, product demands at removal terminals D1−D5 to be satisfied at the end of periods t1−t4 are given in Table 2.

Tmax(2 − wjk − wjk′) + LDkj ≥ LDkj ′, ∀ j , j′ ∈ J , k ∈ K (50)

Tmax(1 − wjk) + LDkj ≥ Lk , ∀ j ∈ J , k ∈ K

(51)

Lk ≥ LDkj , ∀ j ∈ J , k ∈ K

(52)

Constraint (50) states that at any pumping operation, the activity duration of all removal terminals receiving material from the in-transit batches should be of the same length. Moreover, the activity length of each active removal terminal is equal to the length of pumping operation performed at refinery. Such a feature is added to the problem formulation through constraints (51) and (52). 7037



Article

Tables 3 and 4 indicate the performance indicators of Cafaro and Cerda’s model13). 5.2. Example 2 (variant of Example 1). The case study considered in Example 1 is taken again in Example 2, but this time, the admissible flow-rate ranges at the pipeline segments are modified to obtain more realistic operating conditions. The flow rate ranges in pipeline segments are as follows: 700−1200 m3/h for the first segment (refinery−D1), 600−1200 m3/h for segments D1−D2, D2−D3, and D3−D4, and 400−800 m3/h for the last segment, where the pipeline presents a lower diameter. The best pipeline schedule for this example is shown in Figure 11. A total of 26 pumping runs are carried out over the planning horizon, with 22 of them transferring products to more than one distribution center simultaneously. Assuming that the pipeline is initially idle, there are the following: (1) only one flow restart in segment refinery-D1 at operation k1, (2) just one flow resume in segment D1−D2 when executing the same operation k1, (3) three flow restarts in segment D2− D3 at pumping run k2, k10, and k15, (4) four flow restarts in segment D3−D4 at operations k2, k10, k16, and k18, and (5) eight flow resumptions in the last segment at pumping runs k3, k10, k12, k16, k19, k21, k23, and k25. Such flow restarts involve a total activated volume given by (1 × 40000 + 1 × 25000 + 3 × 25000 + 4 × 60000 + 8 × 13500) = 488000 m3. Table S2 (in SI) shows the flow set point for each valve at every pumping operation. From this table, it can be observed that the streamflow rate in each active segment is within the admissible range. Table 4 shows computational results for Example 2. The results in Table 4 reveal that: 1. The pipeline schedule obtained by the CC model includes a higher number of pumping runs; i.e. a total of 40 batch injections against 26 performed at our proposed model. 2. Through the careful coordination among delivering activities at output terminals provided by the proposed MILP model, the activated volume of pipeline is reduced by 16.5%. 3. The stopped volume of pipeline is decreased by 17.18%. 4. The total operating cost of $74,800 is almost 25% less than the one by Cafaro et al.,30 amounting to $98,450.

Table 1. Product Supplies and Interface Costs and Volumes inventory (× 102 m3) P1 P2 P3 P4

interface cost (× 102 US$)/ volume (× 102 m3)

2000 3320 390 1100

P1 − 30/0.30 37/0.37 35/0.35

P1 P2 P3 P4

P2 30/0.30 − X 38/0.38

P3 37/0.37 X − X

P4 35/0.35 38/0.38 X −

The initial condition of the pipeline network is shown at the top of Figure 10. There is a sequence of five old batches i5, i4, i3, i2, and i1 inside the pipeline transporting products P2, P1, P3,P1 and P2, respectively, and featuring the following volumes (400−700−200−200−135), in a hundred cubic meters. The best pipeline schedule found in 5718 s of CPU is depicted in Figure 10. It includes a total of 23 pumping runs injecting a sequence of nine new batches. All of the pumping runs transferred products to more than one removal terminal simultaneously. There is only one temporal stoppage volume in the last pipeline segment during the time interval 23.75−40.42 h. Assuming that the pipeline is initially idle, there are six flow restarts in pipeline segments during the planning horizon distributed as follows: Refinery−D1 (1), D1−D2 (1), D2−D3 (1), D3−D4 (1) and D4−D5 (2). Such flow restarts involve a total activated volume given by (1 × 40000 + 1 × 25000 + 1 × 25000 + 1 × 60000 + 2 × 13500) = 177000 m3. Flow rates at the control valves during each pumping operation are presented in Table S1 (in Supporting Information [SI]). From this table, it can be observed that in many times (15 runs of 23); the pipeline operates at full capacity i.e., 1200 m3/h. In particular, during the pumping run 2 (see Table S1 in SI); the flow rate through the segments refinery−D1, D1−D2 and D2−D3 is 1200 m3/h. It is calculated by dividing the total volume injected (20000 m3) by the pumping run length (16.67 h). For segment D3−D4, the flow rate is reduced to 840 m3/h, because the batch B4 is diverted to terminal D3 at a flow rate of 360 m3/h. Table 3 shows computational results for Example 1. The results in Table 3 reveal that: 1. The pipeline schedule obtained by our proposed model includes a lower number of pumping runs, i.e. a total of 23 batch injections against 33 reported by Cafaro et al.30 2. At the optimum, the activated volume of the pipeline is reduced by 7.08% with respect to the CC model. 3. The stopped volume of the pipeline is decreased by 86.56%. 4. The value of z′ is 21.80% less expensive compared to that of the CC model. 5. Since a two-level planning approach was applied by Cafaro et al.,30 the model dimensions and CPU times are not comparable (the numbers in the parentheses in

6. CONCLUSIONS In this study, an MILP continuous-time formulation for the scheduling of a multiproduct pipeline connecting a unique refinery to several removal terminals has been devised. The proposed formulation allows to consider multiple due dates at period ends. The problem goal is to timely meet all product demands at distribution centers by performing simultaneous deliveries and preventing unnecessary stoppages and restarts of pipeline segments, but taking into account the flow rate

Table 2. Depot Demands at Distribution Centers product demands (× 102 m3) D1

P1 P2 P3 P4

D2

D3

D4

D5

T1

T2

T3

T4

T1

T2

T3

T4

T1

T2

T3

T4

T1

T2

T3

T4

T1

T2

T3

T4

40 100 30 0

30 120 40 0

50 100 30 0

50 120 20 0

100 100 0 0

100 100 0 0

150 100 0 0

120 110 0 0

90 70 20 0

120 80 30 0

100 70 20 0

110 60 30 0

140 200 50 60

180 200 60 80

170 200 50 60

150 220 40 70

100 220 30 70

120 210 20 80

90 250 20 60

100 220 40 90

7038



Article

Figure 10. Optimal pipeline schedule for Example 1.

Table 3. Computational Results for Example 1 no. of pumping runs |K| no. of constraints no. of binary variables no. of continuous variables CPU time (s) interface cost (US$) activated volume (m3) stopped volume (m3) z′ (US$) optimal cost (US$)

Table 4. Computational Results for Example 2

the CC model

our work

33 5651 + (15247) 478 + (1089) 13194 + (10259) 6.04 + (688) 1961500 190500 100500 52050 2013550

23 14848 2207 10480 5718 1961500 177000 13500 40700 2002200

no. of pumping runs |K| no. of constraints no. of binary variables no. of continuous variables CPU time (s) interface cost (US$) activated volume (m3) stopped volume (m3) z′ (US$) optimal cost (US$)

limitations on pipeline segments. To validate the proposed approach, a pair of real-world case studies was solved. Both involve the scheduling of a single unidirectional pipeline transporting four petroleum products from a unique refinery to five offtake stations. Experimental results and comparisons demonstrated the effectiveness of the proposed model in handling large scale systems for long-term horizons at quite acceptable CPU times. Summarizing the results, two interesting

the CC model

our work

40 6735 + (15247) 567 + (1089) 15312 + (10259) 124.9 + (688) 1961500 584500 494500 98450 2059950

26 17454 2489 11742 8899.2 1961500 488000 409500 74800 2036300

conclusions should be mentioned. First, the number of pumping operations (|K|) to discover the optimum significantly decreases and, consequently, substantial savings in pump operating costs are obtained. Second, through the careful coordination among delivering activities at output terminals provided by the proposed model, the number of flow restarts in pipeline segments is considerably deceased. On a future work, 7039



Article

Figure 11. Optimal schedule for Example 2.

■

the proposed model will be generalized to handle multiplesource pipeline networks.

■

Sets

k, k′ ∈ K = set of pumping operations Inew = new batches to be injected into the pipeline during the planning horizon Iold = batches inside pipeline before i, scheduling i, i′ ∈ I = all product batches (I = Iold ∪Inew) j, j′ ∈ J = set of removal terminals and segments p, q, ∈ P = oil derivatives t, t′, ∈ T = set of time periods

ASSOCIATED CONTENT

S Supporting Information *

Tables S1 and S2. This material is available free of charge via the Internet at http://pubs.acs.org.

■

NOMENCLATURE

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]

Parameters

Notes

Tmax = time horizon (hour) EPt = upper extreme of time period t LDmin /LDmax = minimum/maximum activity length of j j removal terminal j during any pumping run vbmin/vbmax = minimum and maximum pumping rate (units/ h)

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors are grateful to anonymous reviewers for their most constructive and insightful comments and suggestions. 7040



Article

vsjmin/vsjmax = minimum/maximum product flow rate (units/ h) in segment j vdjmin/vdjmax = minimum/maximum flow rate (units/h) from pipeline to removal terminal j σj = volumetric coordinate of the removal terminal j from the refinery PV = volume of pipeline from origin to the farthest removal terminal Abmin/Abmax = minimum/maximum batch injection size Admin/Admax = minimum/maximum delivery size to removal terminals Seqp,q = Boolean matrix of possible sequences between subsequent products inside the pipeline, p and q Invp = maximum amount of product p that potentially can be injected to the pipeline during the planning horizon Inp,j = maximum amount of product p that potentially can be transferred to removal terminal j during planning horizon IWi = initial volume of old batch i demp,jt = overall demand of product i to be satisfied by removal terminal j at the time period t f co = fixed cost for fulfilling a pumping operation ca = unit flow restart cost cs = unit flow stoppage cost crp,q = cost of reprocessing a unit interface volume between two products p and q cutp,j = unit backorder penalty cost to tardy meet demand due at period t

■

xki.j = One if and only if a portion of batch i is transferred to the removal terminal j during the pumping run k; otherwise 0 mkt = One if and only if the pumping run k terminates within time period t; otherwise 0 yi,p = One if batch i contains product p; otherwise 0

REFERENCES

(1) Cafaro, D. C.; Cerdá, J. Optimal Scheduling of Multiproduct Pipeline Systems Using a Non-Discrete MILP Formulation. Comput. Chem. Eng. 2004, 28, 2053−2068. (2) Trench C. J. How Pipelines Make the Oil Market Work: Their Networks, Operation and Regulation; Allegro Energy Group: New York, 2001. (3) Milidiu, R. L.; Liporace, F. D. S. Planning of Pipeline Oil Transportation with Interface Restrictions Is a Difficult Problem. 2003 PUC-Rio Inf.MCC56/03. (4) Sasikumar, M.; Prakash, P. R.; Patil, S. M.; Ramani, S. PIPES: A Heuristic Search Model for Pipeline Schedule Generation. Knowl.Based. Syst. 1997, 10, 169−175. (5) Shah, N. Mathematical Programming Techniques for Crude Oil Scheduling. Comput. Chem. Eng. 1996, 20, 1227−1232. (6) Techo, R.; Holbrook, D. L. Computer Scheduling the World’s Biggest Product Pipeline. Pipe. Gas. J. 1974, 4, 27−34. (7) Hane, C. A.; Ratliff, H. D. Sequencing Inputs to MultiCommodity Pipelines. Ann. Oper. Res. 1995, 57, 73−101. (8) Reddy, P. C. P; Karimi, I. A; Srinivasan, R. A New ContinuousTime Formulation for Scheduling Oil Operations. Chem. Eng. Sci. 2004, 59, 1325−1341. (9) Neiro, S. M. S; Pinto, J. M. A General Modeling Framework for the Operational Scheduling of Petroleum Supply Chains. Comput. Chem. Eng. 2004, 28, 871−896. (10) Rejowski, R. J.; Pinto, J. M. Scheduling of a Multiproduct Pipeline System. Comput. Chem. Eng. 2003, 27, 1229−1246. (11) Rejowski, R. J.; Pinto, J. M. Efficient MILP Formulations and Valid Cuts for Multiproduct Pipeline Scheduling. Comput. Chem. Eng. 2004, 28, 1511−1528. (12) Rejowski, R., Jr.; 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−1066. (13) Cafaro, D. C.; Cerdá, J. Dynamic Scheduling of Multiproduct Pipelines with Multiple Delivery Due Dates. Comput. Chem. Eng. 2008, 32, 728−753. (14) Relvas, S.; Matos, H. A.; Barbosa-Povoa, 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−7855. (15) Cafaro, D. C.; Cerdá, J. Efficient Tool for the Scheduling of Multiproduct Pipeline and Terminal Operation. Ind. Eng. Chem. Res. 2008, 47, 9941−9956. (16) Relvas, S.; Matos, H. A.; Barbosa-Povoa, A. P. F. D.; Fialho, J.; Pinheiro, A. S. Reactive Scheduling Framework for a Multiproduct Pipeline with Inventory Management. Ind. Eng. Chem. Res. 2007, 46, 5659−5672. (17) Relvas, S.; Barbosa-Povoa, A. P. F. D.; Matos, H. A. Heuristic Batch Sequencing on a Multiproduct Oil Distribution Center. Comput. Chem. Eng. 2013, 41, 955−968. (18) Magatão, L.; Arruda, L. V. R.; Neves, F. A. Programming Approach for Scheduling Commodities in a Pipeline. Comput. Chem. Eng. 2004, 28, 171−185. (19) Magatão, L.; Arruda, L. V. R.; Neves, F. A. A Combined CLPMILP Approach for Scheduling Commodities in a Pipeline. J. Sched. 2010, 14, 57−87. (20) Mirhassani, S. A.; Ghorbanalizadeh, M. The Multi-Product Pipeline Scheduling System. Comput. Math. Appl. 2008, 56, 891−897. (21) Mirhassani, S. A.; Fani Jahromi, H. Scheduling Multi-Product Tree-Structure Pipelines. Comput. Chem. Eng. 2011, 35, 165−176.

Continuous Variables

Lk = duration of pumping run k LSkj = activity duration of segment j during the pumping run k LDkj = activity duration of removal terminal j during the pumping run k STk = starting time of pumping run k Fki = upper coordinate of batch i at the end of the pumping run k Wki = size of batch i inside the pipeline at the end of the pumping run k Abki = amount of batch i injected into the pipeline through the pumping run k APki,p = amount of batch i containing product p inserted into the pipeline through the pumping run k Adki,j = amount of batch i diverted into the removal terminal j through the pumping run k VPki,p,j = amount of batch i containing product p diverted into removal terminal j through the pumping run k INVi,p,q = interface volume between batches i and (i + 1) if they contain products p and q ACPk = active volume of the pipeline through the pumping run k AVk = pipeline volume activated by run k SVk = pipeline volume stopped at run k undemtp,j = overall unsatisfied demands of product p for removal terminal j at the end of period t to meet at period t +1 Binary Variables

zki = One if and only if a portion of batch i is injected to the pipeline during the pumping run k; otherwise 0 vkj = One if and only if the segment j is active during the pumping run k; otherwise 0 wkj = One if and only if the removal terminal j receives the material from pipeline during the pumping run k otherwise 0 7041



Article

(22) Cafaro, D. C.; Cerdá, J. A Rigorous Mathematical Formulation for the Scheduling of Tree-Structure Pipeline Networks. Ind. Eng. Chem. Res. 2011, 50, 5064−5085. (23) Castro, P. M. Optimal Scheduling of Pipeline Systems with a Resource-Task Network Continuous-Time Formulation. Ind. Eng. Chem. Res. 2010, 49, 11491−11505. (24) Herrán, A.; de la Cruz, J. M.; de Andrés, B. Mathematical Model for Planning Transportation of Multiple Petroleum Products in a Multi-Pipeline System. Comput. Chem. Eng. 2010, 34, 401−413. (25) Herrán, A.; de la Cruz, J. M.; de Andrés, B. Global Search Metaheuristics for Planning Transportation of Multiple Petroleum Products in a Multi-Pipeline system. Comput. Chem. Eng. 2012, 37, 248−261. (26) Cafaro, D. C.; Cerdá, J. Optimal Scheduling of Refined Products Pipelines with Multiple Sources. Ind. Eng. Chem. Res. 2009, 48 (14), 6675−6689. (27) Cafaro, D. C.; Cerdá, J. Operational Scheduling of Refined Products Pipeline Networks with Simultaneous Batch Injections. Comput. Chem. Eng. 2010, 28 (10), 1687−1704. (28) Mirhassani, S. A.; Abbasi, M.; Moradi, S. Operational Scheduling of Refined Product Pipeline with Dual Purpose Depots. App. Math. Model. 2013, 37 (8), 5723−5742. (29) Cafaro, V. G.; Cafaro, D. C.; Mendéz, C. A.; Cerdá, J. Effective Coordination of Simultaneous Delivery Flows into Receipt Terminals of Multiproduct Pipelines. Comput. Aided Process Eng. 2012, 30, 252− 256. (30) Cafaro, V. G.; Cafaro, D. C.; Mendéz, C. A.; Cerdá, J. Detailed Scheduling of Single- Source Pipelines with Simultaneous Deliveries to Multiple Offtake Stations. Ind. Eng. Chem. Res. 2012, 51, 6145−6165.

7042


A General Modeling Framework for the Long-Term Scheduling of

Recommend Documents