ARTICLE pubs.acs.org/IECR
Scheduling Gasoline Blending Operations from Recipe Determination to Shipping Using Unit Slots Jie Li and I. A. Karimi* Department of Chemical and Biomolecular Engineering, National University of Singapore, 4 Engineering Drive 4, Singapore 117576 ABSTRACT: Gasoline can yield 60%70% of the total profit of a refinery. Integrated scheduling of gasoline blending and order delivery operations can increase profit by avoiding ship demurrage, improving customer satisfaction, minimizing quality giveaways, avoiding costly transitions and slop generation, and reducing inventory costs. In this paper, we develop a novel continuous-time formulation using unit slots for an integrated scheduling of gasoline blending operations incorporating features such as multipurpose tanks, parallel nonidentical blenders, minimum run lengths, blender setups, tank changeovers, piecewise-constant profiles for blend component qualities and feed rates, simultaneous receipt/delivery by tanks, etc. We employ and revise the schedule adjustment procedure of Li et al. [AIChE J. 2010, 56, 441465] to avoid solving mixed-integer nonlinear programming (MINLP) problems and illustrate the effectiveness of our new formulation on 14 industry-scale examples.
1. INTRODUCTION The petroleum refining industry constantly faces challenges such as fluctuating product demands, ever-changing crude prices, strict environmental regulations, etc. Optimization can play (and has played) a key role in managing oil refineries, as the refiners seek to exploit all potential cost-saving opportunities to survive and thrive in a low-margin business. Oil refineries have used optimization (especially successive linear programming or SLP) for planning refinery operations for a long time. Refinery operations1 involve three main segments, namely, crude oil storage and processing, intermediate processing, and product blending and distribution. While the first involves upstream operations from ship arrivals to crude distillation,2 the last involves downstream operations from the blending of hydrocarbon fractions to the storage and delivery of final products. Gasoline is one of the most profitable products of a refinery and can account for as much as 60%70% of total profit. A refinery typically blends several gasoline cuts or fractions from various processes to meet its customer orders of varying specifications. However, the large numbers of orders, delivery dates, blenders, blend components, tanks, quality specifications, etc. and nonlinear blending make this problem highly complex and nonlinear. If done manually, it can result in suboptimal schedules and costly quality giveaways. Therefore, the scheduling of gasoline blending and delivery is crucial. Scheduling such operations using advanced optimization techniques such as mixed-integer programming (MIP) can help avoid ship demurrage, improve order delivery and customer satisfaction, minimize quality giveaways, reduce transitions and slop generation, exploit low-quality cuts, and reduce inventory costs. The early work related to this problem focused on the optimal blending of various intermediate fractions and some additives to meet product quality specifications.3,4 Glismann and Gruhn5 developed a two-level decomposition approach to integrate short-term scheduling with blend-recipe optimization. They first solve a nonlinear programming (NLP) problem to obtain blend recipes and product quantities, and then they use mixed-integer r 2011 American Chemical Society
linear programming (MILP) to obtain a schedule for the blending operations. Pinto et al.1 and Joly and Pinto6 developed MIP models for scheduling blending operations involving fuel oil and asphalt. They used linear blending correlations for key elements, but they did not address nonlinear product properties such as sulfur. While they considered transitions in blending pipelines, they did not ensure constant blending rates, nor did they enforce minimum run lengths for blend runs. Jia and Ierapetritou7 proposed a continuous-time event-based MILP formulation for scheduling gasoline blending and order delivery operations simultaneously. While their model allows features such as multipurpose product tanks, one product tank delivering multiple orders, and multiple product tanks delivering one order, it lacks some key features, such as multiple parallel nonidentical blenders, variable recipes, and product specifications. However, as shown by Li et al.,8 it gives infeasible solutions and allows a product tank to hold several products at a time. Mendez et al.9 presented both discrete-time and slot-based continuous-time models for the simultaneous optimization of blending and short-term scheduling. Their model allows parallel identical blenders and determines optimal blend recipes. However, (i) it does not ensure constant blend rates or minimum run lengths; (ii) it does not allow multipurpose product tanks or setup times for blender switching; or (iii) it does not integrate order delivery time windows with blending operations. Furthermore, it uses nonlinear blending correlations for some specifications. Recently, Li et al.8 reported a slot-based continuous-time formulation using process slots10 for the simultaneous treatment of recipe, blending, storage, and order scheduling. They incorporated several operational features such as multiple-purpose product tanks, parallel nonidentical blenders, one blender Received: November 16, 2010 Accepted: May 20, 2011 Revised: April 24, 2011 Published: May 20, 2011 9156
dx.doi.org/10.1021/ie102321b | Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
(3) I component tanks, their initial inventories, limits on their holdups, flow profiles of feeds into the tanks, and limits on the flows out of the tanks; (4) P products and specification limits on their property indices; (5) B blenders, the products that each blender can process, lower limits on the blend times of these products, and limits on their blending rates; (6) J Product tanks, the products that each tank can store, limits on their holdups, the products and holdups at time zero, and delivery (lifting) rates for various products; (7) O orders, their constituent products, amounts, and delivery time windows; and (8) Revenues from product sales, component costs, inventory costs (for components and products), and demurrage costs for orders. Determine: (1) The blenders that each component tank should feed over time, and their feed rates; (2) The products that each blender should produce over time, and their production rates; (3) The products that each product tank should receive over time, from which blender, and at what flow rates; (4) The orders that each product tank should deliver over time and their amounts; and (5) The inventory profiles of various tanks (component and product). Assuming: (1) Flow rate profile of each component from the upstream process is piecewise constant; (2) Component quality profile is also piecewise constant; (3) Mixing in each blender is perfect; (4) Changeover times between products are negligible for product tanks; (5) Changeover times between product runs on blenders are product- and sequence-independent; (6) Each order involves only one product; as discussed earlier, each multiproduct order can be decomposed into several single-product orders; (7) Each order is completed during the scheduling horizon; and (8) Each order is single-product, i.e., multiproduct orders have been broken up into individual one-product orders. Subject to the operating rules: (1) A blender can process, at most, one product at any time. Once it begins processing a product, it must operate for some minimum time before it can switch to another product. (2) A blender can feed, at most, one product tank at any time. In addition to being the industry practice, this helps to decrease the number of tanks in use and increase their utilization. Allowing: (1) A component tank may receive and feed components at the same time; (2) A component tank may feed some or all blenders simultaneously; (3) Multiple component tanks may feed a blender at the same time; (4) A product tank may receive and feed products at the same time;
charging, at most, one product tank at a time, minimum run lengths, etc. They ensured the blending rate to be constant during a run, and they developed a schedule adjustment procedure to avoid solving nonconvex MINLPs. However, the use of process slots makes their model slow on larger examples, because of the low RMIPs and slow improvement of RMIPs. In addition, they assumed unlimited component inventories over the scheduling horizon to avoid more slots. This may not always hold in practice. In this paper, we develop a multigrid continuous-time formulation to improve the efficiency of treating recipe, blending, storage, and order scheduling simultaneously. We achieve this by employing unit slots11 in place of process slots. In addition to including all the operational features of the model of Li et al.,8 we allow limited component inventories, simultaneous receipt/ delivery by product tanks, and blender setup times. We also revise the schedule adjustment procedure of Li et al.8 for unit slots. We use 14 examples from Li et al.8 to evaluate our model and algorithm.
2. PROBLEM STATEMENT Consider Figure 1, which shows a gasoline blending unit (GBU) in a typical refinery. It involves I blend components (i = 1, 2, ..., I) in exactly I dedicated component tanks (i = 1, 2, ..., I), B blenders (b = 1, 2, ..., B), J multipurpose product tanks (j = 1, 2, ..., J), and some lifting ports. It makes P grades of gasoline, which we call products (p = 1, 2, ..., P). At time zero, GBDU has O orders (o = 1, 2, ..., O) with time windows [DDLo , DDU o , where o = 1, 2, ..., O] to fulfill during the scheduling horizon [0, H]. An order may involve multiple gasoline products, but we can break each multiproduct order into several single-product orders. Thus, each order involves a single product. A delivery of order o after DDU o incurs a demurrage (DMo, given in dollars per unit time of delay). The blend components are either gasoline fractions from various processes in a refinery (such as atmospheric distillation, fluid catalytic cracking unit (FCCU), and isoforming unit (IFU)) or additives (such as methyl tert-butyl ether (MTBE) and butane), which enhance the octane rating or act as corrosion inhibitors or lubricators. Each component i has a prefixed, distinct, and known quality or specification, The quality of a blend component is specified in terms of various property indices such as octane number (ON), Reid vapor pressure index (RVPI), sulfur index (SI), viscosity index (VI), etc. The flow profiles, over time, of various blend components into respective component tanks are known a priori. The components feed blend components into the blenders in various proportions to make various products over time. The blenders blend products one at a time in a semicontinuous mode. When a blender switches from one product to another, a minimum setup time (τb) is required. A product tank is assigned to receive product from each blend run. A product tank may hold different products over time. Subsequently, the product tanks load the product orders into vehicles or ships at appropriate times for delivery. We refer to this part as order scheduling. The GBDU operation involves decisions such as recipe determination, allocation of component tanks to blenders, assignment of product tanks to products, and scheduling of blending, transfer, and order loading. Thus, the gasoline blending and scheduling problem addressed in this paper can be stated as follows: Given: (1) A scheduling horizon [0, H]; (2) I components and profiles of their property indices over time; 9157
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Figure 1. Schematic of gasoline blending and distribution.
(5) A product tank may deliver multiple orders at the same time; and (6) Multiple product tanks may deliver an order at the same time. Minimizing: The total operating cost includes material (component), demurrage, transition, and back-order costs. The above problem description differs in three aspects from the one addressed by Li et al.8 First, we allow a product tank to receive and deliver products at the same time; Li et al.8 assumed that a tank cannot feed and deliver at the same time. Second, we do not assume unlimited supply of components from the component tanks over the scheduling horizon. In practice, the scheduling of a particular blend may be delayed due to the unavailability of one or more necessary components or an inability to match the quality from the existing supply of components. We now develop a continuous-time, multigrid11 MILP formulation for the above problem. In addition to the distinguishing problem features mentioned above, the main difference between our new formulation and that of Li et al.8 is the use of unit slots instead of process slots. The advantages of using unit slots are well-established11 in the literature. Our contribution is to extend the multigrid methodology of Susarla et al.,11 which was for multipurpose batch plants, to semicontinuous plants. Susarla et al.11 have detailed the differences and similarities between slotbased and event-based formulations. For the ease of exposition, we make two simplifying assumptions, which we relax later. First, we consider the single-period (SP) scenario, in which the flows of all components into the respective component tanks are constant over the entire scheduling horizon. Second, a product tank cannot receive and deliver products simultaneously. Throughout this paper, each variable is defined with specific ranges of its indices, and each constraint, unless otherwise stated, is written for all valid values of the indices of its constituent variables.
3. MILP FORMULATION We define component tanks (i), blenders (b), and product tanks (j) as units (n). For each unit n, we divide the scheduling horizon [0, H] into K (k = 1, 2, ..., K) ordered contiguous unit slots (see Figure 2). We define the period before time zero as slot zero (k = 0). Let Tnk (k = 0, 1, ..., K; Tn0 g 0, TnK e H) denote the
end time of slot k on unit n, where n becomes i for a component tank, b for a blender, and j for a product tank. Slot k on unit n starts at Tn(k1) and ends at Tnk. In contrast to Li et al.,8 the slots are asynchronous, so Tik, Tbk, and Tjk may vary with the unit and are partially independent.11 We define the following sets in this paper: BP = {(b, p) | blender b can process product p} BJ = {(b, j) | blender b can feed product tank j} JP = {(j, p) | product tank j can hold product p} JO = {(j, o) | product tank j may deliver order o} PI = {(p, i) | product p may use component i} OP = {(o, p) | order o is for product p} Since the slots on each unit n are chronologically ordered, we have Tnk g Tnðk1Þ
ðfor 1 e k e KÞ
ð1Þ
We assume that blending begins always from the start of a slot, but may end at any time within the slot. In other words, the idle time, if any, is always toward the end of the slot. Unless otherwise indicated, an index takes all its legitimate values in all the expressions or constraints in our formulation. 3.1. Blending and Storage. At any time, a blender must be either running or idle. However, we assume that it is always connected to a product tank. To facilitate this, we define a dummy product tank (j = 0). Thus, we have J real product tanks (j = 1, 2, ..., J) and one dummy (j = 0) product tank. If a blender is not running, then we connect it to the dummy tank (j = 0), otherwise it must be connected to a real product tank. Furthermore, we assume that a product tank always has a product. If it is empty, then we assume that it is holding the last product that it had. Chart 1 defines the binary variables used to model the blending and storage operations in the GBDU. Li et al.8 also used the aforementioned variables. They treated vbjk (j = 1, 2, ..., J) as binary and vb0k as 01 continuous, and they proved that xbpk would be binary automatically and, hence, can be treated as 01 continuous. Note that we must assign proper values for xbp0, based on the product that blender b was processing before time zero, and ujp0, based on what was inside tank j before time zero. For a formulation based on unit slots, Susarla et al.11 argued that, if a processing unit and a storage unit are exchanging material, then the unit slots in which this exchange occurs must have the same index on both the unit 9158
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
and the storage. On that basis, we use the following constraints from Li et al.8 in our formulation: vb0k þ
∑
∑ vbjk ¼ 1 j:ðb, jÞ ∈ BJ
ujpk ¼ 1
ðfor 1 e k e KÞ
ðfor 1 e j e J, 1 e k e KÞ
xebk þ xbpk þ xbpðkþ1Þ e 2
ðfor ðb, pÞ ∈ BP, 0 e k < KÞ ð7Þ
ð2Þ
uejk g ujpk ujpðkþ1Þ ðfor ðb, pÞ ∈ JP, 0 e k < K, 1 e j e JÞ
ð3Þ
ð8aÞ
p:ðj, pÞ ∈ JP
vb0k þ
∑ xbpk ¼ 1 p:ðb, pÞ ∈ BP
xbpk g ujpk þ vbjk 1
ðfor 1 e k e KÞ
ðfor ðb, pÞ ∈ JP, 0 e k < K, 1 e j e JÞ
ðfor ðb, pÞ ∈ BP, ðb, jÞ ∈ BJ,
ðj, pÞ ∈ JP, 1 e j e J, 1 e k e KÞ ujpk g xbpk þ vbjk 1
uejk g ujpðkþ1Þ ujpk
ð4Þ
B K1
ð5aÞ
ðfor ðb, pÞ ∈ BP, ðb, jÞ ∈ BJ,
ðj, pÞ ∈ JP, 1 e j e J, 1 e k e KÞ
∑ ∑ xebk g NP B b¼1 k¼1
ð5bÞ
xebk g xbpk xbpðkþ1Þ
ðfor ðb, pÞ ∈ BP, 0 e k < KÞ ð6aÞ
xebk g xbpðkþ1Þ xbpk
ðfor ðb, pÞ ∈ BP, 0 e k < KÞ ð6bÞ
ð8bÞ ð9Þ
where NP is the number of distinct products that must be processed during the scheduling horizon. Since we will impose a penalty for product transitions in a product tank, we do not force uejk = 0. If a blender b is idle during slot k, then it cannot receive from any component tank during k. Otherwise, it must receive at least one component. Therefore, we write vb0k þ yibk e 1 vb0k þ
I
∑ yibk g 1
i¼1
ðfor 1 e k e KÞ ðfor 1 e k e KÞ
ð10aÞ ð10bÞ
3.2. Blend Scheduling and Product Quality. Let BLbk be the time for which blender b processes a product p (p > 0) during its slot k. Li et al.8 did not need this variable. Clearly, BLbk must equal slot length, unless the current blend run is ending during slot k, or the blender is idle. In other words,
BLbk e Tbk Tbðk1Þ
ðfor 1 e k e KÞ
BLbk þ Hðvb0k þ uebk Þ g Tbk Tbðk1Þ
ð11aÞ
ðfor 1 e k e KÞ ð11bÞ
If the blender is idle during slot k, then BLbk must be zero. BLbk e Hð1 vb0k Þ
Figure 2. Schematic of unit-slot design.
ðfor 1 e k e KÞ
ð11cÞ
Chart 1. Binary Variables Used to Model the Blending and Storage Operations in the GBDU
9159
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Next, following Li et al.,8 we define RLbk (0 e k e K) as the length of the current run on blender b at the end of slot k, if the run does not end during slot k; otherwise, RLbk = 0. 8 > > < current run length ðif the current run of blender b does not end during its slot kÞ RLbk ¼ > > :0 ðotherwiseÞ Thus, RLb0 = 0 if a run has ended at time zero; otherwise, it is the current run length at time zero. To compute RLbk, we write RLbk e RLbðk1Þ þ BLbk
ðfor 1 e k e KÞ
RLbk e Hð1 xebk Þ
ðfor 0 e k e KÞ
ð12aÞ
Lastly, Qbk must respect the limits (FLb and FU b ) on the processing rates of blender b. FbL 3 BLbk e Qbk e FbU BLbk
∑
ðfor ðb, pÞ ∈ BP, 1 e k e KÞ I
∑ RLLbpxbpk p:ðb, pÞ ∈ BP
RLbðk1Þ þ BLbk þ RLLb ð1 xebk Þ g
ðfor 1 e k < KÞ RLbðK1Þ þ BLbK g
P
∑
p:ðb, pÞ ∈ BP
where RLLb ¼
max
p:ðb, pÞ ∈ BP
RLLbp xbpK
ð13bÞ
I
J
∑
i¼1
qibk Fi θis g
I
∑
i¼1
I
ðfor 1 e k e KÞ
ð14aÞ
ðfor 1 e k e KÞ
ðfor 1 e k < KÞ
ð14cÞ
Blending requires components and must ensure on-spec products. Let qibk be the volume that blender b receives from component tank i during its slot k; let Gbjk be the volume that it feeds into product tank j; and let Qbk be the volume that it processes. Then, a mass balance around blender b gives us Gbjk ∑ qibk ¼ j:ðb, ∑ i¼1 jÞ ∈ BJ
I
ð19bÞ
!
qibk Fi θLps Mb θLps minðθLps Þ Fmax ð1 xbpk Þ p
ð20aÞ
ðfor ðb, pÞ ∈ BP, 1 e k e KÞ
ð14bÞ
Qbk ¼
p
ðfor ðb, pÞ ∈ BP, 1 e k e KÞ
RLLbp
STbk e STbðk1Þ þ Tbk Tbðk1Þ BLbk
STbk g τb ½xebk vb0ðkþ1Þ
∑
i¼1
ð19aÞ
qibk θis e Qbk θUps þ Mb max ðθUps Þ θUps ð1 xbpk Þ
ð13aÞ
When a blender switches from one product to another, it requires some minimum setup time (τb). Since we allow a blender to be idle across several consecutive slots after it ends a blend run, we must compute the cumulative idle time between successive blend runs. To this end, we define STbk as the accumulated idle time of a blender b by the end of slot k. We compute this using the following constraints: STbk e Hðxebk þ vb0k Þ
p
i¼1
ð12bÞ
P
ð18a; bÞ
The above allows the blending rate to vary from slot to slot during a blend run. This is not realistic, and we follow the procedure of Li et al.8 to address it later. For gasoline quality, Li et al.8 identified nine most commonly used indices and their additive bases. They proposed the following constraints to ensure on-spec products: I L L L qibk θis g Qbk θps Mb θps minðθps Þ ð1 xbpk Þ
Then, to ensure a minimum length (RLLbp) for each blend run, we demand
ðfor 1 e k e KÞ
ðfor 1 e k e KÞ ð15a; bÞ
If blender b receives nothing from component tank i or is idle, or it does not feed product tank j in slot k, then qibk and Gbjk must be zero. qibk e Mb yibk
ðfor 1 e k e KÞ
ð16Þ
Gbjk e Mb vbjk
ðfor 1 e k e KÞ
ð17Þ
I
∑ qibk Fi θis e i∑¼ 1 qibk Fi i¼1
!
θUps þ Mb max ðθUps Þ θUps Fmax ð1 xbpk Þ p
ð20bÞ
ðfor ðb, pÞ ∈ BP, 1 e k e KÞ
Qbk rpiL Mb rpiL minðrpiL Þ ð1 xbpk Þ e qibk e Qbk rpiU p þ Mb max ðrpiU Þ rpiU ð1 xbpk Þ p
ðfor ðb, pÞ ∈ BP, 1 e k e KÞ
ð21a,bÞ
where θis is the known blending index for a property s of component i, Fi is the density of component i, [θLps,θU ps] are the desired limits on property s of product p, and Fmax is the maximum possible density among all products. Equations 19 are for volume-based indices, eqs 20 are for weight-based indices, and the volume fraction of component i in product p satisfies limits [rLpi,rU pi]. 3.3. Order Scheduling. Li et al.8 defined one binary variable (zjok) for order scheduling operations: ( 1 if product tank j delivers order o during its slot k zjok ¼ 0 ðotherwiseÞ ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
They presented the following constraints that are applicable in our new model as well: J
K
∑ ∑ zjok g 1 j:ðj, oÞ ∈ JO k ¼ 1
U U Here, Mb = min[HFU b , max j VPj ] and VPj is the capacity of tank j.
9160
ð22Þ
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Figure 3. An example schedule to illustrate the delivery starts of O1, O2, and O3 by PT-101 at different times in the same slot.
vbjk þ zjok e 1 ðfor ðb, jÞ ∈ BJ, ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
slot k on tank j. ð23Þ
tsjok g Tjðk1Þ ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
zjok e ujpðk1Þ ðfor ðp, jÞ ∈ PJ, ðj, oÞ ∈ JO, ðo, pÞ ∈ OP, 1 e j e J, 1 e k e KÞ
tsjok þ
ð24Þ
ðfor ðp, jÞ ∈ PJ, ðj, oÞ ∈ JO, ðo, pÞ ∈ OP, 1 e j e J, 1 e k e KÞ
ð25Þ zjoðkþ1Þ þ uejk e 1 ðfor ðj, oÞ ∈ JO, 1 e j e J, 0 e k < KÞ
tsjok g DDLo zjok
ð26Þ
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ ð27Þ
J
∑
K
∑ DQjok ¼ TQo
do g tsjok þ
ð29Þ
where DQjok denotes the volume of order o delivered by tank j during slot k, TQo is the required amount for order o, and DRU j is the maximum cumulative rate for all orders. While a product tank may deliver multiple orders during a slot, Li et al.8 required them to begin simultaneously. We allow a delivery to start and end at any time during the slot. For instance, consider Figure 3. PT-101 begins delivering O1 at t = 2 and O2 at t = 4 in the same slot that starts at t = 1 and ends at t = 6. To know these times precisely, we define tsjok as the time at which product tank j begins delivering order o in slot k. Since the delivery rate (DRjo) is constant, the end time for this delivery is tsjok þ DQjok/DRjo. We force this delivery to be entirely within
DQjok DDUo Hð1 zjok Þ DRjo
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
ð28Þ
j:ðj, oÞ ∈ JO k ¼ 1
ð31Þ
If it ends after DDU o , then we have a delivery delay:
DQjok e DRjU ½Tjk Tjðk1Þ ∑ o:ðj, oÞ ∈ JO ðfor 1 e j e J, 1 e k e KÞ
ð30bÞ
Equations 30 allow order delivery to be intermittent from a tank, which is not realistic. We discuss this issue later. In addition to eq 30, we must ensure that a delivery does not begin before its own time window [DDLo , DDU o ].
DQjok e TQo zjok ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
DQjok e Tjk DRjo
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
zjok e ujpk
ð30aÞ
ð32Þ
Finally, the delivery delay (do) must be bounded by (H DDU o ). 3.4. Slot Timings on Units. When using unit slots in dealing with shared resources such as inventories, the main challenge is to monitor the resource usage/production by the various sharing units accurately. The time instances corresponding to the flows in/out of a resource must be ordered chronologically, so that we can get a correct resource profile. To this end, we assume that a blending activity is always left justified (in time) within a slot. In other words, a blend activity must begin at the slot start, but may end before the slot end with an idle time. We also recall that each component tank has a constant inflow of material from upstream units in a singleperiod scenario that we are considering for now. Now, as explained in the Appendix, whenever a blender begins a new blend run in a slot, and a component tank is feeding that blender, then we require that the tank have a matching slot with 9161
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Tbðk1Þ þ BLbk e tik þ H½1 yibk þ yibðkþ1Þ
ðfor 1 e k < KÞ
ð34cÞ Tbðk1Þ þ BLbk g tik H½1 yibk þ yibðkþ1Þ
ðfor 1 e k < KÞ
ð34dÞ While eq 33 applies only when a component tank begins to feed a blender during a blend run or a new run starts, eq 34 is effective only when a component tank ends its feed during a blend run or a blend run ends. If a blend run continues across successive slots with the same component tanks, then no new slots begin. Lastly, the synchronizations in eqs 33 and 34 are necessary, because of the continuous nature of blending operations. This is in contrast to the instantaneous flows that could be reasonably assumed in a batch process.11 For product tanks, the situation is simpler, because they cannot receive and deliver at the same time. So flow into and flow out from a product tank can never occur in the same slot. Therefore, if a blender feeds tank j during slot k, then the start (end) of a slot k on a product tank j must precede (succeed) the start (end) of slot k on the blender. Tjk e Tbk þ H½1 vbjðkþ1Þ ðfor ðb, jÞ ∈ BJ, 1 e j e J, 0 e k < KÞ
ð35aÞ
Tjðkþ1Þ g Tbk þ BLbðkþ1Þ H½1 vbjðkþ1Þ ðfor ðb, jÞ ∈ BJ, 1 e j e J, 0 e k < KÞ
3.5. Inventory Balance. The timing constraints enable us to write the following inventory balances for component and product tanks:
Figure 4. Flowchart for the schedule adjustment procedure.
VCik ¼ VCiðk1Þ þ Fi ½Tik Tiðk1Þ
the same index and start time. Tbk g Tik H½2 xebk yibðkþ1Þ
ðfor 0 e k < KÞ ð33aÞ
Tbk e Tik þ H½2 xebk yibðkþ1Þ
ðfor 0 e k < KÞ ð33bÞ
Similarly, we demand that, if a tank begins feeding a blender during a blend run, then also it should have a matching slot with the same index and start time. In other words, Tbk g Tik H½1 yibðkþ1Þ þ yibk
ðfor 0 e k < KÞ ð33cÞ
Tbk e Tik þ H½1 yibðkþ1Þ þ yibk
ðfor 0 e k < KÞ ð33dÞ
To check the inventory in a tank during or at the end of a blend run, we define an intermediate point tik in each slot, i.e., Ti(k1) e tik e Tik. We then demand that the ends of each blend run and tank feed during a blend run should match the corresponding points on the tank. Tbðk1Þ þ BLbk e tik þ Hð2 xebk yibk Þ
ðfor 1 e k e KÞ
ð34aÞ Tbðk1Þ þ BLbk g tik Hð2 xebk yibk Þ
ð35bÞ
ðfor 1 e k e KÞ
B
∑ qibk b¼1
ðfor 1 e k e KÞ VPjk ¼ VPjðk1Þ þ
B
ð36Þ
O
∑ Gbjk o∑¼ 1 DQjok b¼1
ðfor ðj, oÞ ∈ JO, 1 e k e KÞ
ð37Þ
where Fi is the constant feed rate to tank i from upstream units, VCik (VCLi e VCik e VCU i ) is the inventory of component i at Tik, ) is the holdup in product tank j at Tjk. VPjk (0 e VPjk e VPU j Equation 37 computes inventory levels only at slot points. Since a blend run may end before the slot, we must also compute the inventory level at tik. Therefore, we use VCLi e VCiðk1Þ þ Fi ½tik Tiðk1Þ
B
∑ qibk e VCUi b¼1
ðfor 1 e k e KÞ
ð38a; bÞ
Note that eq 36 neglects the inflows into the component tanks after the end of schedule activities. In other words, eq 36 does not compute the inventory at the end of the horizon if the schedule ends before H. Lastly, a product transition cannot occur on a product tank j, unless the tank holdup at Tjk is zero. Thus, VPjk e VPUj ð1 uejk Þ
ð34bÞ 9162
ðfor 1 e k < KÞ
ð39Þ
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 1. Order Data for (a) Examples 19 and (b) Examples 1014 (a) Order Data for Examples 19
o
p
1
2
Amount
Delivery Rate
( 103 bbl)
( 103 bbl/h)
Delivery Window
3 4, 5 6 7, 8 9 1 2 3 4, 5 6 7, 8 9
1
2
3
4, 5
6
7, 8
9
O1 P1 11 10 10 10 10 10 10 3 5 5
5
5
5
5 [0,12]
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O2 P2 3
3
3
3
3
3
3 3 3 3
3
3
3
3 [0,12]
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O3 P2 60
3
3
3
3
3
3 5 3 3
3
3
3
3 [24,48]
[24,50]
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O4 P1 125 10 10 10 10 10 10 5 5 5
5
5
5
5 [24,72]
[24,50]
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O5 P2 3 O6 P1 -
3 3 3 3 3 3 3 3 3 10 10 10 10 10 10 - 5 5
3 5
3 5
3 5
3 [24,48] 5 -
[48,72] [48,72]
[24,48] [24,48]
[24,48] [24,48]
[24,48] [24,48]
[24,48] [24,48]
[24,48] [24,48]
O7 P2 -
3
- 3 3
3
3
3
3
-
[48,120]
[24,48]
[24,48]
[24,48]
[24,48]
[24,48]
O8 P1 - 120 100 100 90 100 100 - 5 5
5
5
5
5
-
[118,190]
[118,190]
[118,190]
[118,190]
[118,190]
[118,190]
O9 P2 -
- 3 3
3
3
3
3
-
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
O10 P4 - 150 150 150 150 150 150 - 5 5
5
5
5
5
-
O11 P3 -
-
20 20 45 60 60 - - 5
5
5
5
5
-
-
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
O12 P2 -
-
30 30 30 20 20 - - 5
5
5
5
5
-
-
[24,48]
[24,48]
[24,48]
[24,48]
[24,48]
O13 P4 O14 P3 -
-
-
60 45 60 60 - - 10 15 15 20 - - -
5 5
5 5
5 5
5 5
-
-
-
[0,56] [48,72]
[0,56] [48,72]
[0,56] [48,72]
[0,56] [48,72]
O15 P2 -
-
-
20 15 20 20 - - -
4
4
4
4
-
-
-
[0,72]
[0,75]
[0,72]
[0,72]
O16 P2 -
-
-
-
10 20 20 - - -
-
4
5
5
-
-
-
-
[48,72]
[48,72]
[48,72]
O17 P1 -
-
-
-
10 10 10 - - -
-
5
5
5
-
-
-
-
[48,96]
[48,72]
[48,72]
O18 P1 -
-
-
-
10 10 10 - - -
-
5
5
5
-
-
-
-
[48,96]
[48,72]
[48,72]
O19 P2 -
-
-
-
-
60 60 - - -
-
-
5
5
-
-
-
-
-
[0,50]
[0,50]
O20 P2 -
-
-
-
-
40 40 - - -
-
-
5
5
-
-
-
-
-
[144, 168]
[144,168]
O21 P1 O22 P5 -
-
-
-
-
-
30 - - 40 - - -
-
-
-
5 5
-
-
-
-
-
-
[96,120] [144,168]
O23 P3 -
-
-
-
-
-
20 - - -
-
-
-
5
-
-
-
-
-
-
[144,168]
3
3 3
3 3
3 3
3
3
3
3
[150.5,185.5] [150.5,185.5] [150.5,185.5] [150.5,185.5] [150.5,185.5] [150.5,185.5]
(b) Order Data for Examples 1014 Amount ( 103 bbl)
Product o
12 the rest 10
Delivery Rate ( 103 bbl/h)
Delivery Window
11
12
13
14
10
11
12
13
14
10
11
12
13
14
O1 P1
P1
10
10
10
10
10
5
5
5
5
5
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O1 P1 O2 P2
P1 P2
10 3
10 3
10 3
10 3
10 3
5 3
5 3
5 3
5 3
5 3
[0,24] [0,24]
[0,24] [0,24]
[0,24] [0,24]
[0,24] [0,24]
[0,24] [0,24]
O3 P2
P2
3
3
3
3
3
3
3
3
3
3
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O4 P1
P1
10
10
10
10
10
5
5
5
5
5
[0,24]
[0,24]
[0,24]
[0,24]
[0,24]
O5 P2
P2
3
3
3
3
3
3
3
3
3
3
[24,48]
[24,48]
[24,48]
[24,48]
[24,48]
O6 P1
P1
10
10
10
10
10
5
5
5
5
5
[24,48]
[24,48]
[24,48]
[24,48]
[24,48]
O7 P2
P2
3
3
3
3
3
3
3
3
3
3
[24,48]
[24,48]
[24,48]
[24,48]
[24,48]
O8 P1
P1
100 100 100 100 100
5
5
5
5
5
[118,190]
[118,190]
[118,190]
[118,190]
[118,190]
O9 P2 O10 P4
P2 P4
3 3 3 3 3 150 150 150 150 150
3 5
3 5
3 5
3 5
3 5
O11 P3
P3
60
60
60
60
60
5
5
5
5
5
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
O12 P2
P2
20
20
20
20
20
5
5
5
5
5
[24,48]
[24,48]
[24,48]
[24,48]
[24,48]
O13 P4
P4
60
60
60
60
60
5
5
5
5
5
[0,56]
[0,56]
[0,56]
[0,56]
[0,56]
O14 P3
P3
20
20
15
20
20
5
5
5
5
5
[48,72]
[48,72]
[48,72]
[48,72]
[48,72]
O15 P2
P2
20
20
20
20
20
4
4
4
4
4
[0,72]
[0,72]
[0,72]
[0,72]
[0,72]
O16 P2
P2
20
20
20
20
20
5
5
5
5
5
[48,72]
[48,72]
[48,72]
[48,72]
[48,72]
O17 P1 O18 P1
P1 P1
10 10
10 10
10 10
10 10
10 10
5 5
5 5
5 5
5 5
5 5
[48,72] [48,72]
[48,72] [48,72]
[48,72] [48,72]
[48,72] [48,72]
[48,72] [48,72]
O19 P2
P2
60
60
60
60
60
5
5
5
5
5
[0,50]
[0,50]
[0,50]
[0,50]
[0,50]
O20 P2
P2
40
40
40
40
40
5
5
5
5
5
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
[144,168] [144,168] [144,168] [144,168] [144,168] [150.5,185.5] [150.5,185. 5] [150.5,185.5] [150.5,185.5] [150.5,185.5]
9163
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 1. Continued (b) Order Data for Examples 1014 3
o
Delivery Rate ( 103 bbl/h)
Amount ( 10 bbl)
Product 12 the rest 10
11
12
13
14
10
11
12
13
14
Delivery Window 10
11
12
13
14
O21 P5
P1
30
30
30
30
30
5
5
5
5
5
[96,120]
[96,120]
[96,120]
[96,120]
[96,120]
O22 P5
P5
40
40
40
40
40
5
5
5
5
5
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
O23 P3
P3
20
20
20
20
20
5
5
5
5
5
[144,168]
[144,168]
[144,168]
[144,168]
[144,168]
O24 P5 O25 P5
P5 P5
6 20
6 20
6 20
6 20
6 20
3 5
3 5
3 5
3 5
3 5
[96,120] [144,168]
[96,120] [144,168]
[96,120] [144,168]
[96,120] [144,168]
[96,120] [144,168]
O26 P3
P1
-
10
30
10
10
-
4
4
4
4
-
[0,76]
[144,168]
[0,76]
[0,76]
O27 P3
P4
-
20
20
20
20
-
5
4
5
5
-
[120,144]
[72,96]
[120,144]
[120,144]
O28 P4
P1
-
25
3
25
25
-
5
3
5
5
-
[120,144]
[72,96]
[120,144]
[120,144]
O29 P4
P5
-
10
15
10
10
-
5
3
5
5
-
[120,144]
[96,120]
[120,144]
[120,144]
O30 P1
P4
-
15
15
15
15
-
5
3
5
5
-
[120,144]
[96,120]
[120,144]
[120,144]
O31 P2
P1
-
-
15
15
15
-
-
5
5
5
-
-
[96,120]
[120,144]
[120,144]
O32 P5 O33 P1
P1 P4
-
-
20 20
20 20
20 20
-
-
2 5
5 5
5 5
-
-
[96,120] [0,76]
[144,168] [144,168]
[144,168] [144,168]
O34 P3
P4
-
-
20
20
20
-
-
5
5
5
-
-
[120,144]
[168,192]
[168,192]
O35 P3
P5
-
-
30
30
30
-
-
5
5
5
-
-
[120,144]
[168,192]
[168,192]
O36 -
P2
-
-
-
3
3
-
-
-
3
3
-
-
-
[168,192]
[168,192]
O37 -
P1
-
-
-
10
10
-
-
-
5
5
-
-
-
[168,192]
[168,192]
O38 -
P1
-
-
-
40
40
-
-
-
5
5
-
-
-
[168,192]
[168,192]
O39 -
P4
-
-
-
10
10
-
-
-
5
5
-
-
-
[168,192]
[168,192]
O40 O41 -
P5 P1
-
-
-
10 -
10 15
-
-
-
5 -
5 5
-
-
-
[168,192] -
[168,192] [168,192]
O42 -
P2
-
-
-
-
20
-
-
-
-
3
-
-
-
-
[168,192]
O43 -
P3
-
-
-
-
15
-
-
-
-
5
-
-
-
-
[144,168]
O44 -
P5
-
-
-
-
20
-
-
-
-
4
-
-
-
-
[168,192]
O45 -
P4
-
-
-
-
10
-
-
-
-
5
-
-
-
-
[96,120]
3.6. Scheduling Objective. The total operating cost of GBDU from Li et al.8 is given by
TC ¼
I
B
K
B
K1
∑ ∑ ∑ ci qibk þ b∑¼ 1 k∑¼ 1 CBb xebk i¼1 b¼1 k¼1
þ
J
K1
O
∑ ∑ CTj uejk þ o∑¼ 1 DMo do j¼1 k¼1
ð40Þ
where ci is the price of component i ($ per unit volume), CBb is the cost of transition on blender b ($ per occurrence), CTj is the cost of transition in product tank j ($ per occurrence), and DMo is the demurrage cost of order o ($ per unit time). This completes our multigrid formulation for GBDU operations. It comprises eqs 140. It has three main limitations. First, it allows the blending rate to vary from slot to slot and order delivery to be discontinuous; these are undesirable in practice. Second, it assumes constant feed rates for components, or it is a single-period formulation. Third, it does not allow a product tank to receive and deliver product at the same time. For the first limitation, we revise the adjustment procedure of Li et al.8 slightly to accommodate unit slots. The second limitation or the extension to multiple periods is readily addressed as done by Li et al.8 for their model, so we simply state the modified equations without elaboration. However, the third limitation requires a more-detailed consideration. We begin with the first limitation.
4. SCHEDULE ADJUSTMENT A schedule adjustment procedure is necessary to eliminate the slot-to-slot variability of blend rates. The main modification required in the procedure of Li et al.8 concerns the blend length variable (BLbk), which their formulation did not require. Given a solution from our basic model, we compute the correct values for run lengths (CRLbk), volumes (CCQbk), and total volume (TCQbk) processed by a blender in a run as follows: CRLbk ¼ CCQbk ¼ 0
ðif xebk ¼ 1Þ
ð41a,bÞ
CRLbk ¼ CRLbðk1Þ þ ½BLbk
ðif xebk ¼ 0Þ
ð42aÞ
CCQbk ¼ CCQbðk1Þ þ ½Qbk
ðif xebk ¼ 0Þ
ð42bÞ
TCQbk ¼ 0
ðif xebk ¼ 0Þ
TCQbk ¼ CCQbðk1Þ þ ½Qbk
ðif xebk ¼ 1Þ
ð43aÞ ð43bÞ
where [BLbk] and [Q bk] are the values of BL bk and Q bk , repectively. Using the above parameters, the blending rate (R bk) for each blending run at the slot where it ends is given by Rbk ¼ 9164
CCQbðk1Þ þ Qbk CRLbðk1Þ þ BLbk
for k with xebk ¼ 1 and vb0k ¼ 0
ð44Þ
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
1
P3
P3
P2
P2
P1
-
-
-
-
-
C1
C2
C3
C4 C5
C6
C7
C8
C9
tank
PT-101
PT-102
PT-103
PT-104
PT-105
PT-106
PT-107 PT-108
PT-109
PT-110
9165
PT-111
CT-101
CT-102
CT-103
CT-104 CT-105
CT-106
CT-107
CT-108
CT-109
C9
C8
C7
C6
C4 C5
C3
C2
C1
P1
P1
P4
P1 P1
P2
P2
P4
P2
P3
P3
214
Initial Product
13.84
49.34
27.30
19.53
28.29 10.59
59.44
34.58
27.38
-
-
-
-
-
20.20
28.49
14.08
0.00
90.20
1
44.58
49.47
46.91
19.53
44.44 10.59
59.44
67.90
26.46
12.36
60.00
23.96
13.79 12.36
57.59
28.49
25.00
14.08
0.00
30.00
214
Initial Stock ( 103 bbl)
200
100
100
200
250 200
250
200
200
-
-
-
-
-
100
100
150
100
100
1
250
250
250
250
300 200
300
300
250
150
150
200
200 150
150
200
200
150
150
150
214
Capacity ( 103 bbl)
Table 2. Product and Component Tank Data for Examples 114
P2, P3
P2, P3
P2, P3
28
-
-
-
-
-
-
-
-
-
-
-
-
-
P1,P2
-
-
-
-
-
-
-
-
P1,P4
P1,P4
P1,P4
P1,P4 P1,P4
P2,P3
P2, P3
P1,P2 P2P4
P2,P3
P2,P3
P2,P3
1
13, 14
12
-
-
-
-
-
-
-
-
P1,P4
P1,P4
P1,P4
P1,P4 P1,P4
P2,P5
P2,P5
P2P4
-
-
-
-
-
-
-
-
P1,P4
P1,P4
P1,P4
P1,P4 P1,P4
P2,P5
P2,P5
P2P4
-
-
-
-
-
-
-
-
P1,P4
P1,P4
P1,P4
P1,P4 P1,P4
P2,P3,P5
P2,P3,P5
P2P5
P2, P3, P5 P2, P3, P5 P2, P3, P5
P2, P3, P5 P2, P3, P5 P2, P3, P5
P2, P3, P5 P2, P3, P5 P2, P3, P5
911
Storable Products
1.00
0.00
0.00
1.00
1.20 1.00
1.20
1.00
1.00
-
-
-
-
-
15
15
15
15
15
1
1.00
0.00
0.00
0.50
1.00 0.80
1.00
0.50
1.00
20
20
20
20 20
20
20
20
20
20
20
25
1.00
0.00
0.00
0.50
1.00 0.80
1.00
0.50
1.00
25
25
25
25 25
25
25
25
25
25
25
6
1.00
0.00
0.00
0.50
1.00 0.80
1.00
0.50
1.00
20
20
20
20 20
20
20
20
20
20
20
7, 8
1.00
0.00
0.00
0.50
1.00 0.80
1.00
0.50
1.00
30
30
30
30 30
30
30
30
30
30
30
9, 10
1.00
0.00
0.00
0.50
1.00 0.70
1.00
0.50
1.00
30
30
30
30 30
30
30
30
30
30
30
11
1.00
0.00
0.00
0.50
1.00 0.80
1.00
0.50
1.00
30
30
30
30 30
30
30
30
30
30
30
12
Maximum Delivery Rate/Feed Rate ( 103 bbl/h)
1.00
0.50
0.50
0.50
1.00 0.70
1.00
0.50
1.00
30
30
30
30 30
30
30
30
30
30
30
13, 14
Industrial & Engineering Chemistry Research ARTICLE
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 3. Component and Product Property Indices for Examples 114 RBN comp./ product
1
28
RVI
SI
BI
AROI
914
28
914
28
914
28
914
28
914
C1
86.50
86.50
86.50
140.47
140.47
80.00
80.00
0.78
0.78
25.00
25.00
C2
103.66
103.66
103.66
68.92
68.92
40.00
40.00
0.98
0.98
31.70
31.70
C3
111.35
111.35
111.35
87.68
87.68
0.00
0.00
1.20
1.20
48.00
48.00
C4
113.93
113.93
113.93
51.47
51.47
5.00
5.00
0.00
0.00
0.00
0.00
C5
94.50
94.50
94.50
175.59
175.59
0.00
0.00
0.10
0.10
0.00
0.00
C6
118.16
118.16
118.16
19.91
19.91
0.08
0.08
0.01
0.01
0.00
0.00
C7
144.68
144.68
144.68
12.55
12.55
7.50
7.50
0.01
0.01
0.05
0.05
C8 C9
150.66 92.50
150.66 92.50
150.66 92.50
110.59 436.34
110.59 436.34
2.00 30.00
2.00 30.00
0.25 0.09
0.25 0.09
19.20 24.00
19.20 24.00
P1
[110.45,þ¥]
[110.45,þ¥]
[110.45,þ¥]
[15, 170]
[15, 170]
[0,45]
[0,45]
[0,0.86]
[0,0.86]
[0, 35.00]
[0, 35.00]
P2
[111.95,þ¥]
[111.95,þ¥]
[111.95,þ¥]
[15, 170]
[15, 170]
[0, 50]
[0, 50]
[0,0.92]
[0,0.92]
[0, 36.00]
[0, 36.00]
P3
[108.97,þ¥]
[108.97,þ¥]
[108.97,þ¥]
[15, 170]
[15, 170]
[0,44]
[0,42]
[0,0.94]
[0,0.94]
[0, 42.00]
[0, 42.00]
P4
-
[103.24,þ¥]
[103.24,þ¥]
[15, 170]
[15, 170]
[0, 50]
[0, 50]
[0,0.90]
[0,0.90]
[0, 40.00]
[0, 40.00]
P5
-
-
[115.01,þ¥]
-
[15, 170]
-
[0, 48]
-
[0,0.93]
-
[0, 40.00]
OI
DNI
FLI
OXI
4
58
914
48
914
48
914
46
7, 8
911
1214
C1 C2
1.00 23.80
1.00 23.80
1.00 23.80
1.49 1.33
1.49 1.33
3.45 6.25
3.45 6.25
0.25 0.75
0.25 0.75
0.25 0.75
0.25 0.75
C3
0.85
0.85
0.85
1.22
1.22
2.36
2.36
2.00
2.00
2.00
2.00
C4
0.00
0.00
0.00
1.58
1.58
3.56
3.56
1.25
1.25
1.25
1.25
C5
0.40
0.40
0.40
1.50
1.50
1.96
1.96
0.08
0.08
0.08
0.08
C6
0.72
0.72
0.72
1.44
1.44
3.65
3.65
0.00
0.00
0.00
0.00
C7
0.00
0.00
0.00
1.15
1.15
2.96
2.96
0.00
0.00
0.00
0.00
C8
0.15
0.15
0.15
1.35
1.35
5.46
5.46
18.20
18.20
18.20
18.20
C9 P1
0.06 [0, 20.00]
0.06 [0, 20.00]
0.06 [0, 20.00]
1.61 [1.19, 1.67]
1.61 [1.19, 1.67]
7.95 [1.4, 7.60]
7.95 [1.4,7.60]
0.85 [0,1.85]
0.85 [0,2.80]
0.85 [0,2.80]
0.85 [0, 2.80]
P2
[0, 18.00]
[0, 18.00]
[0, 18.00]
[1.20, 1.67]
[1.20, 1.67]
[1.4, 7.25]
[1.4, 7.25]
[0,1.90]
[0,2.75]
[4, 7.25]
[0, 2.75]
P3
[0, 20.00]
[0, 20.00]
[0, 20.00]
[1.18, 1.67]
[1.18, 1.67]
[1.4, 7.20]
[1.4, 7.20]
[0,2.10]
[0,2.90]
[0, 2.90]
[0, 2.90]
P4
[0, 18.00]
[0, 18.00]
[0, 18.00]
[1.19, 1.67]
[1.19, 1.67]
[1.4, 7.50]
[1.4, 7.50]
[0,2.00]
[0,2.70]
[0, 2.70]
[0, 2.70]
P5
-
-
[0, 20.00]
-
[1.20, 1.67]
-
[1.4,7.40]
-
-
[0, 3.00]
[0, 3.00]
Then, we set R bk for all slots within each run to be the above value. Thus, if a run spans slots 36 inclusive, then we set R bk = R b6 for k = 35. Now, to obtain a realistic schedule with the constant blend rates computed above, we fix xbpk, xebk, and vb0k and follow the procedure of Li et al.8 This gives us the following reduced formulation, where [xbpk], [xebk], and [vb0k] respectively denote the optimal values of xbpk, xebk, and vb0k obtained from the full model. Equations 11a11c become BLbk ¼ 0 BLbk ¼ Tbk Tbðk1Þ
for ðb, kÞ with ½vb0 ¼ 1
ð45aÞ
Equations 12 and 13 become RLbk e RLbðk1Þ þ BLbk for ðb, kÞ with ½xebk ¼ ½vb0k ¼ 0 ð46Þ RLbk ¼ 0 for ðb, kÞ with ½vb0k ¼ 1 or ½xebk ¼ 1 and ½vb0k ¼ 0
RLbðk1Þ þ BLbk g
P
∑
p¼1
ð47Þ RLLbp ½xbpk for ðb, kÞ with ½xebk
¼ 1 and ½vb0k ¼ 0, ðb, pÞ ∈ BP, 1 e k < KÞ ð48aÞ RLbðk1Þ þ BLbk g
P
RLLbp ½xbpk ∑ p¼1
for ðb, kÞ with ½vb0k
¼ 0, ðb, pÞ ∈ BP, k ¼ KÞ
for ðb, kÞ with ½xebk ¼ ½vb0k ¼ 0 ð45bÞ
ð48bÞ
Equations 18 becomes BLbk e Tbk Tbðk1Þ
Qbk ¼ 0
for ðb, kÞ with ½xebk ¼ 1 and ½vb0k ¼ 0
Qbk ¼ Rbk 3 BLbk
ð45cÞ 9166
for ðb, kÞ with ½vb0k ¼ 1 for ðb, kÞ with ½vb0k ¼ 0
ð49aÞ ð49bÞ
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 4. Allowable Composition Ranges for Components in Products of Examples 114 C1
C2
product
1
28
914
1
P1
[0,0.22]
[0,0.22]
[0,0.22]
[0.1,1]
P2
[0,0.24]
[0,0.24]
[0,0.24]
[0.1, 1]
P3
[0,0.25]
[0,0.25]
[0,0.25]
[0.1, 1]
P4
-
[0,0.24]
[0,0.24]
P5
-
-
[0, 0.30]
28
1
28
1
28
914
1
24
57
8
912
13, 14
[0.10, 1]
[0.10, 1]
[0,1]
[0,1]
[0,1]
[0,0.40]
[0,0.40]
[0,0.40]
[0,0.40]
[0,0.40]
[0,0.40]
[0.10, 1]
[0.10, 1]
[0,1]
[0,1]
[0,1]
[0,0.45]
[0,0.45]
[0,0.45]
[0,0.45]
[0,0.45]
[0,0.45]
[0.10, 1]
[0.10, 1]
[0,1]
[0,1]
[0,1]
[0,0.43]
[0,0.43]
[0,0.43]
[0,0.43]
[0,0.43]
[0,0.43]
-
[0.10, 1]
[0.10, 1]
-
[0,1]
[0,1]
[0,0.44]
[0,0.44]
[0,0.44]
[0,0.44]
[0,0.44]
-
-
[0.15, 1]
-
-
[0,1]
-
-
-
[0,0.40]
[0,0.40]
C6
914
C4
914
C5 product
C3
1
28
914
-
C7 1
C8
28
914
1
C9
28
914
1
28
914
P1
[0, 0.25] [0, 0.25] [0, 0.25] [0, 0.20] [0, 0.20] [0, 0.20] [0, 0.25] [0, 0.25] [0, 0.25] [0, 0.30] [0, 0.30] [0, 0.30] [0, 0.15] [0, 0.15] [0, 0.15]
P2 P3
[0, 0.25] [0, 0.25] [0, 0.25] [0, 0.22] [0, 0.22] [0, 0.22] [0, 0.25] [0, 0.25] [0, 0.25] [0, 0.30] [0, 0.30] [0, 0.30] [0, 0.18] [0, 0.18] [0, 0.18] [0, 0.25] [0, 0.25] [0, 0.25] [0, 0.18] [0, 0.18] [0, 0.18] [0, 0.25] [0, 0.25] [0, 0.25] [0, 0.30] [0, 0.30] [0, 0.30] [0, 0.20] [0, 0.20] [0, 0.20]
P4
-
P5
-
[0, 0.25] [0, 0.25] -
-
[0, 0.25]
[0, 0.20] [0, 0.20]
-
-
[0, 0.20]
-
[0, 0.25] [0, 0.25] -
∑ vbjk ¼ 1 ½vb0k
j¼1
ðfor ðb, jÞ ∈ BP, 0 < j e J, 0 < k e KÞ ð50Þ
I
∑ qibk θis g
ð51aÞ ðfor ðb, pÞ ∈ BP, ðb, jÞ ∈ BJ, ðp, jÞ ∈ PJ, 0 < j e J, 0 < k e KÞ
ð51bÞ
∑
I
∑
i¼1
[0, 0.17]
Qbk θLps
Mb
ð52aÞ
ð55bÞ
θLps
minðθLps Þ p
ð1 ½xbpk Þ ð56aÞ
qibk θis e Qbk θUps þ Mb max ðθUps Þ θUps ð1 ½xbpk Þ p
I
ð56bÞ
I
ðθLps ÞgFmax ð1 ½xbpk Þ ∑ qibkFi θis g ði∑¼ 1 qibk Fi ÞθLps Mb fθLps min p i¼1
I
∑
i¼1
yibk
ðfor ðb, kÞ with ½vb0k ¼ 0, 1 e k e KÞ
-
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, 0 < k e KÞ
I
i¼1
[0, 0.16] [0, 0.16]
-
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, 0 < k e KÞ
ujpk g ½xbpk þ vbjk 1
1 ½vb0k e
-
[0, 0.30]
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, 0 < k e KÞ
ðfor ðb, pÞ ∈ BP, ðb, jÞ ∈ BJ, ðp, jÞ ∈ JP, 0 < j e J, 0 < k e KÞ
yibk e 1 ½vb0k ðfor ðb, kÞ with ½vb0k ¼ 1, 1 e k e KÞ
-
ðfor ðb, kÞ with ½xebk ¼ 1, 1 e k e KÞ
i¼1
½xbpk g ujpk þ vbjk 1
[0, 0.30] [0, 0.30]
-
Tbðk1Þ þ BLbk g tik Hð1 yibk Þ
Fixing the values of xbpk, xebk, and vb0k, eqs 2, 5, 10, 14a, 14c, 1921, 33a, 33b, 34a, 34b, and 40 become J
-
[0, 0.25]
I
∑
qibk Fi θis e
i¼1
! qibk Fi θUps þ Mb maxðθUps Þ θUps Fmax ð1 ½xbpk Þ p
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, 0 < k e KÞ
ð52bÞ
ð57aÞ
ð57bÞ
Qbk rpiL Mb frpiL minðrpiL Þgð1 ½xbpk Þ e qibk e Qbk rpiU p
STbk ¼ 0 ðfor ðb, kÞ with ½xebk ¼ 0 and ½vb0k ¼ 0, 1 e k e KÞ ð53aÞ
þ Mb fmax ðrpiU Þ rpiU gð1 ½xbpk Þ p ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, 0 < k e KÞ ð58a; bÞ
STbk g τb ðfor ðb, kÞ with ½xebk ¼ 1 and ½vb0ðkþ1Þ ¼ 0, 1 e k < KÞ ð53cÞ
TC ¼
Tbk g Tik H½1 yibðkþ1Þ ðfor ðb, kÞ with ½xebk ¼ 1, 0 e k < KÞ
þ ð54aÞ
Tbk e Tik þ H½1 yibðkþ1Þ ðfor ðb, kÞ with ½xebk ¼ 1, 0 e k < KÞ
ð54bÞ
Tbðk1Þ þ BLbk e tik þ Hð1 yibk Þ ðfor ðb, kÞ with ½xebk ¼ 1, 1 e k e KÞ
I
ð55aÞ
B
K
B
K1
∑ ∑ ∑ ci qibk þ b∑¼ 1 k∑¼ 1 CBb ½xebk i¼1 b¼1 k¼1 J
K1
O
∑ ∑ CTjuejk þ o∑¼ 1 DMo do j¼1 k¼1
ð59Þ
The revised model (RSPM) comprises eqs 1, 3, 8, 14b, 1517, 2232, 33c, 33d, 34c, 34d, 3539, and 4159, whose solution will guarantee constant blend rates with applicable limits and minimum run lengths. As noted by Li et al.,8 the schedule from RSPM may still show intermittent delivery of orders. If a delivery spans contiguous slots, then we can revise the schedule to simply deliver the 9167
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 5. Blender and Economic Data for Examples 114 RLb,0, h and CQb,0 ( 103 bbl) blender
14
5
6, 7
812
Minimum and Maximum Blending Rate ( 103 bbl/h)
13, 14
1
2
36
7
8, 9
1012
Allowable Product
13, 14
1
27
8
912 13, 14
B-101 0 and 0 10 and 150 0 and 0 0 and 0 0 and 0 2.015 1.515 1.520 1.525 1.525 1.530 1.530 P1P3 P1P4 P1P4 P1P5 P1P5 B-102
-
-
-
B-103
-
-
-
0 and 0 0 and 0 -
0 and 0
-
-
-
-
-
-
-
-
1.525 1.530 1.530 -
-
1.525
-
-
-
-
P1P4 P1P5 P1P5 -
-
P1P5
Minimum Blend Run Length (h) P1 17
P2
812
P3
13, 14
17
812
13
14
17
P4
812
13, 14
27
P5
811
12
13, 14
911
12
13, 14
B-101
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
5
5
B-102
-
6
6
-
6
6
6
-
6
6
-
6
6
6
6
5
5
B-103
-
6
-
-
6
6
-
-
6
-
-
-
6
-
-
5
Cost ($/bbl) Example
C1
C2
C3
C4
C5
C6
C7
C8
C9
Example
demurrage cost ( 103 $/h)
scheduling horizon (h)
1 214
20 20
40 24
55 30
45 25
25 22
50 27
58 50
60 50
30 22.5
1 214
2.1 2.5
72 192
Transition Cost ( 103 $/instance) Example PT-101 PT-102 PT-103 PT-104 PT-105 PT-106 PT-107 PT-108 PT-109 PT-110 PT-111 transition cost in blender ( 103 $/instance) 1
9.8
9.8
14.5
9.8
98
-
-
-
-
-
214
14.5
14.5
19
19
19
14.5
19
14.5
19
14.5
complete order at a constant rate.8 Figure 4 shows the complete schedule adjustment procedure.
5. MULTIPERIOD EXTENSION Following the approach of Li et al.,8 we modify eqs 19 and 20 to ensure on-spec products for the piecewise-constant profiles of the multiperiod scenario. I
∑ qibk θist g
i¼1
Qbk θLps
Mb fθLps
minðθLps Þgð1 xbpk Þ p
ðfor ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ
ð60aÞ
I
ðθUps Þ θUps gð1 xbpk Þ ∑ qibk θist e Qbk θUps þ Mb fmax p i¼1 ðfor ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ I
ð60bÞ
ðfor ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ
∑
i¼1
20
Furthermore, eqs 36 and 38 change as follows: VCik ¼ VCiðk1Þ þ
B
∑ Fit ½Tik Tiðk1Þ b∑¼ 1 qibk t ∈ TK ðfor 1 e k e KÞ
ViL
e VCiðk1Þ þ
∑
t ∈ TK
Fit ½tik Tiðk1Þ
B
∑ qibk e
b¼1
ðfor 1 e k e KÞ
qibk Fit θist e ð
I
∑
i¼1
ð61aÞ
qibk Fit ÞθUps þ Mb fmax ðθUps Þ θUps gFmax t ð1 xbpk Þ p
ðfor ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ
ð61bÞ
where TK = {(t, k) | slot k is in period t}, θist denotes the known blending index for a property s of component i during period t, Fit is the density of component i during period t, and Fmax is the t maximum possible density among all products during period t.
ð62Þ
ViU ð63a; bÞ
where Fit is the flow rate of component i in period t. With this, our multiperiod model (MPM) comprises eqs 118, 2135, 37, 3940, and 6063. For revised MPM (RMPM), we only need to change eqs 60 and 61, as follows: I ∑ qibkθist g Qbk θLps Mb θLps minðθLps Þ ð1 ½xbpk Þ p
i¼1
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ
ð64aÞ
I
ðθLps ÞgFmax ∑ qibk Fit θist g ði∑¼ 1 qibk Fit ÞθLps Mb fθLps min t ð1 xbpk Þ p i¼1
I
14.5
I
∑
i¼1
qibk θist e Qbk θUps þ Mb fmax ðθUps Þ θUps gð1 ½xbpk Þ p
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ I
I
∑ qibk Fit θist g i∑¼ 1 qibk Fit
i¼1
!
ð64bÞ
θLps Mb θLps minðθLps Þ Fmax t ð1 ½xbpk Þ p
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ ð65aÞ ! I I qibk Fit θist e qibk Fit θUps þ Mb max ðθUps Þ θUps Fmax t ð1 ½xbpk Þ
∑
i¼1
∑
i¼1
p
ðfor ðb, kÞ with ½vb0k ¼ 0, ðb, pÞ ∈ BP, ðt, kÞ ∈ TK, 0 < k e KÞ ð65bÞ 9168
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 6. Periods, Durations, Slots, and Feed Rates to Component Tanks for Examples 114 Feed Flow Rate to Component Tank ( 103 bbl/h) example
period
duration
slot
CT-101
CT-102
CT-103
CT-104
CT-105
CT-106
CT-107
CT-108
CT-109
1
1
40
1, 2
1.2
0.8
1.2
1.2
0.5
0.8
0.0
0.0
1.0
2
32
3, 4
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
1
60
1, 2
1.2
0.8
1.2
1.2
0.5
0.8
0.0
0.0
1.0
2
132
3
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
1
100
1, 2
1.2
0.8
1.2
1.2
0.5
0.8
0.0
0.0
1.0
2
92
3
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
5
1
100
13
1.2
0.8
1.2
1.2
0.5
0.8
0.0
0.0
1.0
6
2 1
92 120
4, 5 1, 2
0.8 1.2
0.6 0.8
0.6 1.2
0.8 1.2
0.5 0.7
0.6 0.8
0.5 0.0
0.5 0.0
0.0 1.0
2
72
3, 4
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
1
80
13
1.2
0.8
1.2
1.2
0.7
0.8
0.0
0.0
1.0
2
70
46
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3
42
7
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
1
80
1, 2
1.2
0.8
1.2
1.2
0.5
0.8
0.0
0.0
1.0
2
70
3, 4
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3 1
42 80
5, 6 1, 2
0.0 1.0
0.0 0.5
0.0 1.0
0.0 1.0
0.0 0.5
0.0 0.5
0.0 0.0
0.0 0.0
0.0 1.0
2
70
3, 4
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3
42
5, 6
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
1
80
1, 2
1.0
0.5
1.0
1.0
0.8
0.5
0.0
0.0
1.0
2
60
35
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3
52
6, 7
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
1
80
13
1.0
0.5
1.0
1.0
0.7
0.5
0.0
0.0
1.0
2 3
60 52
46 7
0.8 0.5
0.6 0.5
0.6 0.5
0.8 0.5
0.5 0.5
0.6 0.5
0.5 0.0
0.5 0.0
0.0 0.5
1
50
13
1.0
0.5
1.0
1.0
0.8
0.5
0.0
0.0
1.0
2
50
46
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3
50
79
0.5
0.5
0.5
0.5
0.5
0.5
0.0
0.0
0.5
4
42
1012
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
1
50
13
1.0
0.5
1.0
1.0
0.7
0.5
0.5
0.5
1.0
2
50
47
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3 4
50 42
810 1113
0.5 0.0
0.5 0.0
0.5 0.0
0.5 0.0
0.5 0.0
0.5 0.0
0.0 0.0
0.0 0.0
0.5 0.0
1
50
13
1.0
0.5
1.0
1.0
0.7
0.5
0.5
0.5
1.0
2
50
47
0.8
0.6
0.6
0.8
0.5
0.6
0.5
0.5
0.0
3
50
810
0.5
0.5
0.5
0.5
0.5
0.5
0.0
0.0
0.5
4
42
1113
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
2 3, 4
7
8
9
10
11
12
13
14
Thus, RMPM comprises eqs 1, 3, 8, 14b, 1517, 2232, 33c, 33d, 34c, 34d, 35, 37, 39, 4155, 58, 59, and 6265. So far, we do not allow a product tank to receive and deliver simultaneously. Next, we address the aforementioned third limitation.
As we did for xebk and uejk, we write the following to compute and force zejok to be binary:
6. SIMULTANEOUS RECEIPT/DELIVERY BY PRODUCT TANKS If a product tank can receive and deliver simultaneously, then we need one more 01 continuous variable (zejok) to denote the end of the current delivery run as follows:
zejok g zjoðkþ1Þ zjok
(
zejok ¼
1 if tank j ends current delivery run of order o in slot k 0 ðotherwiseÞ
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k < KÞ
zejok g zjok zjoðkþ1Þ
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k < KÞ
ð66aÞ ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k < KÞ
ð66bÞ zejok þ zjok þ zjoðkþ1Þ e 2
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k < KÞ
ð67Þ If tank j does not end delivering order o in slot k, then it must continue in (k þ 1). In other words, the delivery in slot k must end exactly at Tjk. Similarly, if the delivery during slot k does not end, then the delivery in slot (kþ1) must begin 9169
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 7. Solution Statistics of Various Algorithms/Codes for SPM for Examples 114a order Example 1 Example 2 Example 3
5 10 12
number of
number of
number of
cost
total CPU
number of
number of slots
discrete variables
continuous variables
constraints
( 103 $)
time (s)
switching operations
unit slot
4
136
429
1577
5149.73
0.89
2
process slotc
4
100
269
1098
5149.73
0.72
2
unit slot
3
202
733
2743
3658.11
2.12
1
process slotc
4
252
598
2847
3658.11
6.7
1
model
unit slot
3
246
799
3047
3159.12
process slotc
4
315
633
3199
3159.12
18.8
2.62
1
1
3 4
285 372
925 723
3551 3823
4556.67 4556.67
4.20 39.4
1 1
Example 4
15
unit slot process slotc
Example 5
15
unit slot
3
285
925
3563
4556.67
2.41
1
process slotc
4
372
723
3867
4556.67
6.94
1
Example 6
18
unit slot
4
476
1320
5372
5213.88
3024
process slotc
5
576
950
5438
5213.88
10801
2 2
unit slot
5
672
1745
7323
8100.35
10802
3
Example 7
20
process slotc
6
775
1197
7175
8100.35
10801
3
Example 8
20
unit slot process slotc
4 5
595 682
1548 1164
7327 7664
8080.35 8080.35
10819 10801
2 2
Example 9
23
unit slot
4
635
1651
8083
10573.65
10814
4
process slotc
8
1251
1890
13922
10640.05
14400b
6 4
Example 10 Example 11
25 30
unit slot
7
1265
2913
15001
11286.10
14400b
process slotc
9
1503
2209
16492
11326.10
10817
6
unit slot
8
1644
3729
19679
13248.58
46800b
4
process slotc
12
2328
3216
24755
13437.25
46800b
10 5 13
Example 12
35
unit slot process slotc
9 16
2237 3846
4994 4963
26135 39520
14764.86 15475.34
46800b 46800b
Example 13
40
unit slot
12
3378
7222
43139
15646.15
118800b
6
process slot
17
4414
5966
51775
15879.39
108282
13
unit slot
12
3664
7877
46529
17737.39
118800b
9
17
4830
6439
57657
18249.07
118800b
13
c
Example 14
45
c
process slot a
Note: CPU time limit for SPM is set at 10 800 s for Examples 110, 36000 s for Examples 11 and 12, 108000 s for Examples 13 and 14. CPU time limit for RSPM is set at 3600 s for Examples 110, and 10800 s for Examples 1114. b Reached total CPU time limit. c Process-slot from Li et al.8
tsjok e Tbðk1Þ þ Hð2 vbjk zjok Þ ðfor ðj, oÞ ∈ JO, 1 e k e KÞ
exactly at Tjk. Equation 31, along with the following expressions, ensures this. tsjok þ
DQjok g Tjk H½1 zjok þ zejok DRjo
ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
tsjok þ ð68aÞ
tsjoðkþ1Þ e Tjk þ H½1 zjoðkþ1Þ þ zejok ðfor ðj, oÞ ∈ JO, 1 e j e J, 1 e k e KÞ
ð68bÞ
When a product tank can receive and deliver products at the same time, the flows into/out of the tank can occur in the same slot. Then, eq 36 cannot ensure a correct inventory profile. We need additional constraints. Since we define tsjok as the start of delivery of o by tank j in slot k and the end of delivery is given by tsjok þ DQjok/DRjo. Now, if a blender b is feeding a tank j during its slot k, and tank j is delivering order o at the same time, then the start/end times of delivery must match the start/end times of blender during slot k. This gives us the following additional timing constraints for product tank j: tsjok g Tbðk1Þ Hð2 vbjk zjok Þ ðfor ðj, oÞ ∈ JO, 1 e k e KÞ
tsjok þ
ð69aÞ
DQjok g Tbðk1Þ þ BLbk Hð2 vbjk zjok Þ DRjo ðfor ðj, oÞ ∈ JO, 1 e k e KÞ
ð69bÞ
ð70aÞ
DQjok e Tbðk1Þ þ BLbk þ Hð2 vbjk zjok Þ DRjo ðfor j, oÞ ∈ JO, 1 e k e KÞ ð70bÞ
Thus, the model equations for the simultaneous receipt and delivery by product tanks are SPM: eqs 122, 2440, and 6670. RSPM: eqs 1, 3, 8, 14b, 1517, 22, 2432, 33c,33d, 34c34d, 3539, 4159, and 6670. MPM: eqs 118, 21, 22, 2435, 37, 39, 40, 6063, and 6670. RMPM: eqs 1, 3, 8, 14b, 1517, 22, 2432, 33c, 33d, 34c, 34d, 35, 37, 39, 4155, 58, 59, 6265, and 6670. Note that allowing simultaneous receipt and delivery increases the number of slots considerably, and it makes the problem much more difficult. 9170
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Table 8. Solution Statistics of Various Algorithms/Codes for MPM for Examples 114a order Example 1 Example 2 Example 3
5 10 12
Example 4
15
Example 5
15
Example 6 Example 7
18 20
Example 8
20
Example 9
23
Example 10
25
model
30
number of
number of
number of
total CPU
cost
number of
slots
discrete variables
continuous variables
constraints
time (s)
( 103 $)
switching operations
unit slot
4
136
429
1703
1.25
5149.73
2
process slotc
5
130
328
1384
6.89
5149.73
2
unit slot
3
202
733
2683
process slotc
5
329
729
3600
0.69 50.7 1.16
3658.11
1
3658.11
1
unit slot
3
246
799
2987
3159.12
1
process slotc
5
411
769
4048
80.0
3179.12
2
unit slot process slotc
3 5
285 486
925 877
3479 4880
1.97 170
4556.67 4556.67
1 1
unit slot
5
531
1461
6121
46.1
4556.67
1
process slotc
5
486
877
4892
56.9
4556.67
1
unit slot
4
476
1320
5404
1037
5213.88
2
process slotc
6
712
1114
6631
10801
5213.88
2
unit slot
7
986
2381
10337
10814
8100.35
3
process slotc
9
1219
1725
10922
14400b
8682.34
3
unit slot process slotc
6 8
949 1159
2244 1779
11491 12464
14170 10831
8082.85 8123.90
2 4
unit slot
6
1015
2391
12635
10817
10573.65
4
process slotc
9
1423
2107
15708
10833
10636.40
5
unit slot process slotc
Example 11
number of
unit slot
7
1265
2913
15469
11735
11306.10
5
10
1685
2436
18368
10802
11349.54
6
7
1419
3289
17179
36075
13248.58
4
process slotc
14
2742
3720
28957
42342
13348.88
5
Example 12
35
unit slot process slotc
12 17
3050 4099
6572 5258
34151 42026
46800 36467
14809.36 15316.24
7 9
Example 13
40
unit slot
13
3677
7801
46771
118800b
15626.15
5
process slot
18
4686
6302
54857
108005
16187.58
5
unit slot
13
3989
8508
49097
118800b
17737.39
9
process slotc
17
4830
6439
57657
118800b
18678.48
5
c
Example 14
45
a
Note: CPU time limit for SPM is set at 10800 s for Examples 110, 36000 s for Examples 11 and 12, 108000 s for Examples 13 and 14. CPU time limit for RSPM is set at 3600 s for Examples 110, and 10800 s for Examples 1114. b Reached total CPU time limit. c Process-slot-based methodology from Li et al.8
7. NUMERICAL EVALUATION Tables 16 show the data for 14 examples (Examples 114) from Li et al.8 These examples vary in structure, size, scale, and complexity to mimic real-life industry scenarios. Example 1 involves nine components, five product tanks, and a scheduling horizon of 72 h (3 days). Examples 214 have 9 components, 11 product tanks, and a horizon of 192 h (8 days). Examples 17 have one blender, Examples 812 have two, and Examples 13 and 14 have three blenders. Example 1 has three products, Examples 28 have four products, Examples 913 have five products, and Example 14 has six products. The number of orders varies from 10 to 45 in various examples. Li et al.8 used nine product properties: octane number (ON), Reid vapor pressure index (RVPI), sulfur index (SI), benzene index (BI), aromatics index (AI), olefin index (OI), specific gravity index (SGI), flammability index (FI), and oxygenates index (OXI). We considered one property index in Example 1, five in Examples 2 and 3, and nine in Examples 414. The blenders are idle at time zero except in Example 5, where the blender is processing product 1. We used CPLEX 12.1.0 on Dell Precision PWS690 computer (Intel Xeon CPU 3.00 GHz, 16.0 GB memory) running Windows XP to solve the 14 examples.
Tables 7 and 8 present the model statistics for unit-slot and process-slot-based SPM and MPM. Unit-slot-based SPM (USSPM) requires fewer slots in Examples 214, compared to process-slot-based SPM (PSSPM). For instance, Example 12 requires 9 unit slots, compared to 16 process slots. While both USSPM and PSSPM attain the optimal solution in Examples 18, the former solves faster than the latter. For instance, USSPM requires 4.20 CPU s in Example 4 versus 39.4 CPU s for PSSPM. In Example 6, PSSPM fails to prove the optimal solution, even after 3 days of CPU time. While USSPM takes 3024 s to give a solution of $5213.88K, PSSPM requires 10 801 CPU s for the same. Examples 714 are unsolvable to optimality by both models within the allocated CPU times in Table 7. However, USSPM obtains better solutions than PSSPM in all cases. For instance, USSPM gives a cost of $10 573.65K, compared to $10 640.05K for PSSPM in Example 9. In Example 14, USSPM gives $17 737.39K, compared to $18 249.07K for PSSPM. In addition, we also observe that the solutions from USSPM involve, at most, nine product transitions (i.e., on blenders and in product tanks). However, USSPM always gives the same or fewer switches than PSSPM for all examples. For instance, USSPM gives four transitions for Example 10, whereas PSSPM gives six. Similar observations hold true for MPM models 9171
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
Figure 5. (a) Blending schedule for Example 12 (35 orders) from RSPM. (b) Order delivery schedule for Example 12 (35 orders) from RSPM.
as well, except for Example 14, where unit-slot-based (USSPM) gives nine transitions, compared to five for process-slot-based MPM (PSMPM). For illustration, Figure 5 gives the final schedule from SPM for Example 12 (35 orders).
8. CONCLUSIONS We developed a novel multigrid continuous-time formulation for addressing recipe, blending, storage, and order
scheduling simultaneously. Our new formulation improves on our previous work 8 in several respects. By using unit slots instead of process slots, it is able to solve larger problems and is computationally more efficient. In addition, it avoids the assumption of unlimited component inventories, allows blender setup times and simultaneous receipt/delivery by tanks, and revises the schedule adjustment procedure of Li et al.8 for unit slots. 9172
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
to acknowledge Prof. Vladimir Mahalec (McMaster University) for his valuable suggestions on some practical features of real blending systems, which we have tried to include in this paper.
’ NOTATION Sets
BP = set of pairs (blender b, product p) such that b can process p BJ = set of pairs (blender b, product tank j) such that b can feed j PJ = set of pairs (product p, product tank j) such that j can hold p JO = set of pair (product tank j, order o) such that j can deliver o Subscripts
i = component and its tank b = blender p = product o = order j = product tank s = gasoline property specification k = slot t = period Superscripts Figure A1. An example to indicate why, if a blender begins a blend run in a slot, each tank feeding that blender must have a matching slot with the same index and start time. It also illustrates why a check on the tank inventory is needed at the end of the blend run.
’ APPENDIX Consider the example in Figure A1, where CT-101 has a capacity of 1.4 kbbl and feeds B-101 only. The feed rate to CT-101 from upstream units is 1.5 kbbl/h, whereas the consumption rate of B-101 is 3.5 kbbl/h. If BT-101 withdraws from CT-101 during [4, 5] h, then the inventory levels are 0.0 kbbl at 3 h, 1.0 kbbl at 6 h, and 1.5 kbbl at 4 h, and 0.5 kbbl at 5 h. Clearly, we must compute and check the inventory levels at 3, 4, 5, and 6 h. However, this cannot be done by using [3,6] h as one slot on CT-101 and [4,5] h as another on B-101 (Figure A1a). The duration [3, 6] h must be broken up into several slots. [3, 4] h can be merged into the preceding slot. A slot must begin on both B-101 and CT-101 at 4 h. Finally, [4, 6] h can be considered a single slot on both B-101 and CT-101, and we must check the inventory level at 5 h using a continuous variable, as illustrated in Figure A1b. This will guarantee no storage violation for CT-101. From this example, we conclude that whenever a blender begins a blend run in a slot, any component tank feeding that blender must have a matching slot with the same index and the same start time. Moreover, we must check the inventory level at the end of each blend run. ’ AUTHOR INFORMATION Corresponding Author
*Tel.: þ65 6516-6359. Fax: þ65 6779-1936. E-mail: cheiak@ nus.edu.sg.
’ ACKNOWLEDGMENT We would like to acknowledge financial support for this work from The Agency for Science, Technology, and Research (A*Star), under Grant No. 052 116 0074. We would also like
U = upper limit L = lower limit Parameters
Ci = price of component i ($ per unit volume) CBb = cost a transition on blender b ($ per occurrence) CTj = cost of transition in product tank j ($ per occurrence) DDLo = start delivery time of order o DDU o = due date of order o DMo = demurrage cost ($ per unit time) of order o DRjo = delivery rate of product tank j to order o fttc = constant feed rate of component i from upstream units in period t = limits on blending rate of blender b FL/U b Fbk = blending rate of blender b in slot k Fi = constant feed rate of component i from upstream units = limits on the fraction of component i in product p rL/U pi RLLbp = minimum run length of blender b for product p TCQbk = volume of blender b processed in slot k TQo = total amount of order o = capacity limits of component tank i VCL/U i VPU j = maximum capacity of product tank j θL/U ps = limits on property indices of product p θis = blending index for property s of component i Fi = density of component i Fmax = maximum density among all products τb = setup time for blender b Binary Variables
ujpk = 1 if product tank j holds product p during slot k vbjk = 1 if blender b feeds product tank j during its slot k yibk = 1 if blender b receives from component tank i during its slot k zjok = 1 if product tank j is delivering order o during slot k 01 Continuous Variables
uejk = 1 if product tank j switches products at the end of its slot k xbpk = 1 if blender b processes product p during its slot k xebk = 1 if blender b ends its current blend run during its slot k 9173
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174
Industrial & Engineering Chemistry Research
ARTICLE
zejok = 1 if product tank j ends current delivery run of order o during its slot k Continuous Variables
CQbk = volume processed of blender b if it does not end its run during slot k do = demurrage for order o DQjok = volume of order o delivered by product tank j during slot k Gbjk = volume that blender b feeds product tank j during slot k qibk = volume of component i used by blender b during slot k Qbk = volume processed in blender b during slot k RLbk = length of the current run of blender b at the end of slot k STbk = accumulated idle time of a blender b by the end of slot k Tnk = end time of slot k on unit n, where n becomes i for a component tank, b for a blender, and j for a product tank TC = total operating cost ($) VCik = inventory of component tank i at the end of slot k VPjk = inventory of product tank j at the end of slot k
’ REFERENCES (1) Pinto, J. M.; Joly, M.; Moro, L. F. L. Planning and Scheduling Models for Refinery Operations. Comput. Chem. Eng. 2000, 24, 2259–2276. (2) Li, J.; Li, W. K.; Karimi, I. A.; Srinivasan, R. Improving the Robustness and Efficiency of Crude Scheduling Algorithms. AIChE J. 2007, 53, 2659–2680. (3) Dewitt, C. W; Lasdon, L. S.; Waren, A. D.; Brenner, D. A.; Melhem, S. A. OMEGA: An Improved Gasoline Blending System for Texaco. Interfaces 1989, 19, 85–101. (4) Rigby, B.; Lasdon, L. S.; Waren, A. D. The Evolution of Texaco’s Blending Systems: From OMEGA to Starblend. Interfaces 1995, 25, 64–83. (5) Glismann, K.; Gruhn, G. Short-Term Scheduling and Recipe Optimization of Blending Processes. Comput. Chem. Eng. 2001, 25, 627–634. (6) Joly, M; Pinto, J. M. Mixed-Integer Programming Techniques for the Scheduling of Fuel Oil and Asphalt Production. Trans. Inst. Chem. Eng., Part A 2003, 81, 427–447. (7) Jia, Z. Y.; Ierapetritou, M. Mixed-Integer Linear Programming Model for Gasoline Blending and Distribution Scheduling. Ind. Eng. Chem. Res. 2003, 42, 825–835. (8) Li, J.; Karimi, I. A.; Srinivasan, R. Recipe Determination and Scheduling of Gasoline Blending Operations. AIChE J. 2010, 56, 441–465. (9) Mendez, C. A.; Grossmann, I. E.; Harjunkoski, I.; Kabore, P. A. Simultaneous Optimization Approach for Off-Line Blending and Scheduling of Oil-Refinery Operations. Comput. Chem. Eng. 2006, 30, 614–634. (10) Liu, Y.; Karimi, I. A. Scheduling Multistage, Multiproduct Batch Plants with Nonidentical Parallel Units and Unlimited Intermediate Storage. Chem. Eng. Sci. 2007, 62, 1549–1566. (11) Susarla, N.; Li, J.; Karimi, I. A. A New Approach to Scheduling Multipurpose Batch Plants using Unit Slots. AIChE J. 2010, 56, 1859–1879.
9174
dx.doi.org/10.1021/ie102321b |Ind. Eng. Chem. Res. 2011, 50, 9156–9174