Computationally Efficient MILP model for Scheduling a Branched

Mar 11, 2019 - As one essential key performance indicator in continuous-time models, a lower number of event points brings a significant reduction in ...
3 downloads 0 Views 2MB Size
Subscriber access provided by UNIV OF CAMBRIDGE

Process Systems Engineering

Computationally Efficient MILP model for Scheduling a Branched Multiproduct Pipeline System Qi Liao, Pedro Miguel Castro, Yongtu Liang, and Haoran Zhang Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b06490 • Publication Date (Web): 11 Mar 2019 Downloaded from http://pubs.acs.org on March 18, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Computationally Efficient MILP model for Scheduling a Branched Multiproduct Pipeline System Qi Liao , Pedro M. Castro, Yongtu Liang and Haoran Zhang a,b

a

b,1

a

c

China University of Petroleum-Beijing, Beijing Key Laboratory of Urban oil and Gas Distribution Technology, Fuxue Road No.18, Changping District, Beijing 102249, P. R. China b

Centro de Matemática Aplicações Fundamentais e Investigação Operacional, Faculdade de Ciências, Universidade de Lisboa, 1749-016 Lisboa, Portugal

c

Center for Spatial Information Science, The University of Tokyo 5-1-5 Kashiwanoha, Kashiwashi, Chiba 277-8568, Japan

Abstract This paper presents a novel mixed-integer linear programming (MILP) continuous-time formulation for the detailed scheduling of a branched pipeline system with a single refinery and multiple depots, rigorously accounting for forbidden product sequences in every line. The major advantage against previous published work is that the proposed model allows multiple batches to be processed by a node over a single slot, leading to the requirement of fewer event points to represent a schedule. As one essential key performance indicator in continuous-time models, a lower number of event points brings a significant reduction in computational time. Five benchmark problems from the literature involving different pipeline configurations are tested to validate the superiority of the proposed model. In four of them, we were able to find new best solutions.

1

Corresponding author. E-mail: [email protected]. 1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1.

Page 2 of 34

INTRODUCTION A multiproduct pipeline network is a complex and large-scale system composed of refineries,

injection, delivery and pumping stations, storage tanks and connecting pipelines. Pipeline operators should take downstream market demand, refinery production capacity, conveying capacity and storage profiles along the pipeline into account to rationally allocate resources and work out the detailed scheduling plan. The scheduling plan is evaluated and verified from the angle of production safety and feasibility, so as to reach production goals such as minimum makespan , minimal operational cost , 1, 2

3-6

most stable operating condition and the highest market satisfaction . In addition, since the operations 7

8

of multiproduct pipelines are affected by downstream market demand, the scheduling plans need often to be modified during execution, to react to the changing market environment and technology . 9, 10

Moreover, when adjacent batches are conveyed without separation devices, contaminated products are inevitably generated under the diffusion action of convection and turbulence. The tracking and control of contaminant products is one of the key technologies for multiproduct pipeline systems. However, the technical processes and operations management related with contaminated products are complicated, especially in regions with complex topography and large elevation difference. These products cannot be sold as qualified products and need to be processed by blending or refining, which brings extra operation cost for petroleum companies . 11

The scheduling optimization of multiproduct pipelines involves the following subproblems: lot sizing and batch sequence, injection, delivery and pump operations, inventory management and contaminated product control. Each of them involves a large number of factors and even some nonlinear and uncertain items, thus making it challenging to deal with the scheduling optimization . 12

At the same time, a complex pipeline structure and long scheduling horizon greatly increase the model scale, resulting in high computational effort . Therefore, it makes sense for pipeline management to 13, 14

quickly draw the economic and feasible scheduling plans according to the fluctuating market and production environment. Over the past decade, mixed-integer linear (MILP) and nonlinear programming models (MINLP) have been widely used for scheduling. Cafaro et al.

, Castro et al. , Relvas et al. , and MirHassani

3, 15-18

2

ACS Paragon Plus Environment

19, 20

21, 22

Page 3 of 34

Industrial & Engineering Chemistry Research

et al

23, 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

carried out extensive research with continuous-time MILP models for different pipeline

structures. Some studies also presented MINLP models coupled with hydraulic calculation to ensure 25-28

safe and efficient operation. However, MINLP models are nonconvex, so it is difficult to find the optimal solution within an acceptable time range or even a feasible solution. To resolve this problem, some solution strategies gradually appeared, such as decomposition, aggregation and heuristic algorithms. For instance, Rejowski and Pinto

25

and Cafaro et al. used GAMS- DICOPT to solve 27

29

MINLPs accounting for pump energy consumption. By comprehensively considering constraints such as injection, delivery, batch migration, power consumption of pumps and time-dependent electricity price, Zhang et al established a continuous-time MINLP model, then decomposed it into a sequence 28

model and a time-volume model, and finally solved these two models iteratively by a hybrid algorithm to find the optimal or sub-optimal solution of the original problem. Most previous works avoided batch contamination by either optimizing batch sequence or adding constraints to the mathematical models. Research on the former has received a great deal of attention over the last decade . Aiming at batch sequencing and detailed scheduling of branched pipeline 30-32

networks, Mostafaei, et al.

33

presented a unified optimization model considering processing costs

dependent on the volume and type of contamination. Forbidden batch sequences were also taken into account to prevent products with large differences in physical properties to be conveyed adjacently. However, this approach cannot handle multilevel branches. More recently, and based on their previous work , Castro and Mostafaei 34

35

proposed a modular and simpler batch-centric approach for treelike

networks with a single input node that rigorously prevents forbidden sequences. Compared to their product-centric approach , it allows multiple batches to enter/leave a segment during a time slot and 1

finds the same or a better solution for a given number of time slots. In addition, Castro and Mostafaei 35

also contribute to reducing product contamination by avoiding idle segments with multiple products

inside. More sophisticated operations for controlling contamination are considered in Refs. . By regarding 8, 36

the events of batches just arriving at stations as key event points, two significant sets of constraints can be taken into account: (i) more rigorous lower limits are imposed on the flowrate for pipeline 3

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 34

segments with contaminated products; (ii) shutdown operation should be forbidden as far as possible in inclined pipeline segments, especially when heavy oil is higher and light oil is lower. However, more event points are required in the model involving these processing constraints since flowrate constraints for segments change along with the arrival of different mixture interfaces. Aiming at the optimal scheduling of a branched pipeline system with a single refinery, this paper presents a new batch-centric MILP model for systems with product-independent injection and delivery flowrates. Compared with previous work , our model improves in the following aspects: (i) multiple 35

batches can now be injected by an input node and delivered to an output node during a time slot, to generate the optimal solution with fewer time slots; (ii) it tracks the coordinate of batches in lines rather than in segments (a line may include multiple segments), which is a common option in the literature

24, 37

to reduce the problem size. Both aspects lead to an improved computational performance,

with feature (i) requiring the development of novel constraints to rigorously enforce forbidden product sequences. The rest of the paper is organized as follows. Section 2 describes the features of the proposed model. Section 3 gives a short explanation about the main novelty of the proposed model. Section 4 presents the detailed mathematical formulations. Section 5 focuses on five benchmark problems with different pipeline configurations and compares the numerical results with the ones reported in the literature . 35, 37

Finally, the conclusions are provided in Section 6. 2.

PROBLEM STATEMENT This paper deals with a branched multiproduct pipeline system with a single refinery and multiple

depots. Such a system is composed of a set of lines 𝑙 ∈ 𝐿 of volume 𝑣%& connected by junctions 𝑗 ∈ 𝐽. Each line comprises a single input node, |𝑁%+ |=1, and one or more output nodes 𝑛 ∈ 𝑁%- , located at a volume 𝜎/ from the origin. The input node of the mainline (N1) is the refinery 𝑁 0 = {N1} and there is a terminal node on each line (subsets 𝑁%6 ). Output nodes can be depots (subset 𝑁 7 ) or be part of a junction (𝑁8- ) involving an input node (𝑁8+ ), see Figure 1. It should be noted that each node in the system must be an element of one of four subsets 𝑁 0 , 𝑁 7 , 𝑁8+ , and 𝑁8- . If a depot 𝑛 is located at a 4

ACS Paragon Plus Environment

Page 5 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

junction 𝑗, we need to create two output nodes (𝑛′ and 𝑛′′) in the same location, with one representing the depot (𝑛′ ∈ 𝑁 7 ) and the other representing the output node of the junction (𝑛′′ ∈ 𝑁8- ). A dummy segment with a capacity of zero is also generated between these two output nodes, as is the case of EX2 in Section 5.2. On the contrary, if this junction does not include a depot, we can infer that 𝑛 ∈ 𝑁8and need not make a division of this node. Moreover, all nodes have a storage system, except for junction nodes, which do not store products but only transfer products from their connected upstream lines to downstream lines. The refinery node 𝑛 pumps products into the main line at a flowrate belonging to the interval [𝜌/=,?@A , 𝜌/=,?BC ] for the depots to receive them later. Note that different depots can receive products at the same time.

Figure 1. Nodes types in a branched pipeline system. We assume that the initial volumes of all products (𝑝 ∈ 𝑃) inside the system are known. Different batches can initially exist on a given line, with the user allocating the so-called old batches 𝑖 ∈ 𝐼 I%J to &,N &,N products through parameter 𝑦L,M . We can then compute the initial volumes (𝑣%,L ), as well as left (𝑙𝑐%,L ) &,N and right coordinates ( 𝑟𝑐%,L ). Specific product sequences are forbidden to appear in all system

segments 𝑠 ∈ 𝑆 . These are specified in subsets 𝑃M , with 𝑝´ ∈ 𝑃M meaning that 𝑝´ cannot enter a =,?@A =,?BC segment after 𝑝. The inventory in storage at refinery and depot nodes must be within [𝑣/,M , 𝑣/,M ].

The objective is to find the volumes and timing of product injections from the refinery and product =,UAV deliveries to depot nodes that meet the given demand 𝑓/,M in the shortest time possible.

5

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 34

2.1. Representation of Pipeline Network A motivating example is given in Figure 2 to provide a better understanding of the crucial sets and subsets involved in this model. This branched pipeline network transports three products (P1-P3) from a single refinery N1 to five depots (N3, N5, N7, N10 and N12) along a mainline L1 and three branch lines L2-L4. Branch lines L2 and L3 originate from mainline L1 at junctions J1 (i.e. N2 and N6) and J2 (i.e. N4 and N8), respectively. While a level-two branch L4 connects to branch L3 by junction J3 (i.e. N9 and N11). The bottom of Figure 2 lists the sets, subsets and parameters indicating pipeline network topology and node locations. Beyond that, these four lines are divided into a total of eight segments (S1-S8) by twelve nodes (N1-N12). Single-element subsets 𝑆/W (𝑆/7 ) are introduced to specify the adjacent upstream (downstream) segment of node 𝑛 ∈ 𝑁. Notice that any initial node (terminal node) along a line has no upstream (downstream) segment.

Figure 2. Elements of sets and subsets used by our model for a motivating example. 3.

MAIN NOVELTY The major advantage of the batch-centric model in Castro and Mostafaei against their product35

centric model is that multiple products can enter/leave a segment during a time interval. Therefore, 1

the computational performance is improved greatly because fewer event points are required to find the 6 ACS Paragon Plus Environment

Page 7 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

same or a better solution to a problem. Based on these findings, we extend their research by allowing multiple batches to be processed by a node during a slot while guaranteeing the accuracy of the schedules. Let us assume that two batches (I2 and I3) are being diverted to intermediate depot N2 during slot = 𝑡, as shown in Figure 3. In scenario (a), N2 is receiving the batches at a rate of 𝜌=Y,Z and the flowrate [ in its adjacent downstream segment S2 is 𝜌[Y,Z . To generate a valid split between depot N2 and

segment S2 (check Ref. for further details) the transported volumes must satisfy the constraint 35

],a \]^,_^,` b \b^,_^,`

\ ],a

[,J = \]^,_c,` , which brings a nonlinear relationship between the variables (note that variables 𝐹e,L,Z b b^,_c,`

can be generated by writing eq. (29) for every batch). To avoid using this nonlinear constraint, we restrict an intermediate output node with an active downstream segment, to receive at most one batch during a slot. Notice that in scenario (b), segment S2 is idle and so the two batches can be delivered to N2.

Figure 3. A motivating example for multiple batches being delivered to an intermediate node. It should be highlighted that the condition of no flow to a downstream segment of an intermediate node receiving more than one batch per slot, does not lead to suboptimal schedules. This is because one can always increase the number of slots to make it possible for multiple batches to be delivered to the intermediate node while also transferring to the downstream segment. It is just that the different batch deliveries will be assigned to multiple consecutive slots instead of a single one. 7

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 34

Conversely, all refinery nodes with an active adjacent upstream segment can inject at most one batch during a slot. This strategy will be implemented in future work since this paper focuses on a network with a single refinery. 4.

MATHEMATICAL MODEL This model is developed based on the continuous-time representation of Figure 4. Continuous

variable 𝑇Z gives the absolute time of event point 𝑡, while 𝐷Z holds the duration of time slot 𝑡. Slot 𝑡 begins at event point 𝑡 and ends at event point 𝑡+1, leading to eq. (1). The absolute time of all event points is bounded by the given time horizon ℎ, but it is enough to enforce this constraint for the last event point, as seen in eq. (2). 𝑇Zij = 𝑇Z + 𝐷Z ∀𝑡 < |𝑇|

(1)

𝑇|6| ≤ ℎ

(2) slot 1

Time grid

2

1

Variables

T1

time slot 2

D1

event points

T2

D2

slot |T|-2 3

|T|-2

T3

T|T|-2

slot |T|-1 |T|-1

D|T|-2

T|T|-1

t=|T|

T|T| D|T|-1

Timing

h

0

Figure 4. The model relies on a continuous-time representation with a single grid. 4.1. Objective function The objective function discussed in this paper is to minimize the makespan, the absolute time of the last event point in the grid. Completing the transportation task in the shortest possible time helps to achieve the maximum utilization of the pipeline capacity. Another popular alternative

3, 13, 37

is to

minimize the operational cost, which includes pump cost, backorder cost, cost of activated/stopped pipeline volumes, interface cost and inventory cost. min 𝑇|6|

(3)

4.2. Coordinate System One difference with respect to the work of Castro and Mostafaei , is that this paper traces batch 35

migration by means of line coordinates, instead of segment coordinates. Taking line L1 in Figure 2 as 8

ACS Paragon Plus Environment

Page 9 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

an example, we show in Figure 5 its initial left (𝑙𝑐 N ) and right (𝑟𝑐 N ) segment/line coordinates. It can be observed that the number of required parameters for line coordinates is significantly less, and so is the number of corresponding model variables: 𝐿𝐶%,L,Z and 𝑅𝐶%,L,Z ; for left and right coordinates, respectively.

Figure 5. Difference between segment and line coordinates. & Let continuous variable 𝑉%,L,Z hold the volume of batch 𝑖 in line 𝑙 at event point 𝑡. Since we are

dealing with incompressible fluids, the total volume of all batches inside a line must be equal to the line capacity 𝑣%& (eq. (4)). In any line, the right coordinate for batch 𝑖 is computed by summing the volume of all batches 𝑖´ ≥ 𝑖, as stated by eq. (5). Eq. (6) calculates the batch volume inside the line as the difference between its right and left coordinates. The same equation can be used on the initial &,N coordinates to compute the initial volumes 𝑣%,L that appear in eq. (7). The latter is a multiperiod

balance that computes the volume at event point 𝑡 as the volume at the previous event minus the total volume of the batch leaving the line through output nodes during slot 𝑡-1, plus the amount being injected by its input node during that same slot. & ∑L∈+w 𝑉%,L,Z = 𝑣%& ∀𝑙, 𝑡

(4)

& 𝑅𝐶%,L,Z = ∑L´∈+w :L´yL 𝑉%,L´,Z ∀𝑙, 𝑖 ∈ 𝐼% , 𝑡

(5)

& 𝑉%,L,Z = 𝑅𝐶%,L,Z − 𝐿𝐶%,L,Z ∀𝑙, 𝑖 ∈ 𝐼% , 𝑡

(6)

&,N & 𝑉%,L,Z = 𝑣%,L {

Z|j

& = = + 𝑉%,L,Z}j − ∑/∈=~ :L∈+• 𝐹/,L,Z}j + ∑/∈=_ :L∈+• 𝐹/,L,Z}j ∀𝑙, 𝑖 ∈ 𝐼% , 𝑡 w

w

9

ACS Paragon Plus Environment

(7)

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 34

4.3. Nodes =,LJ%• = =1 if node 𝑛 is receiving/sending batch 𝑖 during time slot 𝑡 and 𝑋/,Z =1 Let binary variable 𝑋/,L,Z

otherwise. Movement of the batch requires the transferred volume to belong to interval [𝑓/=,?@A , 𝑓/=,?BC ] and the flowrate to be between the minimum 𝜌/=,?@A and maximum 𝜌/=,?BC values. Conversely, if the node is idle during slot 𝑡, the volume is zero for every batch and the constraint on the slot duration 𝐷Z is relaxed. The two alternatives are modelled through the disjunction

38, 39

in eq. (8).

Note that the disjunction is inclusive over all batches that can appear in the node (𝑖 ∈ 𝐼/ ), thus allowing multiple batches to be transferred during slot 𝑡. This is a major difference compared to the recent model by Castro and Mostafaei , which used an exclusive disjunction to restrict the injection/delivery 35

to a single batch per time slot. = 𝑋/,L,Z =,LJ%• ⎡ ⎤ 𝑋/,Z =,?@A =,?BC = ≤ 𝐹/,L,Z ≤ 𝑓/ ⎢𝑓/ ⎥ = ‚⎢ ‚ ‰𝐹/,L,Z = 0 ∀𝑖 ∈ 𝐼/ ‹ ∀𝑛, 𝑡 < |𝑇| = = 𝐹/,L,Z 𝐹/,L,Z ⎥ L∈+• ⎢ ⎥ 𝐷Z ≤ ℎ =,?BC ≤ 𝐷Z ≤ =,?@A 𝜌/ ⎣ 𝜌/ ⎦

(8) The logic proposition in eq. (8) is converted into eq. (9), while the convex hull reformulation of 34

40

the constraints inside the disjunction generates eqs. (10)-(11). =,LJ%• = 𝑋/,L,Z + 𝑋/,Z ≤ 1 ∀𝑛, 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(9)

= = = 𝑓/=,?@A 𝑋/,L,Z ≤ 𝐹/,L,Z ≤ 𝑓/=,?BC 𝑋/,L,Z ∀𝑛, 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇| b ∑Œ∈_• \•,Œ,` b,Ž•• ••

b ∑Œ∈_• \•,Œ,`

≤ 𝐷Z ≤

b,Ž‘’

••

=,LJ%• + ℎ ∙ 𝑋/,Z ∀𝑛 ∈ 𝑁 0 ∪ 𝑁 7 , 𝑡 < |𝑇|

(10) (11)

4.3.1. Input Node The input node is the left-most node of a line and so lies at coordinate zero. Thus, for a batch to enter the line during slot 𝑡, its left coordinate must be zero at event point 𝑡. Conversely, if the node is idle during slot 𝑡, the constraint on the batch coordinate is relaxed, as seen in eq. (12). ‚• L∈+•

=,LJ%• = 𝑋/,Z 𝑋/,L,Z –‚• – ∀𝑙, 𝑛 ∈ 𝑁%+ , 𝑡 < |𝑇| & 𝐿𝐶%,L,Z = 0 𝐿𝐶%,L,Z ≤ 𝑣% ∀𝑖 ∈ 𝐼/

(12) 10

ACS Paragon Plus Environment

Page 11 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

The convex hull reformulation of the constraint inside the disjunction generates eq. (13). = 𝐿𝐶%,L,Z ≤ 𝑣%& ∙ —1 − 𝑋/,L,Z ˜ ∀𝑙, 𝑛 ∈ 𝑁%+ , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(13)

4.3.1.1. Refinery Node In the network configuration considered in this work, the refinery is the input node of the mainline. Eq. (14) computes the volume of batch 𝑖 leaving the refinery to enter the pipeline, by the sum over all =,J . Eq. (14) also applies to depot products that exist at the node (𝑝 ∈ 𝑃/ ) of disaggregated variables 𝐹/,L,M,Z

nodes. =,J = 𝐹/,L,Z = ∑M∈™• 𝐹/,L,M,Z ∀𝑛 ∈ 𝑁 0 ∪ 𝑁 7 , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(14)

=,N The inventory balance of the refinery’s storage system is given in eq. (15), where parameter 𝑣/,M =,L/ represents the initial inventory of product 𝑝 and 𝜌/,M represents the incoming flowrate of product 𝑝.

The volume at event point 𝑡, is equal to the volume at the previous event point plus the incoming volume during slot 𝑡 − 1 and minus the volume injected in the pipeline during that same slot through all batches 𝑖 that contain product 𝑝. =,N =,L/ =,J = = 𝑉/,M,Z = 𝑣/,M {Z|j + 𝑉/,M,Z}j + 𝜌/,M 𝐷Z}j − ∑L∈+• 𝐹/,L,M,Z}j ∀𝑛 ∈ 𝑁 0 , 𝑝 ∈ 𝑃/ , 𝑡

(15)

At any event point 𝑡, the volume of product 𝑝 in the storage system of a refinery or depot node 𝑛 is subject to given lower and upper bounds. =,?@A =,?BC = 𝑣/,M ≤ 𝑉/,M,Z ≤ 𝑣/,M ∀𝑛 ∈ 𝑁 0 ∪ 𝑁 7 , 𝑝 ∈ 𝑃/ , 𝑡

(16)

4.3.2. Terminal Output Node Output nodes can be distinguished between terminal and intermediate nodes. For terminal output = node 𝑛, batch 𝑖 can be delivered during slot 𝑡 (𝑋/,L,Z =1), if its right coordinate at the end of the slot

(event point 𝑡+1) is equal to the line volume. Eq. (17) can be reformulated into eq. (18). =,LJ%• = 𝑋/,L,Z 𝑋/,Z ‚• – ∀𝑙, 𝑛 ∈ 𝑁%6 , 𝑡 < |𝑇| &– ‚ • & 𝑅𝐶%,L,Zij = 𝑣% 𝑅𝐶%,L,Zij ≤ 𝑣% ∀𝑖 L∈+•

(17) = 𝑅𝐶%,L,Zij ≥ 𝑣%& 𝑋/,L,Z ∀𝑙, 𝑛 ∈ 𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(18) 11

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 34

4.3.2.1. Depot Node The terminal node of a line is also a depot. The inventory balance of a depot’s storage system is like the refinery’s balance but now the pipeline flow is incoming. There is also another term on the right=,•/J at the last event point. hand side of eq (19) to remove the given product demand 𝑓/,M =,N =,J =,•/J = = 𝑉/,M,Z = 𝑣/,M {Z|j + 𝑉/,M,Z}j + ∑L∈+• 𝐹/,L,M,Z}j − 𝑓/,M {Z||6| ∀𝑛 ∈ 𝑁 7 , 𝑝 ∈ 𝑃/ , 𝑡

(19)

4.3.3. Intermediate Output Node The constraints for non-terminal output nodes are more complex. Batch 𝑖 can be delivered to intermediate output node 𝑛 during slot 𝑡 if its right coordinate at the end of the slot has gone past the location of the node 𝜎/ and its left coordinate has not, as seen in the top-left of eq. (20). Then, we need =,j to distinguish if batch 𝑖 is the first to arrive during slot 𝑡 (𝑋/,L,Z =1), or if it arrives in second or later =,Y (𝑋/,L,Z =1). In the former case, we must have the right coordinate at 𝑡, at the node location or to its right.

Otherwise, the constraint on 𝑅𝐶%,L,Z is relaxed to allow the right coordinate of batch 𝑖 to start to the left of the node, as seen in the bottom-left of eq. (20). Note that multiple batches can be received only if the downstream segment is idle, which explains why the left coordinate of the first received batch will not move past 𝜎/ . The right-hand side of eq. (20) relaxes the coordinate values if the node is idle =,LJ%• (𝑋/,Z =1). = 𝑋/,L,Z =,LJ%• ⎡ ⎤ 𝑋/,Z ⎡ ⎤ 𝑅𝐶 ≥ 𝜎 %,L,Zij / ⎢ ⎥ & ⎢ 𝑅𝐶 ≤ 𝑣 ∀𝑖 ⎥ %,L,Zij % 𝐿𝐶%,L,Zij ≤ 𝜎/ 6 ⎥‚ ‚⎢ ⎢ ⎥ ∀𝑙, 𝑛 ∈ 𝑁% \𝑁% , 𝑡 < |𝑇| & =,j =,Y 𝐿𝐶 ≤ 𝑣 ∀𝑖 ⎢ ⎥ % 𝑋/,L,Z 𝑋/,L,Z L∈+• ⎢ %,L,Zij ⎥ –‚• –⎥ & ⎢• ⎣ 𝑅𝐶%,L,Z ≤ 𝑣% ∀𝑖 ⎦ 𝑅𝐶%,L,Z ≥ 𝜎/ 𝑅𝐶%,L,Z ≥ 0 ⎣ ⎦

(20) Eq. (20) features embedded disjunctions on the left-hand side, which can be reformulated following the guidelines given in Vecchietti and Grossmann . The convex hull reformulation of eq. (20) 41

generates eqs. (21)-(24). =,j =,Y = 𝑋/,L,Z = 𝑋/,L,Z + 𝑋/,L,Z ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(21)

= 𝑅𝐶%,L,Zij ≥ 𝜎/ 𝑋/,L,Z ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(22)

12

ACS Paragon Plus Environment

Page 13 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research = = 𝐿𝐶%,L,Zij ≤ 𝜎/ 𝑋/,L,Z + 𝑣%& —1 − 𝑋/,L,Z ˜ ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(23)

=,j 𝜎/ 𝑋/,L,Z ≤ 𝑅𝐶%,L,Z ≤ 𝑣%& ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑡 < |𝑇|

(24)

Another useful constraint is implicit in the disjunction in eq. (20). Eq. (25) states that at most one batch can be received first by the depot during slot 𝑡. =,j ∑L∈+• 𝑋/,L,Z ≤ 1 ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑡 < |𝑇|

(25)

An output node can receive multiple batches only when its adjacent downstream segment 𝑠 ∈ 𝑆/7 is [ idle (𝑋e,Z = 0). This condition is stated in eq. (26), which can be reformulated into eq. (27). =,Y [ 𝑋/,L,Z ⇒ ¬𝑋e,Z ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑠 ∈ 𝑆/7 , 𝑡 < |𝑇|

(26)

=,Y [ 𝑋/,L,Z + 𝑋e,Z ≤ 1 ∀𝑙, 𝑛 ∈ 𝑁%- \𝑁%6 , 𝑖 ∈ 𝐼/ , 𝑠 ∈ 𝑆/7 , 𝑡 < |𝑇|

(27)

4.4. Junction Junction 𝑗 comprises an output node 𝑛 ∈ 𝑁8- that transfers material to an input node 𝑛´ ∈ 𝑁8+ that starts a new line. Eq. (28) simply states that the volume balance for batch 𝑖 and slot 𝑡 is satisfied. = = 𝐹/,L,Z = 𝐹/´,L,Z ∀𝑗, 𝑛 ∈ 𝑁8- , 𝑛´ ∈ 𝑁8+ , 𝑖 ∈ 𝐼/´ , 𝑡 < |𝑇|

(28)

4.5. Segments To enforce minimum 𝜌e[,?@A and maximum 𝜌e[,?BC flowrate constraints on segments, we first need to compute the total volume transported in segment 𝑠 during slot 𝑡. This can be done through eq. (29). More specifically, and for the first node in the line, there are no upstream segments (𝑆/W = ∅) and so = [ the injected volume goes all into its downstream segment, i.e. ∑L 𝐹/,L,Z = 𝐹e´,Z (𝑠´ ∈ 𝑆/7 ). For the

remaining nodes in the line, output nodes, we subtract the volume leaving the node from the volume leaving the upstream segment (𝑠 ∈ 𝑆/W ), to compute the volume entering the downstream segment = [ [ during slot 𝑡, i.e., 𝐹e,Z − ∑L 𝐹/,L,Z = 𝐹e´,Z . = [ = [ 𝐹e,Z {e∈[ ž + ∑L∈+• :/∈=_ 𝐹/,L,Z = 𝐹e´,Z {e´∈[ Ÿ + ∑L∈+• :/∈=~ 𝐹/,L,Z ∀𝑙, 𝑛 ∈ 𝑁% , 𝑡 < |𝑇| •

w

w



(29)

[ Movement through the segment activates binary variables 𝑋e,Z , as seen in eq. (30). Then, eq. (31)

enforces the flowrate constraints.

13

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

[ [ [ 𝑓e[,?@A 𝑋e,Z ≤ 𝐹e,Z ≤ 𝑓e[,?BC 𝑋e,Z ∀𝑠, 𝑡 < |𝑇|



\ ],` ],Ž••

≤ 𝐷Z ≤



\ ],` ],Ž‘’

Page 14 of 34

(30)

[ + ℎ ∙ —1 − 𝑋e,Z ˜ ∀𝑠, 𝑡 < |𝑇|

(31)

4.6. Batch to Product Assignments Let binary variable 𝑌L,M =1 if batch 𝑖 is assigned to product 𝑝. Eq. (32) states that every postulated batch needs to be assigned to exactly one product. Assignment of old batches is done by the user (parameter 𝑦L,M =1) rather than by the optimization, leading to eq. (33). ∑M 𝑌L,M = 1 ∀𝑖 ∈ 𝐼

(32)

𝑌L,M = 𝑦L,M ∀𝑖 ∈ 𝐼 I%J , 𝑝

(33)

Let parameter 𝑓M™,?BC be an upper bound on the maximum amount of product 𝑝 that can be transferred over the time horizon. This bound is used to force to zero all disaggregated variables that do not involve batch 𝑖 and its assigned product 𝑝, for refinery and depot nodes, as seen in eq. (34). =,J ∑/∈=£ ∪=Ÿ :L∈+• ,M∈™• ∑Z¢|6| 𝐹/,L,M,Z ≤ 𝑓M™,?BC 𝑌L,M ∀𝑖, 𝑝

(34)

4.7. Preventing Forbidden Sequences Considering product quality, some products are forbidden to be transported adjacently inside the pipeline. Let subset 𝑃M include products 𝑝´ that cannot follow product 𝑝. Most previous works enforce constraints on the pumping sequence at the left-most input node, where batches are created, i.e., batch 𝑖+1 cannot be assigned to product 𝑝´ if batch 𝑖 is assigned to product 𝑝, as stated by eq. (35). However, forbidden sequences may not be entirely avoided due to the appearance of empty batches in the pipeline. 𝑌L,M + ∑M´∈™¤ 𝑌Lij,M´ ≤ 1 ∀𝑖, 𝑝

(35)

Figure 6 illustrates two scenarios where the formation of an empty batch leads to forbidden sequences in the branched pipeline network. In (a), batch I2 is fully delivered to an intermediate node during slot 𝑡, leading to forbidden sequence P2-P1 in the node’s downstream segment. In scenario (b), when encountering the branch, all I2 remains in line L1, while batch I3 is directed to line L2, leading to the delivery of products P2 (I1) and P1 (I3) in sequence to the terminal node of line L2. 14

ACS Paragon Plus Environment

Page 15 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 6. Appearance of empty batches in the pipeline system leading to forbidden sequences. Castro and Mostafei prevent forbidden sequences after identifying the left-most batch in a segment 35

by means of its non-zero internal volume and zero left coordinates, which works well under their model assumption of a single batch delivery to a depot during a time slot. However, such strategy needs to be adapted to cover the possibility of multiple deliveries. In the bottom row of Figure 6, for scenario (b), we can see that batch I3 is not inside line L2 at event point 𝑡. It then enters line L2 during & & slot 𝑡, leaving it completely at the other end during the same slot. Thus, 𝑉&Y,+¥,Z =𝑉&Y,+¥,Zij =0, even

though batch I3 travels through line L2 during slot 𝑡. In this work, forbidden sequences are prevented through novel constraints involving binary variables &,J 𝑋%,L,M,Z , which are equal to 1 when batch 𝑖 containing product 𝑝 exists or crosses line 𝑙 during slot 𝑡. & Batch 𝑖 is involved in line 𝑙 during slot 𝑡, if it is already there at the start of the slot (𝑉%,L,Z > 0), or if it = has entered the line through one of its input nodes 𝑛 (𝐹/,L,Z > 0) during slot 𝑡. In such cases, eq. (36) &,J together with eq. (41), activates the relevant 𝑋%,L,M,Z binary variable. &,J &,J & = 𝑓%&,?@A ∑M 𝑋%,L,M,Z ≤ 𝑉%,L,Z + ∑/∈=_ :L§+• 𝐹/,L,Z ≤ 𝑓%&,?BC ∑M 𝑋%,L,M,Z ∀𝑙, 𝑖, 𝑡 < |𝑇| w

(36)

If batch 𝑖 of product 𝑝 is involved in line 𝑙 during slot 𝑡 together with a higher-numbered batch 𝑖´ & (auxiliary binary variable 𝑋%,L´,Z = 1), then 𝑖´ is either assigned to a product 𝑝´ with whom 𝑝 forms a

non-forbidden sequence, i.e. 𝑝´ ∉ 𝑃M , or there are other batches in between with assigned products that 15 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 34

form a non-forbidden sequence with 𝑝. This logic proposition can be written as eq. (37), which can be &,J &,J & reformulated into eq. (38), after noting that ∑M´∉™¤ 𝑋%,L´,M´,Z + ∑M´∈™¤ 𝑋%,L´,M´,Z = 𝑋%,L´,Z . 42

&,J &,J &,J & 𝑋%,L,M,Z ∧ 𝑋%,L´,Z ⇒∨M´∉™¤ 𝑋%,L´,M´,Z ∨L«L´´«L´ ∨M´´∉™¤ 𝑋%,L´´,M´´,Z ∀𝑙, 𝑖´, 𝑖 < 𝑖´, 𝑝, 𝑡 < |𝑇|

(37)

&,J &,J &,J ∑L«L´´«L´ ∑M´´∉™¤ 𝑋%,L´´,M´´,Z ≥ 𝑋%,L,M,Z + ∑M´∈™¤ 𝑋%,L´,M´,Z − 1 ∀𝑙, 𝑖´, 𝑖 < 𝑖´, 𝑝, 𝑡 < |𝑇|

(38)

&,J If 𝑋%,L,M,Z =1, then batch 𝑖 must be assigned to product 𝑝. The logic proposition in eq. (39) can be

reformulated into eq. (40). One can then aggregate over 𝑙 and 𝑡 to reduce the number of constraints, leading to eq. (41). Preliminary results have shown that eq. (41) performs better that eq. (40).

5.

&,J 𝑋%,L,M,Z ⇒ 𝑌L,M ∀𝑙, 𝑖, 𝑝, 𝑡 < |𝑇|

(39)

&,J 𝑋%,L,M,Z ≤ 𝑌L,M ∀𝑙, 𝑖, 𝑝, 𝑡 < |𝑇|

(40)

&,J ∑% ∑Z«|6| 𝑋%,L,M,Z ≤ 𝑌L,M ⋅ |𝐿| ⋅ (|𝑇| − 1) ∀𝑖, 𝑝

(41)

COMPUTATIONAL RESULTS In this section, the proposed model is validated by using five benchmark instances from the

=,N literature . Data related to the initial inventory (𝑣/,M ) in the refinery and to the product demand at the 35, 37

=,UAV output nodes (𝑓/,M ) is given in Table 1 and Table 2, respectively. The incoming flowrate to the =,L/ refinery storage system (𝜌/,M ) is zero in EX1-EX4 and 2,000 m /day in EX5. 3

Table 1. Initial inventory in all nodes of benchmark instances. Instance EX1 EX2 EX3 EX4 EX5

Node/Product N1 N1 N1 N1 N1 N2 N3 N4 N6 N7 N9 N11 N13

P1 9,000 6,000 36,000 20,000 200,000 411 209 426 1,628 791 860 1,093 527

P2 8,000 32,000 20,000 200,000 23 31 39 729 140

P3 10,000 10,000 16,000 20,000 200,000 543 411 907 1,364 1,232

388 1,008

1,395 2,103

16

ACS Paragon Plus Environment

P4 5,000 18,000 42,000

P5 8,000 18,000

P6 5,000 18,000

Page 17 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Table 2. Product demand at the end of scheduling horizon. Instance EX1 EX2

EX3

EX4

EX5

Node/Product N2 N4 N6 N3 N5 N6 N8 N10 N11 N3 N4 N6 N8 N9 N11 N3 N5 N7 N10 N12 N2 N3 N4 N6 N7 N9 N11 N13

P1 3,000 3,000 2,000

P2 3,000 2,000

P3 3,000 1,000 2,000 2,000 4,000 1,000 3,500

3,000 1,000

2,000 3,000

1,000 2,000 5,200 4,800 6,300 5,300 6,700 8,800 2,400 4,600 800

2,000 2,000 3,900 3,400 13,200 4,700 2,400 2,900 2,100 4,000 900 1,000

8,215 4,185 8,625 32,550 15,810 17,205 21,855

465 620 775 14,570 2,790

1,800 10,850 8,215 18,135 27,280 24,645

7,750 20,150

42,060

P4 2,000 4,000

2,000

P5 2,000

P6

3,000

4,000

2,000 1,500 2,000

2,300 7,300 3,000 3,000 3,800 3,000 3,500 800

2,000 8,900 3,100 1,200

The model was implemented in GAMS 24.9.2 and solved by CPLEX 12.7.1 running in parallel deterministic mode using up to eight threads. The termination criteria were a relative optimality tolerance of 10 or a maximum wall time limit of 18,000 CPUs. The hardware consisted of a Windows −6

10, 64-bit desktop with an Intel i7-6700K (4.0 GHz) processor and 16 GB of RAM. In the search for the optimal solution, the tuning parameters of the model, the number of event points in the grid |𝑇| and the number of batches |𝐼|, are increased one by one until no better solution is found within 18,000 CPUs. The key performance indicators are given in Table 3 while the computational statistics are listed in Table 4. Table 4 reflects a very similar problem size between the current model and its predecessor in Ref , with the latter typically requiring more equations. Compared to Ref , and for Ex5, our model 35

37

roughly triples the number of binary variables for a similar number of variables and constraints. This 17

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 34

is because our model has the additional feature of being able to generate the detailed delivery schedule, as explained in section 5.5. Table 3. Comparison between the results of this work and the literature. Instance

Source

EX1

Ref , EX1

LP relaxation 135 h

EX2

Ref , EX2

145 h

EX3

Ref , EX3

EX4

EX5

35

35

35

Ref , Motiv. 35

Ref , EX3 37

125.25 h

151.83 h

10.81 days

|𝑇| 5 5 6 6 7 8 9 7 7 8 8 9 10 6 6 7 7 8 8 9 10 6 6 7 7 8 9 10 11 12 6 6 7 7 8 8 9

This work CPUs |𝐼| 7 9.7 8 60.1 7 52.5 8 660 7 33.8 7 34.2 7 94.2 7 51.9 8 829 7 3,319 8 11,452 7 8,783 7 18,000 7 18.7 8 211 7 366 8 2,505 7 2,570 8 18,000 7 18,000 7 140 6 16.7 7 158 6 120 7 1,020 6 316 6 5,931 6 18,000 6 18,000 6 18,000 12 67.4 13 587 12 323 13 11,649 12 675 13 18,000 12 1,345

Makespan Infeasible 195 h 140 h 140 h 135 h 135 h 135 h Infeasible 162.86 h 150.71 h 150.71 h 150.71 h 150.71 h Infeasible Infeasible 150.75 h 147.93 h 127.86 h 127.86 h 125.80 h 125.25 h Infeasible Infeasible 184.69 h 177.30 h 164.45 h 160.95 h 157.44 h 156.67 h 156.67 h Infeasible Infeasible 13.26 d 13.18 d 12.51 d 12.72 d 12.51 d

Makespan 143.33 h 140 h 138.33 h 136.67 h 136.67 h

18

|𝐼| 7 7 7 7 7

170 h 160 h 155.71 h 151.43 h

287 13,613 18,000 18,000

10 11 12 13

7 7 7 7

136.25 h 129.09 h 125.25 h

888 3,392 4,392

11 12 13

7 7 7

190.8 h 170.95 h 166.40 h 178.45 h

4,944 18,000 18,000 18,000

12 13 14 15

6 6 6 6

13.11d

65.08

5

12

a

b c

d e f

g

aGap=3.79%, bGap=2.08%, cGap=0.44%, dGap=2.87%, eGap=0.86%, fGap=3.05%, gGap=6.91%.

ACS Paragon Plus Environment

Literature CPUs |𝑇| 386 9 933 10 10,940 11 18,000 12 18,000 13

Page 19 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Table 4. Statistics related to problem size. Instance

Source

|T|

|𝐼|

Ref , EX1

7 8 9 8 9 10 10 11 11 12 5 6 7 8

7 7 7 7 7 7 7 7 6 6 12 12 12 12

35

EX1 Ref , EX2 35

EX2 Ref , EX3 35

EX3 EX4

Ref , Motiv.

EX5

Ref , EX3

35

37

Binary variables 1,058 1,212 1,368 1,915 2,166 2,417 2,065 2,282 1,872 2,050 1,422 1,752 2,082 2,412

This work Total Equations variables 4,735 4,215 5,062 4,773 5,389 5,331 6,835 7,432 7,369 8,348 7,903 9,264 6,739 7,242 7,238 7,961 5,683 6,991 6,070 7,618 3,183 5,598 3,917 6,851 4,651 8,104 5,385 9,357

Binary variables

Literature Total Equations variables

1,224

3,478

6,108

2,346

6,803

12,103

2,602

6,731

11,124

2,451 515

6,070 3,169

10,748 5,123

5.1. Example 1 In EX1, six products (P1-P6) are conveyed from refinery N1 to three depots, among which N2 and N4 are located in mainline L1 and N6 is located in branched line L2. The branch line L2 is linked to mainline L1 by junction J1 (i.e. N3 and N5). The initial pipeline status, depot locations and network topology can be found in the first line of Figure 7. The allowable operating rate ranges at nodes, given in (m /h) are: [80, 200] for refinery N1 and depot N2, [20, 150] for depot N4, and [20, 100] for depot 3

N6. While the available flowrate ranges at segments S1-S4 are [80, 200], [40, 200], [20, 150], [20, 100] m /h, respectively. The forbidden product sequences for EX1 are P1-P5, P1-P6, P2-P4, P2-P5 3

and vice-versa.

19

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 34

Figure 7. Optimal solution for EX1 (|𝑻|=7 and |𝑰|=7). Figure 7 displays the best-found schedule for |𝑇|=7 and |𝐼|=7, which we know is the optimal solution since it is equal to the LP relaxation (check Table 3). This is because it is possible to operate the refinery pumps at full capacity during the scheduling horizon, i.e. at a maximum rate of 200 m /h. 3

Consequently, the minimum time needed to complete the transportation tasks drops to 135 from 136.67 h (Ref ). 35

In the schedule of Figure 7, four new batches (I4-I7) are being injected. Notice that intervals [5, 25], [25, 55], [80, 115] and [115, 135] h involve more than one batch being pumped from N1. Depot N2 receives at most one batch from the pipeline during a time slot, but depots N4 and N6 sometimes receive more than one. This is because N4 and N6 are terminal nodes of lines L1 and L2 (there are no adjacent downstream segments) and so there is no restriction on the number of batches being discharged to them. More specifically, N4 first receives multiple batches in interval [25, 55], in a total of 3,000 m . Product P2 is the first to arrive (500 m ), since it corresponds to the right-most batch (I1) 3

3

20

ACS Paragon Plus Environment

Page 21 of 34

Industrial & Engineering Chemistry Research

in segment S3 at the 25-h mark. Considering the operating flowrate of 100 m /h after depot N2, we 3

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

can compute that depot N4 receives 1,500 m of P1 (I2) from [30, 45] and finally 1,000 m of P4 (I3) 3

3

from [45, 55] h. The inclusion of multiple batch injections/deliveries during a single time slot, is the major advantage of our model with respect to Ref , and the reason why we can find a better solution 35

within a much shorter time. Looking at junction J1 during intervals [0, 5], [55, 80], [80, 115], [115, 135] h, one can see that there is only one batch crossing it since node N3 is simultaneously feeding segments S3 and S4. 5.2. Example 2 The branched pipeline network studied in EX2 involves three lines (L1-L3) distributing six products (P1-P6) from refinery N1 to six depots. There is a total of three depots (N3, N5 and N6) along the mainline L1. The branch line L2 is linked to mainline L1 by junction J1 (i.e. N2 and N7) and transports products to its terminal depot N8. While the branch line L3 is connected to mainline L1 by junction J2 (i.e. N4 and N9) and transports products to two depots (N10 and N11), with one (N10) being located midway over line L3. The allowable operating rates at nodes, given in (m /h), are as follows: [100, 3

350] for refinery N1 and depot N3, [50, 200] for depot N5, [25, 200] for depot N6, [50, 150] for depot N8, [50, 200] for depot N10, and [25, 100] for depot N11. Moreover, the admissible flowrate ranges at segments S1-S8 are [100, 350], [100, 350], [100, 200], [50, 200], [25, 200], [50, 150], [50, 200], and [25, 100] m /h, respectively. EX2 features as forbidden sequences P1-P3, P2-P5, P3-P5, P4-P6, 3

and vice-versa. Compared with EX1, EX2 features a more complex pipeline system with roughly twice the number of nodes and segments, leading to a significant increase in problem size (see Table 4). From Table 3, we can see that the best-found solution is obtained for |𝑇|=8 and |𝐼|=7. Figure 8 depicts the best-found schedule, where it can be observed that the makespan drops to 150.71 h, compared to the 151.43 h reported in Ref . 35

Three old batches exist initially inside the pipeline, with node N1 beginning operation with the injection of 6000 m of I3 (P3) and then 1000 m of I4 (P4). The remaining injections correspond to 3

3

four new batches (I4-I7), containing a P4-P5-P1-P6 sequence. Note that three intervals, [0, 20], [20, 21

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 34

45] and [65, 105] h, involve the injection of more than one product. Besides N1, there are other occurrences where more than one batch is processed by a node during a slot, such as two batches being dispatched to nodes N5, N6, N8, N10 and N11 during [20, 45], [0, 20] and [65, 105], [130.71, 150.71], [130.71, 150.71], and [65, 105] h, respectively. Nodes N6, N8 and N11 are the terminal nodes of lines L1-L3. Node N5 receives batches I2 and I3 when its downstream segment S5 is idle, while node N10 receives batches I5 and I6 when segment S8 is idle.

Figure 8. Best-found solution for EX2 (|𝑻|=8 and |𝑰|=7). The makespan is higher than the LP relaxation of 145 h in Table 3, which is not surprising considering that the injection flowrate is lower than the maximum of 350 m /h in slots [20, 45] and 3

[65, 105] h. As an example, the reason why slot [20, 45] can only reach 200 m /h, the maximum 3

22

ACS Paragon Plus Environment

Page 23 of 34

Industrial & Engineering Chemistry Research

flowrate in segment S3, can be explained by the node demand data. Since N8 received 1000 m of P3 3

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

during [0, 20] h, its residual demand is 2500 m (check values in Table 2), the volume stored inside 3

segment S6. Thus, no more P3 should be injected to line L2. Similarly, node N3 received the full 2000 m demand of P3 during the first twenty hours of operation and so, all the flow coming out of segment 3

S1 should go to segment S3. 5.3. Example 3 In EX3, four products (P1-P4) are discharged to six depots located along the mainline (L1) and two branch lines (L2 and L3). The network topology of EX3 is similar with that of EX2, except that line L2 holds two depots while line L3 only holds one depot. The allowable operating rate ranges at nodes, given in (m /h), are as follows: [300, 800] for refinery N1, [200, 800] for depot N3, [150, 600] for 3

depot N4, [100, 500] for depot N6, [150, 600] for depot N8, [100, 600] for depot N9, and [100, 300] for depot N11. While the available flowrate ranges at the segments S1-S8 are [300, 800], [200, 800], [150, 600], [150, 500], [100, 500], [150, 600], [100, 600], and [100, 300] m /h, respectively. Product 3

sequences P1-P4 and P4-P1 are forbidden to avoid batch contamination. The first feasible solution appears when selecting 7 event points in the grid and 3 new batches to be injected. It corresponds to a makespan of 150.75 h (check Table 3) that keeps decreasing until finding the global optimal solution of 125.25 for |𝑇|=10 and |𝐼|=7. Note that due to the zero integrality gap, it is much faster to solve such problem to optimality (140.1 CPUs) than the MILP for |𝑇|=9, which reached the specified time limit. The latter is a disadvantage compared to the model by Castro and Mostafaei that was able to reach the global optimal solution for |𝑇|=13 in 4,392 CPUs, without any 35

previous iteration hitting the maximum time limit.

23

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 34

Figure 9. Optimal solution for EX3 (|𝑻|=10 and |𝑰|=7). While featuring the same makespan, the schedule depicted in Figure 9 is not the same as the one reported in Ref . In Ref , there are a total of nine occurrences with multiple batches crossing segments 35

35

during a single interval, an improvement with respect to their previous product-centric model that 1

allowed at most one batch crossing. The formulation proposed in this paper can be seen as another 24

ACS Paragon Plus Environment

Page 25 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

major step in the direction of improved computational performance, since it now also allows multiple batches to be processed by a node during a single interval. Specifically, there are a total of thirteen occurrences, such as batches I3-I4 being allocated to N8 during interval [0, 17.13] h, batches I4-I5 being pumped by N1 during interval [17.13, 40.26] h, etc. As a consequence, we were able to reduce the number of event points required to generate the optimal schedule from 13 to 10, while decreasing the computational time by 96.8%. 5.4. Example 4 EX4 comes from the motivating example in Section 2.1 and features multilevel branching. The allowable operating ranges at the nodes, given in (m /h), are as follows: [90, 164] for refinery N1 and 3

depot N3, [70, 120] for depot N5, [40, 50] for depot N7, [40, 80] for depots N10 and N12. While the available flowrate ranges at the segments S1-S8 are [90, 164], [90, 164], [70, 120], [70, 120], [40, 50], [40, 80], [40, 80], and [40, 80] m /h, respectively. P1-P2 and P2-P1 are forbidden product sequences. 3

By checking Table 3, we can see that the first feasible solution is obtained for |𝑇|=7. Better solutions are then found after raising |𝑇|. However, no further improvements are observed within 18,000 CPUs when |𝑇| is greater than 11. Compared to the best solution found by the model in Ref (166.40 h, for 35

|𝑇|=14), our model found a much better solution up to the maximum time limit (156.67 h, for |𝑇|=11 and |𝐼|=6), as shown in Figure 10. Nevertheless, it should be noted that the computational time for |𝑇|=12 and Ref is well below the maximum limit observed for the current model. However, the 35

current model can obtain much shorter makespans than the 190.8 h from Ref in just 120 CPUs (184.69 35

h, for |𝑇|=7) and 316 CPUs (164.45 h, for |𝑇|=8). In Figure 10, a total of three new batches are pumped from N1 over ten slots, with two slots, [0, 21.83] and [51.34, 64.33] h, involving two batches being injected into the pipeline successively. Apart from N1, there are other nine occurrences of more than one batch being processed by a node over a slot, i.e. two batches crossing junctions J2 and J3 during interval [94.05, 107.19] h, and multiple batches being distributed to terminal nodes N5, N7, N10 and N12 during [0, 21.83], [49.05, 51.34] and [107.19, 132.19], [73.6, 94.05] and [107.19, 132.19], [107.19, 132.19], and [94.05, 107.19] h, respectively. 25

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 34

Interestingly, the makespan is still 4.84 h higher than its LP relaxation in Table 3 although the pumping rate is always kept at the maximum rate of 164 m /h. This is caused by excess P1 initially 3

inside the system. Checking Table 2 and the first line of Figure 10, the total demand for P1 is 7,800 m while 8,593 m of P1 exists inside lines L1 and L2 at the beginning. Hence, 4.84 h is required to 3

3

push the surplus P1 out of the pipeline at a rate of 164 m /h before distributing the new batches to 3

depots.

Figure 10. Optimal solution for EX4 (|𝑻|=11 and |𝑰|=6). 26

ACS Paragon Plus Environment

Page 27 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

5.5. Example 5 The last example deals with a more complicated real-world case , which conveys three products P137

P3 from refinery N1 to eight depots through a mainline (L1) and two branch lines (L2 and L3). The injecting flowrate at the refinery is fixed at 25,200 m /d, and there are no limits for delivery and 3

segment flowrates. During any slot, the minimum volume being transported by a segment is 1,000 m

3

and that being dispatched to a depot is 100 m . Products P1 and P3 cannot be transported adjacently. 3

It can be seen from Table 3 that the best makespan of 12.51 days is obtained for |𝑇|=8 and |𝐼|=12. While Ref considers another objective function, the minimization of total cost, the reported makespan 37

of 13.1 days is 4.7% higher. Note that the model by Cafaro and Cerdá can provide the optimal 37

sequence of batch injections but cannot determine the optimal sequence of product deliveries to depots. It is thus unable to generate a detailed delivery schedule, unlike our model. This is apparent from figure 12 in Ref , where product P2 is delivered to N3 during the first hour (first slot), even though its 37

right coordinate is still distant from the node location at the beginning of the time horizon (check first row of Figure 11). In contrast, our model waits until the second slot to make the required 589 m

3

delivery (620 m of product demand minus the initial inventory of 31 m ; check values in Table 1 and 3

3

Table 2), after all of batch I7 has surpassed N3. Given the ability to generate more detailed schedules, it is not surprising that our model requires more slots (8 vs. 5) and longer computational times (675 vs. 65.08 CPUs). The reduction in makespan is due to the delivery of smaller volumes of product P2 to the depots. From the values in Table 1 and Table 2 one knows that the minimum volumes, given in (m ), are the 3

following: 442 for depot N2, 589 for N3, 736 for N4, 13,841 for N6, 2,650 for N7, 7,362 for N11 and 19,142 for depot N13. Our schedule in Figure 11 delivers the minimum volumes for nodes and N2N4, for a total of 80821 m in all nodes vs. 95579 m for the schedule in Ref . 3

3

27

ACS Paragon Plus Environment

37

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 34

Figure 11. Optimal solution for EX5 (|𝑻|=8 and |𝑰|=12). 6.

CONCLUSIONS In this paper, a novel MILP continuous-time, batch-centric formulation has been developed. The

proposed model can be applied to a branched pipeline system with a single refinery and multiple depots, rigorously preventing forbidden product sequences from being generated in each line at any time. Compared with previous batch-centric models, our model allows multiple batches to be processed by a node over a composite run. Hence, fewer event points are required to discover a highquality solution, resulting in a noticeable improvement in computational time. The superiority of our model has been validated by successfully solving five benchmark problems from the literature . By enabling a lower number of event points, the same or a better schedule is 35, 37

28

ACS Paragon Plus Environment

Page 29 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

found in much less time. More specifically, the solutions returned for EX1 and EX2 by the recent model in Ref at the 18,000 CPUs termination are worse than the ones found by the current model in 35

33.8 CPUs and 3,319 CPUs, respectively. EX4 hits the maximum time limit for both models but the makespan from our best-found schedule is 9.73 h lower than that of Ref . For EX5 from Ref , our 35

37

model can obtain a detailed schedule with a shorter makespan, at the expense of using more computational resources. In all cases, forbidden sequences were avoided despite the occurrence of empty batches in the schedule, which is not the case of most previous models. Overall, the proposed model is highly modular and extensible, and the computational time increases quite reasonably as the complexity grows. It will be generalized in the future, to handle more complicated pipeline configurations, such as mesh-structures with multiple refineries and depots. ACKNOWLEGMENT Financial support from the National Natural Science Foundation of China, Grant No. 51874325, and from

Fundação

para

a

Ciência

e

Tecnologia

through

UID/MAT/04561/2019. NOMENCLATURE Sets/Indices 𝐼/𝑖, 𝑖´ = batches 𝐼 I%J = old batches 𝐼% = batches that can appear in line 𝑙 𝐼/ = batches that can appear in node 𝑛 𝐽/𝑗 = junction 𝐿/𝑙 = lines in the pipeline network 𝑁/𝑛, 𝑛´ = nodes 𝑁 7 = depot node 𝑁8+ = input node of junction 𝑗 𝑁%+ = input node of line 𝑙 𝑁8- = output node of junction 𝑗 𝑁% = Nodes of line 𝑙 𝑁%- = output nodes of line 𝑙 𝑁 0 = refinery node 𝑁%6 = terminal node of line 𝑙 𝑃/𝑝, 𝑝´ = products 𝑃/ = products that can appear in node 𝑛 𝑃M = products that cannot enter a segment after 𝑝 29

ACS Paragon Plus Environment

projects

IF/00781/2013

and

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 34

𝑆/𝑠 = pipeline segments 𝑆/7 = downstream segment to node 𝑛 𝑆/W = upstream segment to node 𝑛 𝑇/𝑡 = event points of the time grid Parameters 𝑓%&,?@A = maximum volume of a batch transported in line 𝑙 during a time slot (m ) 𝑓%&,?BC = minimum volume of a batch transported in line 𝑙 during a time slot (m ) =,UAV 𝑓/,M = demand at the end of scheduling horizon for node 𝑛 and product 𝑝 (m ) =,?BC 𝑓/ = maximum volume transferred from/to node 𝑛 during a time slot (m ) =,?@A 𝑓/ = minimum volume transferred from/to node 𝑛 during a time slot (m ) ™,?BC 𝑓M = maximum volume of product 𝑝 that can be transferred during the time horizon (m ) [,?BC 𝑓e = maximum volume into segment 𝑠 during a time slot (m ) [,?@A 𝑓e = minimum volume into segment 𝑠 during a time slot (m ) ℎ = time horizon (h) &,N 𝑙𝑐%,L = left coordinate of batch 𝑖 in line 𝑙 at the start of the time horizon (m ) &,N 𝑟𝑐%,L = right coordinate of batch 𝑖 in line 𝑙 at the start of the time horizon (m ) =,N 𝑣/,M = initial volume of product 𝑝 in storage system of refinery or depot node 𝑛 (m ) =,?BC 𝑣/,M = maximum volume of product 𝑝 in storage system of node 𝑛 (m ) =,?@A 𝑣/,M = minimum volume of product 𝑝 in storage system of node 𝑛 (m ) & 𝑣% = volume of line 𝑙 (m ) &,N 𝑣%,L = initial volume of batch 𝑖 in line 𝑙 (m ) 𝑦L,M = equal to one when old batch 𝑖 deals with product 𝑝 𝜎/ = volumetric coordinate of node 𝑛 (m ) =,L/ 𝜌/,M = incoming flowrate of product 𝑝 to node 𝑛 (m /h) =,?BC 𝜌/ = maximum flowrate that can be transferred to/from node 𝑛 (m /h) =,?@A 𝜌/ = minimum flowrate that can be transferred to/from node 𝑛 (m /h) [,?BC 𝜌e = maximum flowrate in segment 𝑠 (m /h) [,?@A 𝜌e = minimum flowrate in segment 𝑠 (m /h) 3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

3

Variables & 𝑋%,L,Z = binary variable indicating batch 𝑖 is transported in line 𝑙 during slot 𝑡 &,J 𝑋%,L,M,Z = binary variable indicating batch 𝑖 containing product 𝑝 is transported in line 𝑙 during slot 𝑡

=,@V³U 𝑋/,Z = binary variable indicating node 𝑛 is not processing batches during slot 𝑡 = 𝑋/,L,Z = binary variable indicating node 𝑛 is processing batch 𝑖 during slot 𝑡 =,j 𝑋/,L,Z = binary variable indicating batch 𝑖 is the first to arrive (during slot 𝑡) to node 𝑛 =,Y 𝑋/,L,Z = binary variable indicating batch 𝑖 is not the first to arrive (during slot 𝑡) to node 𝑛 [ 𝑋e,Z = binary variable indicating segment 𝑠 is active during slot 𝑡 𝑌L,M = binary variable indicating batch 𝑖 involves product 𝑝 𝐿𝐶%,L,Z = left coordinate of batch 𝑖 in line 𝑙 at event point 𝑡 (m ) 𝐷Z = duration of time slot 𝑡 (h) = 𝐹/,L,Z = volume of batch 𝑖 transferred to/from node 𝑛 during slot 𝑡 (m ) =,J 𝐹/,L,M,Z = volume of batch 𝑖 of product 𝑝 transferred to/from refinery or depot node 𝑛 during slot 𝑡 (m ) [ 𝐹e,Z = volume entering segment 𝑠 during slot 𝑡 (m ) 𝑅𝐶%,L,Z = right coordinate of batch 𝑖 in line 𝑙 at event point 𝑡 (m ) 𝑇Z = time of event point 𝑡 (h) 3

3

3

3

3

30

ACS Paragon Plus Environment

Page 31 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

= 𝑉/,M,Z = volume of product 𝑝 in storage system of node 𝑛 at event point 𝑡 (m ) & 𝑉%,L,Z = volume in line 𝑙 of batch 𝑖 at event point 𝑡 (m ) 3

3

References 1. Castro, P. M.; Mostafaei, H., Product-centric continuous-time formulation for pipeline scheduling. Computers & Chemical Engineering 2017, 104, 283-295. 2. Castro, P. M., Optimal Scheduling of Multiproduct Pipelines in Networks with Reversible Flow. Industrial & Engineering Chemistry Research 2017, 56, (34), 9638-9656. 3. Cafaro, D. C.; Cerdá, J., Rigorous scheduling of mesh-structure refined petroleum pipeline networks. Computers & Chemical Engineering 2012, 38, 185-203. 4. Cafaro, D. C.; Cerdá, J., Rigorous formulation for the scheduling of reversible-flow multiproduct pipelines. Computers & Chemical Engineering 2014, 61, 59-76. 5. Liao, Q.; Liang, Y.; Xu, N.; Zhang, H.; Wang, J.; Zhou, X., An MILP approach for detailed scheduling of multi-product pipeline in pressure control mode. Chemical Engineering Research and Design 2018, 136, 620-637. 6. Zhou, X.; Zhang, H.; Qiu, R.; Liang, Y.; Wu, G.; Xiang, C.; Yan, X., A hybrid time MILP model for the pump scheduling of multi-product pipelines based on the rigorous description of the pipeline hydraulic loss changes. Computers & Chemical Engineering 2019, 121, 174-199. 7. Chen, H.; Zuo, L.; Wu, C.; Wang, L.; Diao, F.; Chen, J.; Huang, Y., Optimizing detailed schedules of a multiproduct pipeline by a monolithic MILP formulation. Journal of Petroleum Science and Engineering 2017, 159, 148-163. 8. Zhang, H.; Liang, Y.; Liao, Q.; Shen, Y.; Yan, X., A self-learning approach for optimal detailed scheduling of multi-product pipeline. Journal of Computational and Applied Mathematics 2018, 327, 41-63. 9. Harjunkoski, I.; Maravelias, C. T.; Bongers, P.; Castro, P. M.; Engell, S.; Grossmann, I. E.; Hooker, J.; Méndez, C.; Sand, G.; Wassick, J., Scope for industrial applications of production scheduling models and solution methods. Computers & Chemical Engineering 2014, 62, 161-193. 10. Zhang, H.; Liang, Y.; Xiao, Q.; Wu, M.; Shao, Q., Supply-based optimal scheduling of oil product pipelines. Petroleum Science 2016, 13, (2), 355-367. 11. Zhang, H.; Liang, Y.; Liao, Q.; Ma, J.; Yan, X., An MILP approach for detailed scheduling of oil depots along a multi-product pipeline. Petroleum Science 2017, 14, (2), 434-458. 12. Relvas, S.; Matos, H. A.; Barbosa-Póvoa, A. P. F. D.; Fialho, J.; Pinheiro, A. S., Pipeline Scheduling and Inventory Management of a Multiproduct Distribution Oil System. Industrial & Engineering Chemistry Research 2006, 45, (23), 7841-7855. 13. Mostafaei, H.; Ghaffari Hadigheh, A., A General Modeling Framework for the Long-Term Scheduling of Multiproduct Pipelines with Delivery Constraints. Industrial & Engineering Chemistry Research 2014, 53, (17), 7029-7042. 14. Liao, Q.; Zhang, H.; Wang, Y.; Zhang, W.; Liang, Y., Heuristic method for detailed scheduling of branched multiproduct pipeline networks. Chemical Engineering Research and Design 2018, 140, 82-101. 15. Cafaro, V. G.; Cafaro, D. C.; Méndez, C. A.; Cerdá, J., Detailed Scheduling of Operations in Single-Source Refined Products Pipelines. Industrial & Engineering Chemistry Research 2011, 50, (10), 6240-6259. 16. Cafaro, V. G.; Cafaro, D. C.; Méndez, C. A.; Cerdá, J., Detailed Scheduling of Single-Source Pipelines with Simultaneous Deliveries to Multiple Offtake Stations. Industrial & Engineering Chemistry Research 2012, 51, (17), 6145-6165. 17. Cafaro, D. C.; Cerdá, J., Optimal scheduling of multiproduct pipeline systems using a nondiscrete MILP formulation. Computers & Chemical Engineering 2004, 28, (10), 2053-2068. 18. Cafaro, D. C.; Cerdá, J., Optimal Scheduling of Refined Products Pipelines with Multiple Sources. Industrial & Engineering Chemistry Research 2009, 48, (14), 6675-6689. 19. Castro, P. M.; Barbosa-Póvoa, A. P.; Matos, H. A.; Novais, A. Q., Simple Continuous-Time Formulation for Short-Term Scheduling of Batch and Continuous Processes. Industrial & Engineering Chemistry Research 2004, 43, (1), 105-118. 31

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 34

20. Castro, P. M., Optimal Scheduling of Pipeline Systems with a Resource−Task Network Continuous-Time Formulation. Industrial & Engineering Chemistry Research 2010, 49, (22), 1149111505. 21. Relvas, S.; Boschetto Magatão, S. N.; Barbosa-Póvoa, A. P. F. D.; Neves, F., Integrated scheduling and inventory management of an oil products distribution system. Omega 2013, 41, (6), 955-968. 22. Relvas, S.; Barbosa-Póvoa, A. P. F. D.; Matos, H. A.; Fialho, J.; Pinheiro, A. S., Pipeline scheduling and distribution centre management—A real-world scenario at CLC. In Computer Aided Chemical Engineering, 2006; Vol. 21, pp 2135-2140. 23. MirHassani, S. A.; Abbasi, M.; Moradi, S., Operational scheduling of refined product pipeline with dual purpose depots. Applied Mathematical Modelling 2013, 37, (8), 5723-5742. 24. MirHassani, S. A.; Fani Jahromi, H., Scheduling multi-product tree-structure pipelines. Computers & Chemical Engineering 2011, 35, (1), 165-176. 25. Rejowski, R.; Pinto, J. M., A novel continuous time representation for the scheduling of pipeline systems with pumping yield rate constraints. Computers & Chemical Engineering 2008, 32, (4), 1042-1066. 26. Rejowski, R.; Pinto, J. M., A rigorous MINLP for the simultaneous scheduling and operation of multiproduct pipeline systems. In Computer Aided Chemical Engineering, 2005; Vol. 20, pp 10631068. 27. Cafaro, V. G.; Cafaro, D. C.; Méndez, C. A.; Cerdá, J., MINLP model for the detailed scheduling of refined products pipelines with flow rate dependent pumping costs. Computers & Chemical Engineering 2015, 72, 210-221. 28. Zhang, H.; Liang, Y.; Liao, Q.; Wu, M.; Yan, X., A hybrid computational approach for detailed scheduling of products in a pipeline with multiple pump stations. Energy 2017, 119, 612-628. 29. Viswanathan, J.; Grossmann, I. E., A combined penalty function and outer-approximation method for MINLP optimization. Computers & Chemical Engineering 1990, 14, (7), 769-782. 30. Relvas, S.; Barbosa-Póvoa, A. P. F. D.; Matos, H. A., Heuristic batch sequencing on a multiproduct oil distribution system. Computers & Chemical Engineering 2009, 33, (3), 712-730. 31. MirHassani, S. A.; BeheshtiAsl, N., A heuristic batch sequencing for multiproduct pipelines. Computers & Chemical Engineering 2013, 56, 58-67. 32. Polli, H. L.; Magatão, L.; Magatão, S. N. B.; Neves, F.; Arruda, L. V. R., Collaborative Approach Based on Heuristic Algorithm and MILP Model To Assignment and Sequencing of Oil Derivative Batches in Pipeline Networks. Industrial & Engineering Chemistry Research 2017, 56, (9), 2492-2514. 33. Mostafaei, H.; Castro, P. M.; Ghaffari-Hadigheh, A., A Novel Monolithic MILP Framework for Lot-Sizing and Scheduling of Multiproduct Treelike Pipeline Networks. Industrial & Engineering Chemistry Research 2015, 54, (37), 9202-9221. 34. Mostafaei, H.; Castro, P. M., Continuous-time scheduling formulation for straight pipelines. AIChE J 2017, 63, 1923-1936. 35. Castro, P. M.; Mostafaei, H., Batch-centric scheduling formulation for treelike pipeline systems with forbidden product sequences. Computers & Chemical Engineering 2018, https://doi.org/10.1016/j.compchemeng.2018.04.027. 36. Liao, Q.; Zhang, H.; Xu, N.; Liang, Y.; Wang, J., A MILP model based on flowrate database for detailed scheduling of a multi-product pipeline with multiple pump stations. Computers & Chemical Engineering 2018, 117, 63-81. 37. Cafaro, D. C.; Cerdá, J., A Rigorous Mathematical Formulation for the Scheduling of TreeStructure Pipeline Networks. Industrial & Engineering Chemistry Research 2011, 50, (9), 5064-5085. 38. Raman, R.; Grossmann, I. E., Modelling and computational techniques for logic based integer programming. Computers & Chemical Engineering 1994, 18, (7), 563-578. 39. Balas, E., Disjunctive Programming. Annals of Discrete Mathematics 1979, 5, (2), 3-51. 40. Balas, E., Disjunctive Programming and a Hierarchy of Relaxations for Discrete Optimization Problems. SIAM Journal on Algebraic Discrete Methods 1985, 6, (3), 466-486. 41. Vecchietti, A.; Grossmann, I. E., Modeling issues and implementation of language for disjunctive programming. Computers & Chemical Engineering 2000, 24, (9), 2143-2155. 32

ACS Paragon Plus Environment

Page 33 of 34 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

42. Castro, P. M.; Grossmann, I. E., Generalized Disjunctive Programming as a Systematic Modeling Framework to Derive Scheduling Formulations. Industrial & Engineering Chemistry Research 2012, 51, (16), 5781-5792. 43. Castro, P. M.; Mostafaei, H., New continuous-time scheduling formulation for multilevel treelike pipeline systems. In Computer Aided Chemical Engineering, Friedl, A.; Klemeš, J. J.; Radl, S.; Varbanov, P. S.; Wallek, T., Eds. Elsevier: 2018; Vol. 43, pp 973-978.

33

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

Table of Contents Graphic

I3 ! %∈./

I2

I2

Input node I2

Terminal output node

I1

!

I1

%∈./

I1

' 𝑋#,%,& 𝐿𝐶*,%,& = 0

' 𝑋#,%,& 𝑅𝐶*,%,&12 = 𝑣*4

I1

I2 I2

Intermediate output node I2

I3

I2 I3

',2 𝑋#,%,&

I3

I2

',5 𝑋#,%,&

I2

I2 I3 9 Active downstream segment (𝑋8,& =1)

',2 𝑋#,%,& I2 I3 9 =0) Idle downstream segment (𝑋8,&

Processing of multiple batches

',5 9 𝑋#,%,& ⇒ ¬𝑋8,&

34

ACS Paragon Plus Environment

I3 I2

I3

I2

I4

I2 I1

I2

I4 I3

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 34