ARTICLE pubs.acs.org/IECR
Dynamic Scheduling for Ethylene Cracking Furnace System Chuanyu Zhao, Chaowei Liu, and Qiang Xu* Dan F. Smith Department of Chemical Engineering, Lamar University, Beaumont, Texas 77710, United States ABSTRACT: The cracking furnace system is crucial for an olefin plant. Its operation needs to follow a predefined schedule to process various feeds continuously, meanwhile conducting a periodically decoking operation for each furnace when its performance apparently decreases. In practice, because the feed supply changes dynamically, the routine furnace scheduling is better performed in a dynamic and reactive way, through which the furnace operations can be smartly rescheduled with respect to any delivery of new coming feeds. Thus, the feeds from the new delivery and the leftover inventories can be timely, feasibly, and optimally allocated to different furnaces for processing to obtain the maximum average net profit per time. Facing this challenge, this paper develops a new MINLP-based reactive scheduling strategy, which can dynamically generate reschedules based on the new feed deliveries, the leftover feeds, and current furnace operating conditions. It simultaneously addresses all the major scheduling issues of a cracking furnace system, such as semicontinuous operation, nonsimultaneous decoking, secondary ethane cracking, and seamless rescheduling. The efficacy of the study and its significant economic potential are demonstrated by a comprehensive case study.
1. INTRODUCTION Ethylene and propylene are the most widely produced organic compounds in the world and are critically important for the chemical/petrochemical industries and daily life. Because of the rapid growth of the worldwide economy and population, the global consumptions of ethylene and propylene in 2008 were approximately 113 million metric tons and 73 million metric tons, respectively. It is predicted that ethylene consumption will grow by an average of 3.7% per year from 2008 to 2013, and that propylene consumption will increase by about 4.9% per year during the same period.1 To meet the increasing global demand for ethylene and propylene, technology advancement for improving olefin plant productions is being enthusiastically pursued. Certainly, the technological innovation for the cracking furnace operation is one of the most important tasks. The cracking furnace system is crucial for an olefin plant because major product and intermediate yields are mainly determined by its pyrolysis operation, and it consumes about half the entire plant energy supplies. The performance of a cracking furnace depends on its design, feed characteristics, and operation. In general, the unique operational features of a cracking furnace system include the following: (i) The cracking furnace operation is a semicontinuous and performance-decaying exercise. This is because the coke is formed, deposited, and accumulated on the reaction coils during thermal cracking, which increases its heat-transfer resistance and pressure drop, resulting in lower energy efficiency, reaction selectivity, and productivity. Thus, the product yield will dynamically change with respect to time, and a cracking furnace must be shut down for decoking (or cleanup) after continuously running for some time so that the performance of cracking furnaces can be restored. Usually, the cleanup takes 13 days. The run length between two adjacent cleanups, the so-called batch processing time, is normally 3090 days. Figure 1 gives an illustrative example of the related concepts. (ii) When multiple furnaces are simultaneously employed for production, the decoking operation is subject to an important r 2011 American Chemical Society
constraint that simultaneous furnace cleanups are generally not allowed. This is because an olefin plant usually has only one decoking facility, which allows cleaning one furnace at a time; meanwhile, simultaneously shutting down two or more cracking furnaces for cleanup causes significant upsets to the downstream process.2 Therefore, this constraint requires that any two cleanups for any two furnaces should have no time overlap (see Figure 1). (iii) The cracking operation produces multiple products. Thus, the furnace economic performance should consider all of them. Furthermore, the pyrolysis mixture contains an amount of ethane. The industrial practice for handling ethane is to recycle it as an internal feed and conduct the secondary cracking in a specific furnace. Therefore, an ethane furnace is usually designated in an olefin plant for specifically cracking the recycled ethane feed. Because of the importance of the cracking furnace, many studies have been conducted for improving cracking furnace operations. However, the major focuses are in the areas of pyrolysis simulation, control, and operational optimization of a single furnace.310 When multiple furnaces are employed to process multiple feeds with the consideration of various product prices and manufacturing costs, the scheduling for the entire furnace system, i.e., optimally coordinating information elements of feed, furnace, time, quantity, and sequence for the maximum net profit, becomes critically important for an olefin plant. The study of production scheduling for a cracking furnace system is still lacking because the thermal cracking process involves highly complex reaction kinetics. It is very hard to integrate such reaction complexity into the scheduling activity for the entire furnace system. For the simplified multifurnace scheduling, Jain and Grossmann initiated an MINLP (mixedinteger nonlinear programming) model for the cyclic scheduling Received: February 14, 2011 Accepted: September 21, 2011 Revised: August 25, 2011 Published: September 21, 2011 12026
dx.doi.org/10.1021/ie200318p | Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Figure 1. Illustrative example for the dynamic scheduling.
of multiple feeds cracked by multiple furnaces and developed a solving algorithm for the global optimality.11 However, the developed MINLP model did not give considerations to the nonsimultaneous cleanups and the secondary ethane cracking. Schulz et al. studied the scheduling for ethane-fed ethylene plants with the consideration of recycled ethane by taking into account the interactions between the entire plant operation and furnace performance.12 Based on that, they developed an MINLP model to study the optimal performance of cyclic furnace shutdowns and downstream separation systems.13 These models could capture the decaying performance throughout furnace operations and were claimed to increase the overall plant net profit. However, they still allowed simultaneous cleanup scenarios and only addressed the cracking of one type of feed. In 2006, Lim et al. developed a neural network based scheduling model, which employed online information of ethylene and propylene yields, coil skin temperature, and pressure drop to support the scheduling decisions.14 To tackle the computational complexity of the developed large-scale MINLP model, they also developed three alternative solution strategies. Unfortunately, this study did not consider the situation of cracking multiple feeds. To further handle the unexpected uncertainties, they developed an MINLP scheduling model to determine the optimal decoking policy.15 The scheduling activity is conducted in a dynamic scheme, and it is triggered when the gaps between the model predictions and the measurements of the tube skin temperature and the pressure drop exceed the threshold values. The model can determine appropriate rescheduling points before actual operational problems arise. Recently, Liu et al. developed a cyclic scheduling model by considering multiple feeds, multiple cracking furnaces, and nonsimultaneous cleanup constraints.2 However, the model did not consider the issue of the secondary ethane cracking. To overcome this deficiency, Zhao et al. developed a new MINLP model, which enhanced the applicability of the scheduling results.16 Based on the literature survey, it is noticed that dynamic feed supply issues for the cracking furnace scheduling are completely
neglected. As known, in order to enhance the production flexibility of an olefin plant, multiple cracking feeds are adopted by olefin plants. Because the feedstock market is volatile, there are dynamic changes in the feed supply, such as the changes of shipping time and quantity of the new coming feeds. Under such situations, it is better to perform the routine furnace scheduling in a dynamic and reactive way with respect to any finalized new coming feeds. The furnace operations should be dynamically rescheduled once the supply scenarios are determined so that all the feeds of the new delivery and the leftover inventories can be timely, feasibly, and optimally allocated to different furnaces for processing to obtain the maximum average net profit per time (see Figure 1). Studies on reactive rescheduling methods and applications have been extensively reported,1721 where decisions must be dynamically adjusted when updated changes of uncertainties/input scenarios are available. They are mainly used in batch processes instead of semicontinuous and performancedecaying processes such as cracking furnace systems. Conceivably, the handling of such a dynamic scheduling problem presents a big challenge. Facing this challenge, this paper develops a new MINLP-based reactive scheduling strategy, which can dynamically and seamlessly generate reschedules based on the new feed deliveries, the leftover feeds, and the current furnace operation conditions with every rescheduling initialization. This is an important improvement compared with previous studies, because the scheduling is never cyclic and can be dynamically rescheduled at any time with respect to a feed supply change. It simultaneously addresses all the major scheduling issues of a cracking furnace system, such as semicontinuous operation, nonsimultaneous decoking, secondary ethane cracking, and seamless rescheduling. The efficacy of the study and its significant economic potential are demonstrated by a case study.
2. DYNAMIC SCHEDULING FRAMEWORK An ethylene cracking furnace system consists of several parallel furnaces for cracking multiple feeds and one ethane 12027
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Figure 2. Scope of the studied cracking furnace system.
Figure 3. Methodology framework.
furnace specifically used for cracking the recycled ethane. The scope of the system is shown in Figure 2. In the system, every feed including the secondary ethane produces multiple products (e.g., ethylene, propylene, ethane, and other products), which will be accounted for in the operating profit. The new purchased feeds are delivered to the olefin plant in batches. However, the arrival date, feed type, and quantity are dynamically changing. When a
supply scenario is determined, the plant will reactively reschedule the furnace system operation based on the new coming feeds and the leftover inventories. The objective of each schedule is to maximize the average net profit per time for the furnace system operation. Figure 3 shows the framework of the developed dynamic scheduling strategy. First, an MINLP-based model is developed for each rescheduling activity. With the input of the initial data, the model is solved by GAMS22 to obtain the solution of the current scheduling problem. Then, the solution feasibility will be checked before the implementation. If it is satisfied, the current schedule will be continuously executed until the delivery of the coming new feeds has been determined. When a delivery scenario is determined, the new feed types and amounts are also known. The data plus the current leftover inventory and the current furnace operation conditions will be utilized to generate the initial conditions for rescheduling. Note that the delivery may occur when a furnace is still in cracking operation. If the new schedule starts right after the new feed delivery, it could cause unfeasible solutions such as one furnace cracking two feeds in just one batch processing, or induce significant upsets for the cracking furnace operation. To smoothly and seamlessly connect the current schedule, the starting time of the new schedule should be appropriately determined with respect to each furnace. Thus, this paper developed a methodology to handle this issue. The general idea is to maintain the ongoing batch operation of each cracking furnace based on the current schedule, and the new schedule will start right after the cleanups. Therefore, the starting time information of each furnace will be a part of the initial condition used for the rescheduling. Once all the initial conditions are ready, the scheduling model will be solved again to obtain the new schedule. Iteratively, the cracking furnace system performs rescheduling for any new delivery of feeds, and thus scheduling activities are accomplished dynamically and reactively. 12028
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
3. DYNAMIC SCHEDULING MODEL In this section, the dynamic scheduling strategy for an ethylene cracking furnace system is introduced, which includes the MINLP scheduling model and the initial condition identification for rescheduling. For illustration, some equations and auxiliary figures are employed to explain the scheduling model. All the variables and parameters are explained both in context and in the Nomenclature. 3.1. Scheduling Model. 3.1.1. Objective Function. The objective in this dynamic scheduling model is to maximize the average net profit per time during the current schedule, which takes into account the sale incomes of products subtracting the raw material cost, operational cost, and cleanup cost as shown in eq 1. max J ¼
∑ ∑ ∑f∑ i ∈ SI i ∈ SJ k ∈ SK i ∈ SL
Z Pl ð
ti, j, k
Di, j ðci, j, l þ ai, j, l ebi, j, l t Þ
0
dtÞ ðCri þ Cv i, j ÞDi, j ti, j, k Csi, j yi, j, k g Tt ð1Þ where Tt is the total time length of the current schedule. ci,j,l + ai,j,lebi,j,l t characterizes the dynamic change of product l's yield with respect to time, when feed i is cracked in furnace j during one batch processing. The negative or positive sign of bi,j,l determines the increasing or decreasing trends of a product yield with time. Because only one furnace (denoted as the NCth furnace) is used to crack the ethane (usually denoted as the NFth feed), and all the other furnaces can process all the feeds except ethane, some associated particular parameters of ci,j,l and ai,j,l are set to zero for this purpose. Di,j is the feed flow rate when feed i is processed in furnace j. Pl is the unit price of product l. The cleanup cost is controlled by a binary value yi,j,k. yi,j,k is 1 if feed i is processed in the kth batch in furnace j; otherwise, it is 0. Because of the binary control, the batch processing time, ti,j,k, is subject to yi,j,k (if yi,j,k is 0, ti,j,k will be also 0), which will be presented in the model constraints. To calculate the net profit of one batch operation, the material and operational costs, (Cri + Cvi,j)Di,jti,j,k, and the cleanup cost, Csi,jyi,j,k, need to be deducted. Thus, the numerator of eq 1 is the total net profit during the current schedule. It is divided by the total time, suggesting that the objective function is to maximize the average net profit per time (day) of the furnace system during the current schedule. Note that the ethane product from all of the cracking operations will be recycled to the ethane furnace for secondary cracking. To account the total ethane revenue (sale income minus material cost), the actual calculation is represented by eq 2. Since PNP = CrNF, eq 2 can be completely embedded into eq 1 for the total profit calculation. Also note that the amount of consumed ethane should be less than the total ethane inventory. ENF ¼
∑ ∑ ∑ ½PNPð i ∈ SI i ∈ SJ i ∈ SK
Z t i, j, k 0
CrNF DNF, NC tNF, NC, k
Figure 4. Illustrative examples for binary variable of xj,k,j0 ,k0 2. (a) xj,k,j0 ,k0 = 1 (the kth cleanup in furnace j is no-overlap behind the k0 th cleanup in furnace j0 ). (b) xj,k,j0 ,k0 = 0 (the kth cleanup in furnace j is no-overlap ahead of the k0 th cleanup in furnace j0 ).
of eq 1 is formulated as below. λi, j, k, 1 e zi, j, k, 1 , " i ∈ SI;
λi, j, k, n e zi, j, k, n1 þ zi, j, k, n , " k ∈ SK;
ð3Þ
" i ∈ SI; " j ∈ SJ;
"n ∈ SN and 1 < n < NN
λi, j, k, NN e zi, j, k, NN1 ,
" i ∈ SI;
ð4Þ
" j ∈ SJ;
" k ∈ SK
ð5Þ
∑
n ∈ SN
zi, j, k, n ¼ 1,
" i ∈ SI;
" j ∈ SJ;
" k ∈ SK; "n ∈ SN and n < NN
∑
n ∈ SN
λi, j, k, n ¼ 1,
" i ∈ SI;
" j ∈ SJ;
ð6Þ " k ∈ SK ð7Þ
ti, j, k ¼
∑
n ∈ SN
λi, j, k, n PTn ,
" i ∈ SI;
" j ∈ SJ;
" k ∈ SK
ð8Þ PFi, j, n ¼
Di, j ðci, j, NP þ ai, j, NP ebi, j, NP t Þ dtÞ
R ti,j,k
" j ∈ SJ; " k ∈ SK
∑ Pl ð l ∈ SL
Z PT n 0
Di, j ðci, j, l þ ai, j, l ebi, j, l t Þ dtÞ
ðCri þ Cvi, j ÞDi, j PTn
ð2Þ
Overall, it can be seen that ∑i∈SI ∑l∈SL Pl( 0 Di,j(ci,j,l + ai,j,le ) dt) in eq 1 is the total sale income for all the products during one batch operation. The numerator of eq 1 is a highly nonlinear function. It can be approximated by piecewise linearization to facilitate solution identifications according to the application requirement of online reactive scheduling. By selecting n points, n 1 linear sections could be imposed on a nonlinear function. The derived piecewise linear approximation23 of the numerator bi,j,l t
pi, j, k ¼
∑
n ∈ SN
ð9Þ
λi, j, k, n PFi, j, n , " i ∈ SI; " j ∈ SJ; " k ∈ SK
ð10Þ
max J ¼ 12029
∑ ∑ ∑
i ∈ SI i ∈ SJ i ∈ SK
ðpi, j, k Csi, j yi, j, k Þ
ð11Þ
Tt dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
In piecewise linear approximation, a positive variable λi,j,k,n and a binary variable zi,j,k,n are introduced to ensure the linear combination of two successive values. From eqs 36, because of the binary variable zi,j,k,n, and ∑n∈SN zi,j,k,n = 1 in eq 6, at most two successive λi,j,k,n values are strictly positive. Also λi,j,k,n is a variable with a value from 0 to 1 (eq 7). One batch time length is divided into n 1 segments by n points PTn, and the batch processing time ti,j,k is the linear combination of two successive time points PTn
(as shown in eq 8). PFi,j,n is the value of net profit during a batch time PTn before subtracting the cleanup cost (eq 9), and the approximated net profit before subtracting the cleanup cost from a single batch pi,j,k is the linear combination of two successive PFi,j,n values (as shown in eq 10). The objective function can be simplified by eq 11. Although the reduced objective function in eq 11 is still a nonlinear function, it is simplified a lot by removal of the
Table 1. Parameter Values for the Cracking Furnace System furnace j parameter τi,j (day)
Di,j (tons/day)
tloi,j (day)
tupi,j (day)
ai,j,Pa
ai,j,Pb
ai,j,Pd
bi,j,Pa (1/day)
bi,j,Pb (1/day)
bi,j,Pd (1/day)
ci,j,Pa
feed i
1
2
3
4
5
6
Fa
2
2
2
2
2
0
Fb
2
2
3
2
2
0
Fc
3
2
3
3
2
0
Fpd Fa
0 204.4
0 214.9
0 188.2
0 202.3
0 212.7
3 0
Fb
202.2
213.7
187.2
199.6
210.6
0
Fc
201.0
211.1
185.8
197.6
209.8
0
Fpd
0
0
0
0
0
210.1
Fa
20
23
18
20
23
0
Fb
19
21
17
18
22
0
Fc
18
19
16
17
20
0
Fpd Fa
0 70
0 72
0 65
0 69
0 73
59 0
Fb
53
54
49
52
56
0
Fc
41
43
39
43
45
0
Fpd
0
0
0
0
0
101
Fa
0.0539
0.0500
0.0574
0.0510
0.0585
0
Fb
0.0539
0.0560
0.0534
0.0528
0.0584
0
Fc
0.0129
0.0140
0.0139
0.0125
0.0134
0
Fpd Fa
0 0.0267
0 0.0243
0 0.0288
0 0.0242
0 0.0278
0.0542 0
Fb
0.0267
0.0287
0.0269
0.0244
0.0257
0
Fc
0.0149
0.0149
0.0143
0.0146
0.0148
0
Fpd
0
0
0
0
0
0.0264
Fa
0.2445
0.1675
0.3466
0.2318
0.2901
0
Fb
0.0783
0.0855
0.0940
0.0723
0.0671
0
Fc
0.0436
0.0485
0.0525
0.0420
0.0394
0
Fpd Fa
0 0.0050
0 0.0054
0 0.0053
0 0.0049
0 0.0050
0.1659 0.0010
Fb
0.0060
0.0062
0.0054
0.0058
0.0057
0.0010
Fc
0.0070
0.0070
0.0069
0.0073
0.0068
0.0010
Fpd
0.0010
0.0010
0.0010
0.0010
0.0010
0.0051
Fa
0.0041
0.0043
0.0038
0.0037
0.0038
0.0040
Fb
0.0041
0.0042
0.0038
0.0038
0.0043
0.0040
Fc
0.0105
0.0100
0.0115
0.0112
0.0100
0.0040
Fpd Fa
0.0040 0.0007
0.0040 0.0010
0.0040 0.0006
0.0040 0.0007
0.0040 0.0006
0.0039 0.0010
Fb
0.0012
0.0011
0.0010
0.0013
0.0014
0.0010
Fc
0.0014
0.0013
0.0012
0.0015
0.0016
0.0010
Fpd
0.0010
0.0010
0.0010
0.0010
0.0010
0.0010
Fa
0.5900
0.6378
0.5637
0.5987
0.6432
0
Fb
0.3600
0.3698
0.3548
0.3626
0.3628
0
Fc
0.3000
0.3109
0.2979
0.3011
0.2985
0
Fpd
0
0
0
0
0
0.6067
12030
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Table 1. Continued furnace j feed i
parameter ci,j,Pb
ci,j,Pd
1
2
Fc
0.1600
0.1523
0.1631
0.1572
0.1513
0
Fpd
0
0
0
0
0
0.0666
Fa
0.5701
0.4430
0.6944
0.5359
0.5733
0
Fb
0.1136
0.1208
0.1293
0.1076
0.1024
0
Fc
0.0686
0.0735
0.0775
0.0670
0.0644
0
Fpd
0
0
0
0
0
0.4732
" i ∈ SI
" j ∈ SJ
ð13Þ
∑
yi, j, 1 ¼ 1,
" j ∈ SJ
∑
yi, j, k e 1,
" j ∈ SJ;
i ∈ SI
ð14Þ " k ∈ SK, k > 1
0.0599 0.1524
0 0
Table 2. Operational Cost and Cleanup Cost Values furnace j parameter
feed i
1
2
3
4
5
6
Cvi,j ($/ton) Fa
16.55
17.74
15.55
18.18
16.66
0
Fb
16.80
17.61
18.27
17.91
18.25
0
Csi,j ($)
Fc
17.20
15.50
16.06
17.29
16.18
0
Fpd
0
0
0
0
0
16.94
Fa
86,700.0 92,758.3 81,960.6 94,973.8 86,733.3 0
Fb
87,361.0 90,964.4 79,695.6 80,087.9 82,357.1 0
Fc
89,154.0 87,879.6 84,239.9 81,375.3 92,398.7 0
Fpd
0
0
0
0
0
88,625.2
constraints can be modeled as shown in eqs 16 and 17.16
ð12Þ
∑ ðDi, j k ∈∑SK ti, j, kÞ e Di, jtloi, j , j ∈ SJ
0.0724 0.1620
6
0.0747 0.1663
3.1.3. Integrality Constraints. Prior to the solution identification, it is unknown how many batches will be needed for each furnace during a schedule. Thus, the total number of batches, NB, will be given a heuristic upper limit. Some batches may not be utilized for cracking. Some logic constraints borrowed from Liu et al.2 may help reduce the solution searching space as shown in eqs 14 and 15. Equation 14 indicates that for each furnace the first batch slot should always be used for cracking a feed. Equation 15 means that each batch slot could only be permitted for cracking one feed at most. i ∈ SI
5
0.0610 0.1581
∑ ðDi, j k ∈∑SK ti, j, kÞ e Lni þ Lli , j ∈ SJ " i ∈ SI;
4
0.0650 0.1600
exponential part from the equation. The nonlinearity of eq 11 is only from the division. Because both numerator and denominator are linear and differentiable, the objective function is a pseudoconvex or pseudoconcave function. 3.1.2. Mass Balance Constraints. For each cracking feed, the total amount consumed by all cracking furnaces should not be over the inventory limit. The feed inventory except ethane comes from two sources: the new coming feed supply (Lni) and the leftover inventory (Lli) for the current schedule. For the ethane feed, although there is no ethane supply coming from the purchase, its inventory is also from two sources: the recycled ethane (LnNF) produced from cracking operations for the current schedule, and the leftover ethane inventory (LlNF) for the current schedule. Thus, no matter which feed, the consumed amount should be less than the total inventory (Lni + Lli), which is shown in eq 12. Equation 13 suggests that, after a schedule is completely executed, any leftover feeds will not have enough amounts to comprise a minimum batch in any furnace.
Lni þ Lli
3
Fa Fb
yi, NC, k ¼ 0,
" i ∈ SI, i 6¼ NF,
" k ∈ SK
ð16Þ
yNF, j, k ¼ 0,
" j ∈ SJ, j 6¼ NC,
" k ∈ SK
ð17Þ
3.1.5. Timing Constraints. Timing constraints are used to define and restrict the time variables of ti,j,k, Sj,k, and Ej,k. Equation 18 suggests that the batch processing time, ti,j,k, should be within a time window, where tloi,j and tupi,j are the lower and upper bounds of a batch run length. Sj,l, as shown in eq 19, is designated as the starting time of the first batch in furnace j. It is equal to the leftover batch processing time (Tlj), which is determined based on the delivery time of the new feed and the operational conditions of this furnace. Tlj will be calculated in the section of the initial condition identification. Equations 20 and 21 present the relations between starting and ending time points for every batch processing. In eq 22, the total time length of a schedule minus the leftover batch processing time calculated from the previous schedule should be exactly equal to the summation of all the batch processing times plus the cleanup time. Equation 23 gives the upper limit of the starting time of a furnace. Note that the starting and ending times of a batch are arranged in such a way that if the kth batch is not actually utilized, its starting and ending time points will be exactly the same.
ð15Þ
tloi, j yi, j, k e ti, j, k e tupi, j yi, j, k , " j ∈ SJ,
3.1.4. Feed Allocation Constraints. Because only one cracking furnace is specifically used to process ethane, and all the other furnaces can process all the feeds except ethane, these feed allocation
Sj, 1 ¼ Tlj , 12031
" k ∈ SK
" j ∈ SJ
" i ∈ SI, ð18Þ ð19Þ
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Table 3. Additional Parameter Values of Feeds and Products furnace j
feed i
Tlj (day)
1
16
2
2
3
Cri ($/ton)
Lni (tons)
product l
Lli (tons)
Pl ($/ton)
Fa (gas)
444.13
22 833.7
20 899.6
Pa (ethylene)
Fb (naphtha)
415.49
89 456.3
40 722.7
Pb (propylene)
718.53
27
Fc (light diesel)
411.01
16 313.3
1 005.0
Pc (others)
356.67
4
21
Fpd (ethane)
444.13
27 102.9
10 925.2
Pd (ethane)
444.13
5
0
6
24
866.67
Figure 5. Optimal scheduling results for ethylene production: (a) ethylene yield in current schedule; (b) ethylene yield in current schedule and reschedule.
Sj, k ¼ Ej, k1 þ
∑
i ∈ SI
τi, j yi, j, k1 ,
" k ∈ SK, k 6¼ 1 Ej, k ¼ Sj, k þ
∑
i ∈ SI
ti, j, k ,
∑ ∑
" j ∈ SJ;
" j ∈ SJ;
i ∈ SI k ∈ SK
ðti, j, k þ τi, j yi, j, k Þ ¼ Tt Tlj ,
ð22Þ
ð20Þ " k ∈ SK
" j ∈ SJ
ð21Þ
Sj, NB e Tt, 12032
" j ∈ SJ
ð23Þ
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Table 4. Detailed Mass Flow Based on the Optimal Current Schedule furnace 1
2
3
4
5
batch
feed
batch starting
batch processing
number
type
time (day)
time (day)
Fa (tons)
Fb (tons)
Fc (tons)
Fpd (tons)
Pa (tons)
Pb (tons)
Pc (tons)
Pd (tons)
1
Fb
16
49
9907.8
2941.5
1344.8
5249.6
371.9
2
Fb
67
46
9301.2
2767.4
1261.2
4924.7
347.9
3
Fb
115
45.1
9119.2
2715.0
1236.1
4827.4
340.7
4
Fb
162.1
46.9
9483.2
2819.7
1286.3
5022.1
355.1
41.9
8852.5
2607.5
1239.1
4773.2
232.7
36
7693.2
2360.3
1011.0
4037.6
284.3
7707.1
8466.0
2364.1 2494.1
1012.7 1183.9
4045.5 4565.9
284.8 222.1 2693.8
1
Fc
2
2
Fb
45.9
3 4
Fb Fc
83.9 122
36.1 40.1
5
Fa
164.1
44.9
9648.2
5605.5
374.7
974.2
1
Fb
27
35
6552.0
1938.4
924.3
3447.4
241.9
2
Fb
65
26
4867.2
1447.4
684.6
2557.5
177.7
3
Fb
94
36
6739.2
1992.6
951.0
3546.5
249.1
4
Fb
133
36
6739.2
1992.6
951.0
3546.5
249.1
5
Fb
172
36
6739.2
1992.6
951.0
3546.5
249.1
1 2
Fb Fb
21 69
46 46
9181.6 9181.6
2770.3 2770.3
1281.4 1281.4
4786.5 4786.5
343.4 343.4
3
Fb
117
47.1
9401.2
2834.5
1312.4
4902.3
352.0
4
Fb
166.1
42.9
8562.8
2588.9
1194.0
4460.8
319.1
1
Fa
0
3430.3
2
Fb
58
42.8
3
Fa
102.8
48.2
56
4
Fa
153
56
6
1 2
Fpd Fpd
24 120
93 88
total
6851.9
414.1
1214.9
9003.2
2668.0
1160.1
4839.7
335.4
10 262.8
5918.7
353.3
1041.9
2948.9
11 911.2
11 911.2
43 733.4
130 178.9
3.1.6. Nonsimultaneous Cleanup Constraints. Cleanup time intervals among different furnaces should have no overlap. Through the study of the nonoverlap scenarios, as shown in Figure 4, eqs 24 and 25 can be utilized to linearly represent such nonoverlap situations. Here, another binary variable, xj,k,j0 ,k0 , is introduced,2 which is designated as 1 if the kth cleanup in furnace j is no-overlap behind the k0 th cleanup in furnace j0 ; otherwise, if the kth cleanup in furnace j is no-overlap ahead of the k0 th cleanup in furnace j0 , it is 0. As shown in eqs 24 and 25, xj,k,j0 ,k0 equals 1 and 0 resulting in the nonoverlap scenarios of parts a and b, respectively, of Figure 4. Note that this constraint cannot be applied to the last cleanup of each furnace because all the furnaces have to stop running thereafter. However, this situation is seldom attained as the new coming feeds will supposedly trigger a new schedule before the plant inventory is depleted. ðxj, k, j0 , k0 1ÞM e Ej, k Sj0 , k0 þ1 e xj, k, j0 , k0 M, " j, j0 ∈ SJ, j < j0 ;
" k, k0 ∈ SK, k, k0 ¼ 6 NB
ð24Þ
xj, k, j0 , k0 M e Ej0 , k0 Sj, kþ1 e ð1 xj, k, j0 , k0 ÞM, " j, j0 ∈ SJ, j < j0 ;
" k, k0 ∈ SK, k, k0 ¼ 6 NB
ð25Þ
17 318.5
6851.9
414.1
1214.9
3430.3
19 539.3 18 488.8
10 474.1 9931.3
863.9 814.1
2052.9 1932.6
6148.4 5810.8
38 028.1
89 698.6
23 500.5
86 297.6
29 762.2
cleanup in furnace j0 (k00 > k), and vice versa.
∑
i ∈ SI
yi, j, k g
∑
i ∈ SI
yi, j, k0 , " k, k0 ∈ SK and k < k0 , " j ∈ SJ
ð26Þ xj, k, j0 , k0 e xj, k00 , j0 , k0 , " j, j0 ∈ SJ, j < j0 ; " k, k0 , k00 ∈ SK, k < k00
ð27Þ
xj, k, j0 , k0 g xj, k, j0 , k00 , " j, j0 ∈ SJ, j < j0 ; " k, k0 , k00 ∈ SK, k < k00
ð28Þ
Equation 29 represents the logic relations characterized by binary variables of x and y: if two consecutive x values change from 1 to 0, i.e., xj,k,j0 ,k0 xj,k,j0 ,k0 +1 = 1, then there must be a batch slot being used between these two cleanups; i.e., ∑i∈SI yi,j0 ,k0 +1 = 1. The algebraic constraint of eq 29 is shown in eq 30. xj, k, j0 , k0 xj, k, j0 , k0 þ1 ¼ 1 w
3.1.7. Additional Logic Constraints. Three types of logic constraints can be employed to further reduce the solution space. Equations 26 through 28 are borrowed from a previous study.16 Specifically, eq 26 shows that if one batch slot of a furnace is not utilized, the followed batches of that furnace should not be utilized either. Equations 27 and 28 show that if the kth cleanup in furnace j is after the k0 th cleanup in furnace j0 , then the k00 th cleanup in furnace j should also be after the k0 th
6 j0 ; " j, j0 ∈ SJ, j ¼
∑
i ∈ SI
yi, j0 , k0 þ1 ¼ 1,
" k, k0 ∈ SK, k0 < NB
ð29Þ
yi, j0 , k0 þ1 g xj, k, j0 , k0 xj, k, j0 , k0 þ1 , 6 j0 ; " j, j0 ∈ SJ, j ¼
12033
∑
i ∈ SI
" k, k0 ∈ SK, k0 < NB
ð30Þ
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Table 5. Detailed Cash Flow Based on the Optimal Current Schedule batch furnace number 1
2
3
4
5
feed type
batch
batch
ethylene
propylene
ethane
other
feed
starting time (day)
processing time (day)
revenue (M$)
revenue (M$)
revenue (M$)
product revenue (M$)
cost (M$)
operation cost (M$)
cleanup
batch net
cost (M$)
profit (M$)
1
Fb
16
49
2.55
0.97
0.17
1.87
4.12
0.17
0.09
1.18
2
Fb
67
46
2.40
0.91
0.15
1.76
3.86
0.16
0.09
1.10
3
Fb
115
45.1
2.35
0.89
0.15
1.72
3.79
0.15
0.09
1.08
4
Fb
162.1
46.9
2.44
0.92
0.16
1.79
3.94
0.16
0.09
1.13
1
Fc
2
41.9
2.26
0.89
0.10
1.70
3.64
0.14
0.09
1.09
2 3
Fb Fb
45.9 83.9
36 36.1
2.05 2.05
0.73 0.73
0.13 0.13
1.44 1.44
3.20 3.20
0.14 0.14
0.09 0.09
0.92 0.92
4
Fc
122
40.1
2.16
0.85
0.10
1.63
3.48
0.13
0.09
1.04
5
Fa
164.1
44.9
4.86
0.27
1.20
0.35
4.29
1 2
Fb Fb
27 65
35 26
1.68 1.25
0.66 0.49
0.11 0.08
1.23 0.91
2.72 2.02
0.17 0.12
0.09 0.08
2.13 0.76
3
Fb
94
36
1.73
0.68
0.11
1.26
2.80
0.09
0.08
0.55
0.12
0.08
0.78
4
Fb
133
36
1.73
0.68
0.11
1.26
2.80
0.12
0.08
0.78
5 1
Fb Fb
172 21
36 46
1.73 2.40
0.68 0.92
0.11 0.15
1.26 1.71
2.80 3.81
0.12 0.16
0.08 0.08
0.78 1.12
2
Fb
69
3
Fb
117
4
Fb
166.1
1
Fa
0
46
2.40
0.92
0.15
1.71
3.81
0.16
0.08
1.12
47.1
2.46
0.94
0.16
1.75
3.91
0.17
0.08
1.15
42.9
2.24
0.86
0.14
1.59
3.56
0.15
0.08
1.04
56
5.94
0.30
1.52
0.43
5.29
0.20
0.09
2.61
2
Fb
58
42.8
2.31
0.83
0.15
1.73
3.74
0.16
0.08
1.04
3
Fa
102.8
48.2
5.13
0.25
1.31
0.37
4.56
0.17
0.09
2.25
6
4 1
Fa Fpd
153 24
56 93
5.94 9.08
0.30 0.62
1.52 2.73
0.43 0.73
5.29 8.68
0.20 0.33
0.09 0.09
2.61 4.06
2
Fpd
120
88
8.61
0.58
2.58
0.69
8.21
0.31
0.09
3.85
total
77.74
16.89
13.22
30.78
97.52
3.95
2.06
35.09
3.1.8. Variable Bounds. All the continuous variables have a lower bound of 0, and all the starting times, ending times, batch processing times, and the total time scale should be less than the upper bound. x and y are defined as binary variables. Ej, k , Sj, k , ti, j, k , Tt g 0, " j ∈ SJ, " k ∈ SK
tR 150
" i ∈ SI,
Ej, k , Sj, k , ti, j, k , Tt e M, " j ∈ SJ, " k ∈ SK
ð31Þ " i ∈ SI,
yi, j, k , xj, k, j0 , k0 , zi, j, k, n ∈ f0, 1g, " j, j0 ∈ SJ,
Table 6. Initial Conditions for the Optimal Rescheduling
" k, k0 ∈ SK,
ð32Þ
furnace j TlRj (days) Kj 1
12.1
3
Fa (gas)
40 781.1
21 560.2
2 3
14.1 22
4 4
Fb (naphtha) Fc (light diesel)
44 919.5 13 247.6
24 822.8 0.0
4
16.1
3
Fpd (ethane)
22 718.4
0.0
5
3
3
6
61
2
TlRj ¼ Ej, Kj t R þ
" i ∈ SI, " n ∈ SN and n < NN
∑
τi, j yi, j, Kj ,
∑
ðDi, j
i ∈ SI
LlRi ¼ Lni þ Lli
ð33Þ
LnRi (tons) LlRi (tons)
feed i
j ∈ SJ
∑
k ∈ Kj
" j ∈ SJ ti, j, k Þ,
Other than a cyclic scheduling, the furnace system can be dynamically rescheduled at any time with respect to a feed supply change. When the delivery scenario of new feeds is finalized, the delivery time tR, the new feed types, and quantities are all fixed, and the rescheduling for the cracking furnace system will be triggered. To seamlessly connect the new schedule to the previous one, the initial conditions for the rescheduling, which need to be identified, are shown in eqs 3438. k
" j ∈ SJ
ð34Þ
ð35Þ
" i ∈ SI ð36Þ
3.2. Initial Condition Identification for Rescheduling.
Kj ¼ fkj min ðt R sj, k Þ g 0g,
LlRNF ¼ LnNF þ LlNF DNF, NC
LnRNF ¼
∑ ∑ ∑ð i ∈ SI j ∈ SJ k ∈ K j
Z ti, j, k 0
∑
k ∈ Kj
tNF, NC, k
ð37Þ
Di, j ðci, j, NP þ ai, j, NP ebi, j, NP t Þ dtÞ
ð38Þ Here, eq 34 is used to identify the index of the processing batch where the new feeds are timely allocated based on the previous 12034
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Table 7. Detailed Mass Flow Based on the Optimal Reschedule furnace 1
2
3
4
5
6 total
batch
feed
batch starting
batch processing
number
type
time (day)
time (day)
Fa (tons)
Fb (tons)
Fc (tons)
Fpd (tons)
Pa (tons)
Pb (tons)
Pc (tons)
Pd (tons)
1
Fb
12.1
34.3
6929.8
2078.5
935.8
3660.0
255.5
2
Fb
48.4
36
7279.2
2180.9
983.7
3845.6
269.0
3
Fb
86.4
36
7279.2
2180.9
983.7
3845.6
269.0
4
Fb
124.4
36
7279.2
2180.9
983.7
3845.6
269.0
1
Fa
14.1
34.3
7365.1
4291.5
282.5
741.4
2049.7
2
Fa
50.4
36
7736.4
4506.0
297.4
778.7
2154.3
3 4
Fa Fa
88.4 126.4
36 34
7736.4 7306.6
4506.0 4258.0
297.4 280.2
778.7 735.0
2154.3 2033.4
1
Fb
22
30.4
5685.6
1686.5
800.9
2989.5
208.7
2
Fc
55.4
35.3
6558.7
1850.1
991.6
3545.9
171.1
3
Fc
93.7
4
Fb
132.7
36
6688.8
1886.5
1011.5
3616.2
174.6
26.7
4998.2
1485.8
703.2
2626.6
182.6
1
Fb
16.1
34.3
6840.7
2079.4
951.4
3557.7
252.2
2
Fb
52.4
36
7185.6
2182.1
999.9
3738.1
265.5
3 4
Fb Fb
90.4 128.4
36 32
7185.6 6387.2
2182.1 1944.4
999.9 887.8
3738.1 3320.1
265.5 234.9
1
Fa
3
36
7657.2
4433.0
259.4
772.5
2192.3
2
Fa
41
36
7657.2
4433.0
259.4
772.5
2192.3
3
Fa
79
4 1
Fa Fpd
117 61
36 43.4 98.4
7657.2
67 050.3
9225.2 62 341.3
13 247.5
4433.0
259.4
772.5
2192.3
20 668.0 20 668.0
5328.4 11 053.9 71 160.9
315.6 917.8 14 402.2
934.4 2184.6 50 799.3
2646.8 6511.7 26 944.7
feed
Table 8. Detailed Cash Flow Based on the Optimal Reschedule
furnace
batch
batch
ethylene
propylene
ethane
other
batch
feed
starting
processing
revenue
revenue
revenue
product
cost
operation
cleanup
batch net
number
type
time (day)
time (day)
(M$)
(M$)
(M$)
revenue (M$)
(M$)
cost (M$)
cost (M$)
profit (M$)
1
Fb
12.1
34.3
1.80
0.67
0.11
1.31
2.88
0.12
0.09
0.81
2
Fb
48.4
36
1.89
0.71
0.12
1.37
3.02
0.12
0.09
0.85
3
Fb
86.4
36
1.89
0.71
0.12
1.37
3.02
0.12
0.09
0.85
4
Fb
124.4
36
1.89
0.71
0.12
1.37
3.02
0.12
0.09
0.85
1
Fa
14.1
34.3
3.72
0.20
0.91
0.26
3.27
0.13
0.09
1.61
2
Fa
50.4
36
3.91
0.21
0.96
0.28
3.44
0.14
0.09
1.69
3
Fa
88.4
36
3.91
0.21
0.96
0.28
3.44
0.14
0.09
1.69
4 1
Fa Fb
126.4 22
34 30.4
3.69 1.46
0.20 0.58
0.90 0.09
0.26 1.07
3.25 2.36
0.13 0.10
2 3
Fc Fc
55.4 93.7
35.3 36
1.60 1.63
0.71 0.73
0.08 0.08
1.26 1.29
2.70 2.75
0.11 0.11
0.09 0.08 0.08
1.59 0.65 0.78
0.08
0.79
4
Fb
132.7
26.7
1.29
0.51
0.08
0.94
2.08
0.09
0.08
0.56
1
Fb
16.1
34.3
1.80
0.68
0.11
1.27
2.84
0.12
0.08
0.82
2
Fb
52.4
36
1.89
0.72
0.12
1.33
2.99
0.13
0.08
0.87
3
Fb
90.4
36
1.89
0.72
0.12
1.33
2.99
0.13
0.08
0.87
4 1
Fb Fa
128.4 3
32 36
1.69 3.84
0.64 0.19
0.10 0.97
1.18 0.28
2.65 3.40
0.11 0.13
0.08 0.09
0.76 1.66
2
Fa
41
36
3.84
0.19
0.97
0.28
3.40
0.13
0.09
1.66
3
Fa
79
36
3.84
0.19
0.97
0.28
3.40
0.13
0.09
1.66
4
Fa
117
43.4
4.62
0.23
1.18
0.33
4.10
0.15
0.09
2.01
6
1
Fpd
61
98.4
9.58
0.66
2.89
0.78
9.18
0.35
0.09
4.29
total
61.67
10.35
11.97
18.12
70.17
2.81
1.81
27.32
1
2
3
4
5
schedule. The index of Kj may be different with respect to the furnace index of j. Based on the identified Kj, the leftover batch
processing time for each furnace can be calculated by eq 35. It represents the starting time for the furnace j to conduct the 12035
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Figure 6. Heuristic scheduling results for ethylene production: (a) ethylene yield in current schedule; (b) ethylene yield in current schedule and reschedule.
rescheduling activity. As aforementioned, the reason for this arrangement is that it is better to postpone the rescheduling activity until a furnace completes its current batch processing and the following cleanup (see Figure 1). Equation 36 calculates the leftover inventory for the rescheduling. Since ethane feed comes from the recycled ethane, its calculation is given by eqs 37 and 38. Specifically, eq 37 calculates the amount of leftover ethane; eq 38 gives the generated ethane from the completed cracking operation based on the current schedule, which is designated as LnNF for the rescheduling. In summary, the static scheduling model including eqs 328 and 3033 and the rescheduling initialization described by eqs 3438 constructs a completed reactive scheduling MINLP model. Since the nonlinearity of each static scheduling model comes from the pseudoconvex or pseudoconcave objective function, each reactive scheduling activity actually handles a pseudoconvex or pseudoconcave optimization problem.
4. CASE STUDY To demonstrate the efficacy of the developed dynamic scheduling strategy, a case study with industrial data has been conducted. In this
case study, an ethylene plant processes four types of feeds (NF = 4): liquefied petroleum gas (LPG), naphtha, light diesel, and ethane (denoted as Fa, Fb, Fc, and Fpd, respectively). Among them, ethane feed comes from the recycled self-cracking product. Thus, it generally works in the way of self-producing and self-consuming. The feeds of LPG, naphtha, and light diesel will be purchased intermittently, but their exact delivery date and quantity are considered unpredictable. Four types of products (NP = 4) are considered: ethylene, propylene, ethane, and other products (denoted as Pa, Pb, Pd, and Pc, respectively). The plant has six cracking furnaces in total (NC = 6). Among them, five are ordinary furnaces, which can process LPG, naphtha, and light diesel; one furnace is an ethane furnace, processing the recycled ethane specifically. The product yields for each cracking furnace are nonlinearly changing with respect to time. For instance, the ethylene yield is decreasing with batch processing time, while the propylene yield is increasing with batch processing time. Generally, the product yield change can be approximated by an exponential function as shown in eq 1. For this case study, the yield model parameters for the major products such as ethylene, propylene, and ethane are obtained through the regression of industrial historical data, 12036
dx.doi.org/10.1021/ie200318p |Ind. Eng. Chem. Res. 2011, 50, 12026–12040
Industrial & Engineering Chemistry Research
ARTICLE
Figure 7. Optimal scheduling results for propylene production: (a) propylene yield in current schedule; (b) propylene yield in current schedule and reschedule.
which is shown in Table 1. The yield of the left products is equal to 1 minus the summation of the above three product yields. Obviously, the product yield varies with feeds and furnaces. Table 1 also gives the cleanup time and batch processing time limits for each furnace. The data of raw material cost, operational cost, cleanup cost, and the product prices are listed in Tables 2 and Table 3 based on the industrial practice. In the scheduling model, the upper bound M for the total schedule time is 365 days. Based on the given data and some initial inventory data, the developed MINLP dynamic scheduling model is solved by GAMS with the solver of DICOPT.24 MILP and NLP subproblems are solved by CPLEX and CONOPT, respectively.25,26 As shown in the elaborated methodology framework, the MINLP model will be reactively resolved whenever a new delivery of feeds has been determined. For comparison, a heuristic scheduling method based on the industrial expertise is also employed to give the schedules for the furnace system. The comparison is based on the same initial conditions specified in Table 3. Both methods will give a current schedule and a reschedule. The developed MINLP model using piecewise linearization method for the current schedule has 1575 binary variables, 3081
continuous variables, and 3883 constraints. The solving time for the optimal scheduling of the piecewise linearized model is
