- Email: [email protected]
Supply chain optimization of large-scale continuous processes
Computers and Chemical Engineering 29 (2005) 1305–1316
Supply chain optimization of large-scale continuous processes E.P. Schulz, M.S. Diaz ∗ , J.A. Bandoni Planta Piloto de Ingenier´ıa Qu´ımica—PLAPIQUI (Universidad Nacional del Sur—CONICET), Camino La Carrindanga Km 7, 8000 Bah´ıa Blanca, Argentina Available online 11 April 2005
Abstract In this work, a supply chain model for a petrochemical complex is formulated as a multiperiod mixed integer nonlinear program (MINLP). The model includes production, product delivery, inventory management and decisions such as individual production levels for each product, as well as operating conditions for each plant in the complex. Two multiperiod MINLP are formulated in order to compare models with different levels of rigorousness. The first one includes linear mathematical models derived for the ethylene plants while the second considers semi-rigorous nonlinear models. Simplified models based on rigorous existing models tuned with actual plant data are considered for the natural gas plants, taking into account variations in production with key plant operating variables. Data for chemical transformations and utilities consumption from the literature have been used to model the remaining plants. © 2005 Published by Elsevier Ltd. Keywords: Multiperiod MINLP; Short term planning; Petrochemical complex
1. Introduction Petrochemical industries are within wider supply chains and involve a number of players, which span sourcing, manufacturing and distribution. Every day, much money is spent on feedstocks, chemicals, utilities, equipment, labor and effluent treatment in order to produce marketable process-industry products on time. In a petrochemical complex there may be several plants operating under different conditions producing several categories of petrochemicals (basic, intermediate and end-product). The responses to forecast market demands have to be coordinated and simultaneous planning of production and each plant production distribution has to be undertaken. Decisions have to be taken at the different stages of the chain, i.e. operating conditions of individual processes, intermediate product handling and end-product storage and distribution. A model that can account for the different interactions among the plants to provide a suitable production and operation plan for the multi-plant facility is a valuable tool for the decision-making process. Turkay, Asakura, Fujita, Hui, and Natori (1998) applied logic-based approaches proposed by Turkay and Grossmann ∗
Corresponding author.
0098-1354/$ – see front matter © 2005 Published by Elsevier Ltd. doi:10.1016/j.compchemeng.2005.02.025
(1996) to solve a mixed integer linear problem for the total site optimization of a petrochemical complex. Furlonge and Young-Hoon (2000) presented a methodology to address long-term planning for development of the industry in Trinidad and Tobago, selecting an initial structure and policy for the development of an olefin-based complex. Al-Sharrah, Alatiqi, Elkamel, and Alper (2001) proposed an environmental objective in the production planning in the petrochemical industry. Bok, Grossman, and Parks (2000) proposed a bilevel decomposition strategy for the solution of a multiperiod MILP model for the supply chain of continuous processes taking into account sales, intermittent deliveries, production shortfalls, delivery delays, inventory profiles and changeovers. As stated by Shah (2004), supply chain planning and scheduling is a tool to optimize production, distribution and storage resources of a fixed network to respond properly to the external conditions, such as orders and demands forecasts, over a short to medium term horizon. Production process representation depends on the gross margin of the business. Business with low margin facilities, such as refining and petrochemicals, require “property based” representations, which use process conditions and process models. An example of property based planning is provided by Jackson and
1306
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Nomenclature (a) Sets h furnace i plant (natural gas plants: NGI and NGII, ethylene plants: EPI and EPII) in inlet stream to a tank, mixer or splitter j component k intermediate tank m mixer out outlet stream from a tank, mixer or splitter poly polyethylenes (HDPE: high density polyethylene, LDPE: low density polyethylene, LLDPE: linear low density polyethylene) s splitter t, n time intervals (from 1 to H, number of time periods) v {vC3 , vC4 , vC5+ } tanks and ships for propane (vC3 ), butane (vC4 ), gasoline (vC5+ ), ammonia (vNH3 ) or urea (vurea) delivery (b) Parameters ai,j (K−1 ), bi,j (bar−1 ) and ci,j correlation coefficients for natural gas plants recoveries inventory cost (US$/kmol) Civv 0 Cpoly initial inventory of polyethylene poly (t) demvv ship v capacity t,U fit,L i,j –fii,j minimum–maximum individual molar flow of component j in plant i in period t (kmol/d) t,U Fvvt,L –Fvv out,v out,v minimum–maximum transfer rate of propane, butane or gasoline from tank v to ship v in period t (kmol/d) IOpoly storage target of polyethylene poly (t) Su,i separation factor in unit u in the separation train in plant i max minimum–maximum loading duration (d) tvmin –tv v v period when ship v starts loading its product TAv (d) v0k,j initial moles of component j in tank k (kmol) VkL –VkU minimum–maximum levels in intermediate tank k (kmol) Vv0 initial moles in tank k (only one component) (kmol) xii,j inlet molar fraction of component j in plant i Greek letters αi,j , βi,j correlation coefficients for ethylene plants simplified models ξj penalty coefficient for product p unmet demand t,U ηt,L –η i,j i,j minimum–maximum component j recovery in demethanizer bottom stream in plant i in period t
θ poly
penalty coefficient for polyethylene poly unmet inventory target (US$/t)
(c) Binary variables yatv denotes when ship v starts loading its product ydtv denotes when ship v finishes loading its product (d) Continuous variables Cdem penalty for not meeting forecast demand Cinv total inventory cost (US$) CIO penalty on deviations from inventory targets for polyethylenes Convth,i ethane conversion in furnace h in period t t inventory of polyethylene poly in period t (t) Cpoly deltav shortfall in ship v capacity deltatpoly demand shortfall for polyethylene poly in period t (t) demtpoly forecast demand for polyethylene poly in period t (t) ffoti,h j production of component j in furnace h in plant i in period t (kmol/d) fiti,j individual molar inlet molar flow of component j in plant i in period t (kmol/d) fmxtin,m,j –fmxtout,m,j individual molar inlet–outlet molar flow of component j in mixer m in period t (kmol/d) individual molar outlet molar flow of compofoti,j nent j in plant i in period t (kmol/d) Fpotpoly sales of polyethylene poly in period t (t/d) fsptin,s,j –fsptout,s,j individual molar inlet–outlet molar flow of component j in splitter s in period t (kmol/d) ftbti,u,j molar flowrate of component j as bottom product in unit u in the separation train in plant i (kmol/d) fttti,u,j molar flowrate of component j as top product in unit u in the separation train in plant i (kmol/d) fvtin,k,j –fvtout,k,j individual molar inlet–outlet molar flow of component j in intermediate tank k in period t (kmol/d) Ffiti total furnaces charge in plant i in period t (kmol/d) Fiti total molar inlet molar flow in plant i in period t (kmol/d) Fmxtin,m –Fmxtout,m total molar inlet–outlet molar flow in mixer m in period t (kmol/d) Fpitpoly production of polyethylene poly in period t (t/d) Frti ethane recycle stream in plant i in period t (kmol/d)
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Frgti
residual gas total flow in plant i in period t (kmol/d) Fsptin,s –Fsptout,s total molar inlet–outlet molar flow in splitter s in period t (kmol/d) Fvtin,k –Fvtout,k total molar inlet–outlet molar flow in intermediate tank k in period t (kmol/d) Fvttout,v outlet molar flow for train or truck continuous delivery of product v (only one component in the stream) in tank v in period t (kmol/d) Fvvtin,v –Fvvtout,v inlet–outlet molar flow for ship discontinuous delivery of product v (only one component in the stream; v = vNH3 , vurea) in tank v in period t (kmol/d) Pit demethanizer top pressure in plant i in time t (bar) t PCGCi cracked gas compression suction pressure in plant i in period t PENPOLYtpoly unmet storage target for polyethylene poly in period t (t) Pfith,i furnace h in plant i inlet pressure in period t (bar) Pfoth,i furnace h in plant i outlet pressure in period t (bar) Rdti hydrocarbon dilution ratio with steam in plant i in period t Rlti ethylene/ethane ratio at the entrance of the separation train in plant i in period t Tit high pressure separation tank in plant i in time t (K) TDv period when ship v finishes loading its product vtk,j moles of component j in intermediate tank k in period t (kmol) Vkt total moles in intermediate tank k in period t (kmol) xspts,j molar fraction of component j in splitter s in period t xvtk,j molar fraction of component j in intermediate tank k in period t xwtv 0–1 continuous variable to denote if ship v is loading its product at time t Greek letter ηti,j component j recovery in demethanizer bottom stream in plant i in period t
1307
benefits of solving for the complex together than for the individual refineries separately. This work presents an integrated multiperiod optimization model for the supply chain of an existing large-scale petrochemical complex with short grid spacing. The mathematical model accounts for supply, production and product multimodal delivery. The petrochemical complex under study comprises two natural gas liquid processing (NGL) plants, two ethylene plants, a chlorine plant, a VCM plant, a PVC plant, two polyethylene plants, an ammonia and an urea plant. Mathematical models have been derived for the NGL and ethylene plants, based on rigorous existing models tuned with actual plant data (Diaz & Bandoni, 1996; Diaz, Serrani, de Beistegui, & Brignole, 1995). Simplified models take into account variations in production with key plant operating variables, such as temperature and pressure in separation units or conversion in reactors. An alternative nonlinear semi-rigorous model for ethylene plants was formulated in order to show the benefits in solving more rigorous models. Available yield data for chemical transformations have been used to model the rest of the petrochemical complex. Multiperiod modeling allows economic parameters, such as demands, to vary with time. The objective function is the maximization of total profit for the entire site, considering economic issues and drivers such as contractual obligations, operating costs and inventory costs subject to constraints on mass balances, bounds on product demands, equipment capacities and conditions, yields and intermediate and final product storage tanks limitations. The most profitable set of product-stock to produce and sell is determined. Final products are intermittently distributed by ships and continuously by trains and trucks while satisfying storage tanks capacities. Discontinuous product delivery is represented by the use of binary variables. Bilinear nonconvex constraints arise from products between inlet flows and product calculated recovery in natural gas plants and from blending and storage of multicomponent streams. The resulting nonconvex large scale mixed integer nonlinear models have been solved with GAMS using an initial point provided by the solution of a linear model (MILP) obtained applying a reformulation–linearisation technique which provides a valid relaxation of the original feasible region. Numerical results provide a realistic insight on supply chain management in the petrochemical complex.
2. Petrochemical complex description Grossmann (2003) who proposed Lagrangean decomposition (spatial and temporal) techniques for the solution of a non-linear programming (NLP) problem that models multisite production planning of production, transportation and sales in a chemical company. Neiro and Pinto (2003) describe a planning model of a refinery complex as a mixed integer nonlinear programming (MINLP) problem, including process unit models and showing significant cost
The petrochemical complex subject of this study has an annual production of 4,800,000 t. Fig. 1 shows its model representation. There are two natural gas processing plants, whose main objective is to extract ethane from natural gas to use it as raw material in the ethylene plants. natural gas plant I, next to the complex, is fed with 24 MM m3 /d of natural gas. Residual gas (mainly methane) is recompressed to pipeline pressure and part of it is taken as feed for the ammonia plant
1308
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Fig. 1. Petrochemical complex model representation.
and the rest is delivered as sales gas. Pure ethane, propane, butane and gasoline are plant products. Natural gas plant II has its cryogenic sector (referred to as demethanizing plant) several kilometers away from the conventional separation train (NGL fractionation plant). The demethanizing plant is fed with 36 MM m3 /d of natural gas. Light gases are separated from the heavy ones (ethane, propane, butanes and gasoline) and injected to the natural gas pipeline. The rich liquefied mixture (5 MM m3 /d of heavy gases) is stored in thermal vessels and pumped to the petrochemical complex where it is fed to containers to equalize the charge, i.e. to damp any pulsation or flow changes that may occur anywhere along the pipeline. The mixture undergoes a conventional distillation train in the NGL fractionation plant to obtain LPG (liquefied petroleum gas: propane, butane and gasoline) and ethane. Ethylene plants process 2300 t/y of pure ethane; ethylene is provided as raw material to polyethylene and VCM plants and the rest is exported. Ethylene plant II is more efficient than ethylene plant I. The ammonia plant production is 2050 t/d, most of which is fed to the urea plant to produce 3250 t/d of urea. Urea, ammonia, polyethylenes and PVC are delivered by ship, train and trucks.
3. Mathematical model The optimization of the supply chain of a petrochemical complex is formulated as a multiperiod MINLP problem resulting from a uniform time discretisation with 1-day time intervals. Binary variables are associated to intermittent product delivery. Two alternative MINLP models have been formulated, with different degree of rigorousness. The first one is formulated with linear models representing the plants, where nonlinearities arise from bilinear equations standing for mixers and splitters. The second one adds nonlinear process models for some plants and operating conditions. The objective is to compare their performance and decide their application depending on the amount of information required. An approximate mixed integer linear model (MILP) has also been formulated applying linearisation techniques to bilinear equations (Quesada & Grossmann, 1995) to obtain a valid initial point required to solve the MINLPs. The objective function is the maximization of total profit, defined as the difference between sales revenue and the total operating cost plus penalties for not meeting demands and inventory targets during the given time horizon. A brief
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316 Table 1 Feed gas composition for natural gas plants Component
Natural gas plant I
Natural gas plant II
CO2 N2 CH4 C2 H6 C3 H8 C4 H10 C4 H10 (iso) C5 H12 C5 H12 (iso) C6 H14
0.0065 0.0144 0.9043 0.0461 0.0176 0.0033 0.0044 0.0015 0.0009 0.0010
0.0212 0.0082 0.9032 0.0381 0.0161 0.0030 0.0049 0.0013 0.0015 0.0025
Table 2 Polyethylenes’ yields Plant
t of product/kmol of C2 H4
LDPE LLDPE HDPE
0.0272 0.0275 0.0273
description of plants and their mathematical models is given below. 3.1. Plants description and models 3.1.1. Natural gas separation plants In both natural gas processing plants, pipeline gas composition has been considered as constant (Tables 1 and 2). The total molar flowrate is a variable and the individual flows have been calculated as (NGI and NGII stand for natural gas plants I and II, respectively): fiti,j = Fiti xii,j
i = NGI, NGII; ∀j; t = 1, . . . , H
(1)
1309
Natural gas plants are fed from several pipelines. The main difference in the composition of the natural gas plants is the carbon dioxide content, as it is shown in Table 1, requiring basic (NGI) and gas subcooled (NGII) process turboexpansion technologies. The cryogenic sector is the core of these plants. In the basic turboexpansion process, inlet gas is cooled by heat exchange with residual gas and demethanizer side and bottom reboilers. The partially condensed gas feed is then sent to a high-pressure separator. The vapor is expanded through a turboexpander to obtain the low temperatures required for high ethane recovery and is then fed to the top of a demethanizer column. The liquid from the high-pressure separator enters the demethanizer at its lowest feed point. Methane and lighter components, the top product, constitute the residual gas, which are re-injected to the pipeline. Carbon dioxide and ethane distribute between top and bottom streams. Heavier hydrocarbons are obtained as bottom product and are further fractionated to obtain pure ethane, propane, butane and gasoline. Fig. 2 shows the cryogenic sector of a basic turboexpansion plant. In a gas subcooled turboexpansion process, a fraction of the vapor from the cold tank is condensed and subcooled by heat exchange with residual gas coming from the demethanizer. The remaining vapor is expanded through the turboexpander and fed to the middle of the column. The subcooled liquid is flashed and fed to the top of the demethanizer column. Linear models have been derived for ethane (ηti,C2 H6 ) and carbon dioxide (ηti,CO2 ) recovery in demethanizer bottom stream, as function of main operating variables in the cryogenic sector: high pressure separation tank temperature and demethanizing column top pressure. These correlations have been obtained based on simulations with a plant rigorous model (Diaz et al., 1995) for conditions not far from actual operating points, with correlation coefficients ranging from
Fig. 2. Cryogenic sector in a basic turboexpansion process.
1310
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
fsptout,s,j = Fsptout,s xspts,j
0.96 to 0.99 ηti,j
=
ai,j Tit
+ bi,j Pit
+ ci,j ,
ai,j , bi,j < 0;
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
(13) (2)
Product flows (foti,j ) are calculated as products between inlet flowrate (fiti,j ) and recovery (ηti,j ), resulting in bilinear
equations
foti,j = fiti,j ηti,j ,
∀j, s, out ∈ Sout ; t = 1, . . . , H
i = NGI, NGII;
j = CO2 , C2 H6 ; t = 1, . . . , H
(3)
In the MILP model, these equations have been replaced by (Quesada & Grossmann, 1995): t,L t t,L t,L t foti,j ≥ fit,L i,j ηi,j + ηi,j fii,j − fii,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
(4)
(5)
(6)
t,U t t,L t,U t foti,j ≤ fit,L i,j ηi,j + ηi,j fii,j − fii,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
(14)
where Sout is the set of outlet streams and component j represents carbon dioxide, ethane, propane, butanes, pentanes and hexane. In the linear model, Eqs. (13) and (14) have been replaced by valid linearisations analogous to Eqs. (4)–(7) (Quesada & Grossmann, 1995). Mass balances in mixers with multicomponent streams are modeled as follows: fmxtin,m,j = fmxtout,m,j ∀j, m; t = 1, . . . , H (15)
Fmxtin,m = Fmxtout,m
∀j, m; t = 1, . . . , H
(16)
in ∈ Si
t,L t t,U t,L t foti,j ≤ fit,U i,j ηi,j + ηi,j fii,j − fii,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
∀j, s; t = 1, . . . , H
in ∈ Si
t,U t t,U t,U t foti,j ≥ fit,U i,j ηi,j + ηi,j fii,j − foi,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
fsptin,s,j = Fsptin,s xspts,j
(7)
where Sin is the set of inlet streams to mixer m and j represent the same components as for the splitter equations. There are multicomponent storage tanks in natural gas plant II and charge equalizers at the entrance of the fractionation plant. The material balance equation for each intermediate container at each discrete interval is described as the initial volume plus the summation of inflows minus the summation of outflows up to each time interval (Lee, Pinto, Grossman, & Park, 1996): t
t
For propane, butanes, pentanes and hexane, 100% recovery is assumed. Residual gas flowrate (mainly methane and nitrogen) is calculated as the difference between feed gas flowrate and the summation over the heavier components flowrates
vtk,j
foti,j = fiti,j ,
Lee, Pinto, Grossmann, & Park, 1996 modeled crude oil blending operations considering key component balances. In the present work, additional equations to enforce that the composition inside the tank is equal to the outlet stream composition have been added, Eqs. (18)–(21) (Li & Hui, 2002; Reddy, Karimi, & Srinivasan, 2004). All the tanks have a minimum-security inventory level and a maximum volume, Eq. (22) Vkt = vtk,j ∀k; t = 1, . . . , H (18)
i = NGI, NGII;
j = C3 H8 , C4 H10 , C4 H10 (iso), C5 H12 , C5 H12 (iso), C6 H14 ; t = 1, . . . , H Frgti = Fiti −
(8) j=N2 ,CH4
foti,j ,
i = NGI, NGII; t = 1, . . . , H (9)
The stream of liquids to the NGL fractionation plant is composed of carbon dioxide, ethane, propane, butanes, pentanes and hexane. Mass balances in each multicomponent streams splitter s are calculated as follows: fsptout,s,j = fsptin,s,j ∀j, s; t = 1, . . . , H (10) out ∈ Sout
Fsptout,s
=
Fsptin,s
∀s; t = 1, . . . , H
(11)
out ∈ Sout
j
xspts,j = 1 ∀s; t = 1, . . . , H
(12)
=
vok,j
+
fvnin,k,j
−
n=1
fvnout,k,j
∀j, k;
n=1
t = 1, . . . , H
(17)
j
vtk,j = Vkt xvtk,j ∀j, k; t = 1, . . . , H fvtout,k,j ∀k; t = 1, . . . , H Fvtout,k =
(19) (20)
j
fvtout,k,j = Fvtout,k xvtout,k,j VkL ≤ Vkt ≤ VkU
∀j, k; t = 1, . . . , H
∀k; t = 1, . . . , H
(21) (22)
After fractionation, there are storage tanks for natural gas plants products: propane (vC3 ), butane (vC4 ) and gasoline (vC5+ ). Each product is stored in a different tank v (i.e., there is one tank per individual product) therefore Eqs. (18)–(21)
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
are not required and Eq. (17) is replaced by (23) for only one component. These products are intermittently delivered by ship, one ship per product. These tanks have also capacity lower and upper limits Vvt
=
Vvo
t
+
Fvvnin,v
−
n=1
t
Fvvnout,v ,
n=1
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H
(23)
Ship delivery for these products is modeled as shown in Eqs. (20)–(28) (Lee et al., 1996) yatv = 1, v = vC3 , vC4 , vC5+ (24) t
t
t
t
ydtv = 1,
v = vC3 , vC4 , vC5+
(25)
t yatv = TAv ,
v = vC3 , vC4 , vC5+
(26)
t ydtv = TDv ,
v = vC3 , vC4 , vC5+
(27)
TDv − TAv ≥ tvmin v ,
v = vC3 , vC4 , vC5+
(28)
TDv − TAv ≤ tvmax v ,
v = vC3 , vC4 , vC5+
(29)
xwtv =
t n=1
yanv − ydnv ,
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H (30)
t Fvvtout,v ≤ Fvvt,U out,v xwv ,
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H
(31)
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H demvv ≥
t
from the storage tank v to the corresponding ship v. Total product delivery in the ship v cannot exceed its capacity (demvv ), Eq. (33), and penalties have been added to the objective function to account for the shortfall in the ship capacity (deltavv ), Eq. (34). 3.1.2. Ethylene plants Ethane from natural gas plants is fed to both ethylene plants, where ethylene is obtained by thermal cracking. Each plant can be fed from either natural gas plant, fulfilling contract conditions of type TOP (take or pay) with both plants. Fig. 3 shows the flowsheet of one of the ethylene plants. There are eight furnaces in parallel whose feed also includes an important ethane recycle stream. Main products from furnaces are hydrogen, methane, acetylene, ethylene, ethane, propylene, propane, butanes, butylene and pentanes. Cracked gas is compressed and acetylene is converted to ethylene in hydrogenation reactors prior to its entrance to the separation train. The first column in the separation train is the demethanizer: the top product stream is off-gas, which is used as fuel gas in furnaces, and the bottom stream is further separated to obtain pure ethylene, propane, propylene, butane, gasoline and an ethane stream that is partly recycled to the furnaces and partly sent to storage (ethane sphere). The two MINLP models in the present work only differ on the ethylene plant models, as described below. MINLP1 includes simplified linear models and MINLP2, semi-rigorous nonlinear models of the ethylene plants. • Linear model (included in MINLP1) Linear models have been derived for products as function of raw material flowrate, based on simulations with a semi-rigorous plant model (Schulz, Diaz, & Bandoni, 2000) (EPI and EPII stand for ethylene plants I and II, respectively): foti,j = αi,j Fiti + βi,j ,
t Fvvtout,v ≥ Fvvt,L out,v xwv ,
(32)
deltavv = demvv −
t
v = vC3 , vC4 , vC5+ Fvvtout,v ,
(33)
v = vC3 , vC4 , vC5+ (34)
Discontinuous delivery is modeled with binary variables yatv and ydtv to denote the time period t when that ship v starts or finishes loading product its product, respectively. Each ship loads the corresponding product v only once throughout the horizon, Eqs. (24) and (25), and it starts loading at t = TAv and finishes at t = TDv , Eqs. (26) and (27). Loading time is limited, Eqs. (28) and (29); xwtv is a 0–1 continuous variable that is equal to 1 if ship v is loading its product at time t, Eq. (30). Eqs. (31) and (32) are operating constraints on t,U t product transfer rate Fvvtout,v (Fvvt,L out,v ≤ Fvvout,v ≤ Fvvout,v )
i = EPI, EPII;
j = CH4 , C2 H4 , C3 H8 , C3 H6 , C4 H10 , C5 H12 ; t = 1, . . . , H
Fvvtout,v ,
1311
(35)
This model comprises six linear equations to calculate product flowrates per time period, per plant. • Nonlinear model (included in MINLP2) To obtain the profile of the operating variables at the ethylene plants when solving the production chain model, a semi-rigorous model which includes nonlinear correlations for each cracking furnace production (ffoti,h,j ) has been considered, Eq. (36); ethane recycle (Frti ), Eq. (37); furnaces h inlet pressure (Pfith,i ), Eq. (38); and simplified models for the distillation columns in the separation train, Eqs. (39) and (40). These variables are calculated as function of main plant optimization variables: hydrocarbon dilution ratio with steam (Rdti ), ethane conversion in furnace h in period t (Convti,h ), cracked gas compressor suction t ), and other variables as total furnaces pressure (PCGCi t charge, Ffii , ethylene/ethane ratio at the entrance of sepa-
1312
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Fig. 3. Typical ethylene plant (EPI).
ration train Rlti,j . Therefore, this process model comprises 139 equations per ethylene plant per time period (considering eight furnaces in each plant) while the linear model only comprises six equations per plant per time period
as follows:
ffoti,h,j = f (Ffiti , Rdti , Pfith,i , Convti,h ),
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
i = EPI, EPII; ∀j, h; t = 1, . . . , H
(36)
frti = f (ffoti,h,j , PCGC , Rlti , Su,i ), i = EPI, EPII; t = 1, . . . , H
(37)
t 0 = Cpoly + Cpoly
t
Fpinpoly −
t
n=1
Fponpoly ,
n=1
(41)
When inventory levels do not meet given storage targets (IOpoly ) an economic penalty (PENPOLYpoly ) is added to the objective function (Jackson et al., 2003) t PENPOLYtpoly ≥ IOpoly − Cpoly ,
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
(42)
Pfith,i = f (Ffiti , Rdti , Pfoth,i , Convti,h ), i = EPI, EPII; t = 1, . . . , H
(38)
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
t ffttu,i,j = f (PCGCi , Rlti , Su,i ),
i = EPI, EPII; ∀j, u; t = 1, . . . , H
(39)
(43)
Sales cannot exceed the daily forecast demand (demtpoly ) and the difference between the demand and the sales (deltatpoly ), Eq. (45) is penalized in the objective function demtpoly ≥ Fpotpoly ,
t Fbtu,i,j = f (PCGCi , Rlti , Su,i ),
i = EPI, EPII; ∀j, u; t = 1, . . . , H
t − IOpoly , PENPOLYtpoly ≥ Cpoly
(40)
Ethylene plant II is more efficient than ethylene plant I. Therefore ethylene plant II correlations have been adjusted to obtain higher production levels. 3.1.3. Polyethylene plants Ethylene is raw material for polyethylene plants: linear low density polyethylene (LLDPE), low density polyethylene (LDPE) and high density polyethylene (HDPE). Linear inlet–outlet models have been included, based on the information given in Table 2. Polyethylenes bulk storage is modeled
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
(44)
deltatpoly = demtpoly − Fpotpoly , poly = HDPE, LDPE, LLDPE;
t = 1, . . . , H
(45)
3.1.4. VCM and PVC plants These plants also use ethylene as raw material. Linear models based on the information given in Table 3 have been considered. Product daily demands have been considered.
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316 Table 3 NaOH, VCM, EDC and PVC yields
Table 5 Economical data
Product
t of product/t of C2 H4
NaOH VCM EDC PVC
0.0573 0.0093 0.0087 0.0653
Table 4 Ammonia and urea yields Product
kmol of product/kmol of methane
NH3 (NH2 )2 CO
0.2329 0.9506
3.1.5. Ammonia and urea plants Residual gas (99% methane) from natural gas plant I is used as raw material for the ammonia plant. Methane and air in stoichiometric ratio are used to obtain nitrogen and hydrogen for ammonia synthesis. Carbon dioxide is also obtained and is later used in the urea production. Ammonia is used in urea production and it is also sold as final product. Linear correlations have been derived for both plants production, as function of natural gas consumption based on data from Table 4. Discontinuous delivery of urea and ammonia by ship has been modeled with equations analogous to Eqs. (20)–(28). They are delivered by train and truck according to a forecast daily demand for each product Vvt = Vvo +
t n=1
Fvvnin,v −
t
1313
Fvvnout,v −
n=1
t
Fvtnout,v ,
Natural gas (US$/kmol) Residual gas (US$/kmol) C2 H4 (US$/kmol) C2 H6 (US$/kmol) C3 H8 (US$/kmol) C4 H10 (US$/kmol) C5 H12 (US$/kmol) NH3 (US$/kmol) (NH2 )2 CO (US$/kmol) LDPE (US$/t) HDPE (US$/t) LLDPE (US$/t) NaOH (US$/t) VCM (US$/t) EDC (US$/t) PVC (US$/t) Inventory cost (US$/kmol) Electricity (kWh)
Ethylene plants: Linear correlations have been derived from rigorous plant simulations for required fuel gas for furnaces (both plants) and boilers (only ethylene plant I) as function of ethane flowrate. Remaining plants: Utilities consumption for polyethylenes, VCM, PVC, ammonia and urea plants have been taken from the literature. They are shown in Table 5. Tank inventory costs (ethane, propane, butane, gasoline, ammonia and ethylene) are calculated according to the trapezoidal area, where Civv is the cost coefficient. Urea storage is also handled this way
n=1
v = vNH3 , vurea; t = 1, . . . , H
(46)
There are capacity lower and upper limits for ammonia tank and urea storage room. • Objective function: The objective function is the maximization of net profit along a 20 days horizon, defined as the difference of the incomes from product sales and the costs (operating costs, raw materials costs and tank inventory costs) and the penalties for not satisfying forecast demands, ship capacities and inventory targets profit = sales revenue − raw materials costs −operating costs − Cinv − CIO − Cdem
CInv =
• Operating costs: Natural gas plants: the main operating costs are associated to compression loads. Compressors are driven by gas turbines. In natural gas plant I, compression (33–60 bar) and recompression (18–35 bar) are required. In natural gas plant II, not far from gas wells, no inlet gas compression is required. Linear models for the gas consumption in the compressor turbines have been derived from rigorous simulations as function of raw material flow.
Civv
v
V t + V t−1 k k 2 t
(48)
• Penalty terms: Penalties on deviations from inventory targets for polyethylenes and for not meeting forecast demand (or ships capacity) for the remaining products have been included t θpoly (49) PenPolypoly CIO = t
poly
Cdem = (47)
1.8123 1.63863 11.73 4.1253 7.7625 10.23236 12.70224 2.805 9.9 640 720 600 230 400 240 1080 0.05 0.09454
j
ξj
t
deltatj
(50)
4. Discussion of results Two multiperiod mixed integer nonlinear programming mathematical models have been formulated for the production chain optimization of a petrochemical complex: MINLP1 includes Eqs. (1)–(3), (8)–(35) and (41)–(50) and MINLP2 comprises Eqs. (1)–(3), (8)–(34) and (36)–(50). The main difference between these models is the ethylene plant
1314
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Fig. 4. Optimal ethylene production profile in ethylene plant II with alternative MINLP models.
model which is simplified and linear in MINLP1 and semirigorous nonlinear in MINLP2. An approximate MILP model has been also formulated replacing bilinear constraints in the MINLPs by valid linearisations (McCormick, 1976), Eqs. (4)–(7), to provide a valid initial point for the MINLPs. Main variables profiles of the MINLPs have been compared, as it is shown in Figs. 4 and 5. Although they show a similar trend, the inclusion of semi-rigorous models for both ethylene plants (MINLP2) renders more uniform production profiles than the obtained with the linear models (MINLP1). In ethylene plant II, total ethylene production linear model has an average error of 6.6% (Fig. 4). Average errors of 8% can be seen in Fig. 5 for the ammonia plant. Linear models for plants have been derived for operating conditions not far from actual ones. Consequently, errors introduced when using linear models are mainly due to process flowrates that are out of range. This could have been improved by replacing linear models by piecewise linear functions. However, as nonlinear models have shown a good performance, further analysis has been made with these ones. Regarding natural gas plant II, propane, butanes and gasoline tanks profiles and the arrival days for the corresponding ships are shown in Fig. 6. After product delivery by ship, the plant production decreases. Natural gas plant I operates at its highest production level throughout the entire time horizon (demethanizer top pressure at its lower bound, 17 bar, and highest ethane recovery). Both ethylene plants require maximum allowed ethane amount from this plant, mainly due to ethane lower price.
Fig. 6. Propane, butanes and gasoline tanks profiles in natural gas plant II.
Figs. 7 and 8 show total furnaces charge and conversion profiles for three of the eight furnaces in ethylene plant II. A slight increase in conversion during the last 4 days is associated to the decrease in total furnaces charge, as it is always seen in actual plant behaviour. As an example of polyethylene plants, Fig. 9 shows optimal LLDPE production profile and inventory level, which most of the time meets the inventory target (360 t, horizontal line). Urea production profile, inventory level and amount daily delivered by trains and trucks are plotted in Fig. 10. Maximum stored amount is 45,150 t and the ship arrival day is 18. Urea amount daily delivered by trains and trucks always meets constant daily demands (900 t). Finally, a different scenario has been analysed with MINLP2 model (Case 2). A 20% increase in HDPE has been specified. Fig. 11 shows that a slight increase in ethylene production throughout the time horizon is enough to fulfil
Fig. 7. Furnaces total feed and ethane recycle profiles in ethylene plant II.
Fig. 5. Optimal urea production profile with alternative MINLP models.
Fig. 8. Furnaces conversion profiles in ethylene plant II.
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
1315
Table 6 Computational details
Constraints Continuous variables Binary variables CPU (s) Solver Objective function (US$)
MINLP1
MINLP2
6932 6123 200 724 DICOPT++ (CONOPT2/OSL) 28624443
16373 18408 200 5002 DICOPT++ (CONOPT2/OSL) 30014774
the increased raw material demand. HDPE forecasted daily demand and HDPE production profile are plotted in Fig. 12, showing that the new demand is satisfied most of the days with an overall increase in HDPE production of 14%.
Fig. 12. Comparison between HDPE production in original Case 1 and Case 2 (20% increase in HDPE demand).
5. Computational details Fig. 9. LLDPE production and inventory profiles. Inventory target at 360 t (horizontal line).
Multiperiod MILP and MINLPs models are coded in GAMS 2.25 (Brooke, Kendrick, & Meeraus, 1992) modelling environment and solved with OSL and DICOPT++ (CONOPT2 and OSL), respectively. MINLP models statistics are shown in Table 6 and both models have required three major iterations to achieve convergence. The approximate MILP model, solved to obtain a valid initial point for MINLP, comprises 13,376 constraints with the same number of continuous and binary variables as MINLP1 model. The problems have been solved in a 1 GHz Pentium III PC, with 256 MB of RAM.
6. Conclusions Fig. 10. Urea production profile, inventory level and amount delivered by trains and trucks.
Fig. 11. Comparison between ethylene production in ethylene plant II in original MINLP2 problem (Case 1) and with a 20% increase in HDPE daily demand (Case 2).
Current petrochemical sites tend to be single integrated facilities to optimize total processing economics for converting natural gas to final industrial products. A model that represents the entire supply chain can provide a better insight on global profit. Two production chain optimization models (MILP1 and MINLP2) of increasing degree of complexity have been formulated. Their solutions provide coordination of responses to demands, production planning, production distribution and inventory level management. The most detailed MINLP solution also determines main plants operating conditions and optimization variables such as conversion in each cracking furnace. Two scenarios corresponding to different product demands have been analysed. An optimization approach in which all the processes are highly integrated has proved to be a useful tool to detect
1316
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
operational and economical improvements. The inclusion of more detailed process models, the exploration of solution methods that can cope with the resulting large-scale MINLP and the inclusion of uncertainty in several process parameters are part of current work.
References Al-Sharrah, G. K., Alatiqi, I., Elkamel, A., & Alper, E. (2001). Planning an integrated petrochemical industry with an environmental objective. Industrial and Engineering Chemical Research, 40, 2103–2111. Brooke, A., Kendrick, D., & Meeraus, A. (1992). GAMS—A user’s guide, release 2.25. Danvers, MA: Boyd and Frase Publ. Co. Bok, J., Grossman, I. E., & Parks, S. (2000). Supply chain optimisation in continuous flexible process networks. Industrial and Engineering Chemical Research, 39, 1279–1290. Diaz, M. S., & Bandoni, A. (1996). A mixed integer optimization strategy for a real large scale chemical plant. Computers and Chemical Engineering, 20(5), 531–545. Diaz, M. S., Serrani, A., de Beistegui, R., & Brignole, E. (1995). An MINLP strategy for the debottlenecking problem in an ethane extraction plant. Computers and Chemical Engineering, 19, 175– 178. Furlonge, H. I., & Young-Hoon, A. (2000). Towards the development of an optimal long-term structure and policy for the development of Trinidad and Tobagoˇıs petrochemical industry. Part I. The methanebased complex. West Indian Journal of Engineering, 22(2). Jackson, J. R., & Grossmann, I. E. (2003). Temporal decomposition scheme for nonlinear multisite production planning and distribution models. Industrial and Engineering Chemical Research, 42, 3045–3055. Jackson, J., Hofmann, J., Wassick, J., & Grossmann, I. E. (2003). A nonlinear multiperiod process optimization model for production planning in multi-plant facilities. In I. E. Grossmann & C. M. McDonald (Eds.).
Fourth International Conference on Foundations of computer Aided Process Operations (FOCAPO), 281–284. Li, W., & Hiu, C.-W. (2002). Scheduling crude oil unloading. storage and processing. Industrial and Engineering Chemical Research, 41, 6723–6734. Lee, H., Pinto, J. M., Grossmann, I. E., & Park, S. (1996). Mixed-integer linear programming model for refinery short-term scheduling of crude oil unloading with inventory management. Industrial and Engineering Chemical Research, 35, 1630–1641. McCormick, G. P. (1976). Computability of global solutions to factorable nonconvex programs. Part I. Convex underestimating problems. Mathematical Programming, 10, 146–175. Neiro, S. M. S., & Pinto, J. M. (2003). Supply chain optimization of petroleum refinery complexes. In I. E. Grossmann & C. M. McDonald (Eds.), Foundations of computer-aided process operations (pp. 59–72). Florida: CACHE AICHE INFORMS. Quesada, I., & Grossmann, I. E. (1995). Global optimization of bilinear process networks with multicomponent flows. Computers and Chemical Engineering, 19(12), 1219–1242. Reddy, P. C. P., Karimi, I. A., & Srinivasan. (2004). A new continuoustime formulation for scheduling crude oil operations. Chemical Engineering Science, 59, 1325–1341. Schulz, E., Diaz, M. S., & Bandoni, A. (2000). Interaction between process plant operation and cracking furnaces maintenance policy in an ethylene plant. Computer Aided and Process Engineering, 8, 487–492. Shah, N. (2004). Process industry supply chains: Advances and challenges. In A. Barbosa-P´ovoa & H. Matos (Eds.). Computer Aided and Process Engineering, 18, 123–138. Turkay, M., & Grossmann, I. E. (1996). Disjunctive programming techniques for the optimization for process systems with discontinuous investment costs—Multiple size regions. Industrial and Engineering Chemical Research, 35, 2611–2623. Turkay, M., Asakura, T., Fujita, K., Hui, C. W., & Natori Y. (1998). Total site optimization of a petrochemical complex. In J. F. Pekny, & G. Blau (Eds.), Foundations of computer-aided process operations, AIChE symposium series no. 320: vol. 94 (pp. 185–189) Production editor, CACHE Publications, B. Carnahan.
Supply chain optimization of large-scale continuous processes E.P. Schulz, M.S. Diaz ∗ , J.A. Bandoni Planta Piloto de Ingenier´ıa Qu´ımica—PLAPIQUI (Universidad Nacional del Sur—CONICET), Camino La Carrindanga Km 7, 8000 Bah´ıa Blanca, Argentina Available online 11 April 2005
Abstract In this work, a supply chain model for a petrochemical complex is formulated as a multiperiod mixed integer nonlinear program (MINLP). The model includes production, product delivery, inventory management and decisions such as individual production levels for each product, as well as operating conditions for each plant in the complex. Two multiperiod MINLP are formulated in order to compare models with different levels of rigorousness. The first one includes linear mathematical models derived for the ethylene plants while the second considers semi-rigorous nonlinear models. Simplified models based on rigorous existing models tuned with actual plant data are considered for the natural gas plants, taking into account variations in production with key plant operating variables. Data for chemical transformations and utilities consumption from the literature have been used to model the remaining plants. © 2005 Published by Elsevier Ltd. Keywords: Multiperiod MINLP; Short term planning; Petrochemical complex
1. Introduction Petrochemical industries are within wider supply chains and involve a number of players, which span sourcing, manufacturing and distribution. Every day, much money is spent on feedstocks, chemicals, utilities, equipment, labor and effluent treatment in order to produce marketable process-industry products on time. In a petrochemical complex there may be several plants operating under different conditions producing several categories of petrochemicals (basic, intermediate and end-product). The responses to forecast market demands have to be coordinated and simultaneous planning of production and each plant production distribution has to be undertaken. Decisions have to be taken at the different stages of the chain, i.e. operating conditions of individual processes, intermediate product handling and end-product storage and distribution. A model that can account for the different interactions among the plants to provide a suitable production and operation plan for the multi-plant facility is a valuable tool for the decision-making process. Turkay, Asakura, Fujita, Hui, and Natori (1998) applied logic-based approaches proposed by Turkay and Grossmann ∗
Corresponding author.
0098-1354/$ – see front matter © 2005 Published by Elsevier Ltd. doi:10.1016/j.compchemeng.2005.02.025
(1996) to solve a mixed integer linear problem for the total site optimization of a petrochemical complex. Furlonge and Young-Hoon (2000) presented a methodology to address long-term planning for development of the industry in Trinidad and Tobago, selecting an initial structure and policy for the development of an olefin-based complex. Al-Sharrah, Alatiqi, Elkamel, and Alper (2001) proposed an environmental objective in the production planning in the petrochemical industry. Bok, Grossman, and Parks (2000) proposed a bilevel decomposition strategy for the solution of a multiperiod MILP model for the supply chain of continuous processes taking into account sales, intermittent deliveries, production shortfalls, delivery delays, inventory profiles and changeovers. As stated by Shah (2004), supply chain planning and scheduling is a tool to optimize production, distribution and storage resources of a fixed network to respond properly to the external conditions, such as orders and demands forecasts, over a short to medium term horizon. Production process representation depends on the gross margin of the business. Business with low margin facilities, such as refining and petrochemicals, require “property based” representations, which use process conditions and process models. An example of property based planning is provided by Jackson and
1306
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Nomenclature (a) Sets h furnace i plant (natural gas plants: NGI and NGII, ethylene plants: EPI and EPII) in inlet stream to a tank, mixer or splitter j component k intermediate tank m mixer out outlet stream from a tank, mixer or splitter poly polyethylenes (HDPE: high density polyethylene, LDPE: low density polyethylene, LLDPE: linear low density polyethylene) s splitter t, n time intervals (from 1 to H, number of time periods) v {vC3 , vC4 , vC5+ } tanks and ships for propane (vC3 ), butane (vC4 ), gasoline (vC5+ ), ammonia (vNH3 ) or urea (vurea) delivery (b) Parameters ai,j (K−1 ), bi,j (bar−1 ) and ci,j correlation coefficients for natural gas plants recoveries inventory cost (US$/kmol) Civv 0 Cpoly initial inventory of polyethylene poly (t) demvv ship v capacity t,U fit,L i,j –fii,j minimum–maximum individual molar flow of component j in plant i in period t (kmol/d) t,U Fvvt,L –Fvv out,v out,v minimum–maximum transfer rate of propane, butane or gasoline from tank v to ship v in period t (kmol/d) IOpoly storage target of polyethylene poly (t) Su,i separation factor in unit u in the separation train in plant i max minimum–maximum loading duration (d) tvmin –tv v v period when ship v starts loading its product TAv (d) v0k,j initial moles of component j in tank k (kmol) VkL –VkU minimum–maximum levels in intermediate tank k (kmol) Vv0 initial moles in tank k (only one component) (kmol) xii,j inlet molar fraction of component j in plant i Greek letters αi,j , βi,j correlation coefficients for ethylene plants simplified models ξj penalty coefficient for product p unmet demand t,U ηt,L –η i,j i,j minimum–maximum component j recovery in demethanizer bottom stream in plant i in period t
θ poly
penalty coefficient for polyethylene poly unmet inventory target (US$/t)
(c) Binary variables yatv denotes when ship v starts loading its product ydtv denotes when ship v finishes loading its product (d) Continuous variables Cdem penalty for not meeting forecast demand Cinv total inventory cost (US$) CIO penalty on deviations from inventory targets for polyethylenes Convth,i ethane conversion in furnace h in period t t inventory of polyethylene poly in period t (t) Cpoly deltav shortfall in ship v capacity deltatpoly demand shortfall for polyethylene poly in period t (t) demtpoly forecast demand for polyethylene poly in period t (t) ffoti,h j production of component j in furnace h in plant i in period t (kmol/d) fiti,j individual molar inlet molar flow of component j in plant i in period t (kmol/d) fmxtin,m,j –fmxtout,m,j individual molar inlet–outlet molar flow of component j in mixer m in period t (kmol/d) individual molar outlet molar flow of compofoti,j nent j in plant i in period t (kmol/d) Fpotpoly sales of polyethylene poly in period t (t/d) fsptin,s,j –fsptout,s,j individual molar inlet–outlet molar flow of component j in splitter s in period t (kmol/d) ftbti,u,j molar flowrate of component j as bottom product in unit u in the separation train in plant i (kmol/d) fttti,u,j molar flowrate of component j as top product in unit u in the separation train in plant i (kmol/d) fvtin,k,j –fvtout,k,j individual molar inlet–outlet molar flow of component j in intermediate tank k in period t (kmol/d) Ffiti total furnaces charge in plant i in period t (kmol/d) Fiti total molar inlet molar flow in plant i in period t (kmol/d) Fmxtin,m –Fmxtout,m total molar inlet–outlet molar flow in mixer m in period t (kmol/d) Fpitpoly production of polyethylene poly in period t (t/d) Frti ethane recycle stream in plant i in period t (kmol/d)
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Frgti
residual gas total flow in plant i in period t (kmol/d) Fsptin,s –Fsptout,s total molar inlet–outlet molar flow in splitter s in period t (kmol/d) Fvtin,k –Fvtout,k total molar inlet–outlet molar flow in intermediate tank k in period t (kmol/d) Fvttout,v outlet molar flow for train or truck continuous delivery of product v (only one component in the stream) in tank v in period t (kmol/d) Fvvtin,v –Fvvtout,v inlet–outlet molar flow for ship discontinuous delivery of product v (only one component in the stream; v = vNH3 , vurea) in tank v in period t (kmol/d) Pit demethanizer top pressure in plant i in time t (bar) t PCGCi cracked gas compression suction pressure in plant i in period t PENPOLYtpoly unmet storage target for polyethylene poly in period t (t) Pfith,i furnace h in plant i inlet pressure in period t (bar) Pfoth,i furnace h in plant i outlet pressure in period t (bar) Rdti hydrocarbon dilution ratio with steam in plant i in period t Rlti ethylene/ethane ratio at the entrance of the separation train in plant i in period t Tit high pressure separation tank in plant i in time t (K) TDv period when ship v finishes loading its product vtk,j moles of component j in intermediate tank k in period t (kmol) Vkt total moles in intermediate tank k in period t (kmol) xspts,j molar fraction of component j in splitter s in period t xvtk,j molar fraction of component j in intermediate tank k in period t xwtv 0–1 continuous variable to denote if ship v is loading its product at time t Greek letter ηti,j component j recovery in demethanizer bottom stream in plant i in period t
1307
benefits of solving for the complex together than for the individual refineries separately. This work presents an integrated multiperiod optimization model for the supply chain of an existing large-scale petrochemical complex with short grid spacing. The mathematical model accounts for supply, production and product multimodal delivery. The petrochemical complex under study comprises two natural gas liquid processing (NGL) plants, two ethylene plants, a chlorine plant, a VCM plant, a PVC plant, two polyethylene plants, an ammonia and an urea plant. Mathematical models have been derived for the NGL and ethylene plants, based on rigorous existing models tuned with actual plant data (Diaz & Bandoni, 1996; Diaz, Serrani, de Beistegui, & Brignole, 1995). Simplified models take into account variations in production with key plant operating variables, such as temperature and pressure in separation units or conversion in reactors. An alternative nonlinear semi-rigorous model for ethylene plants was formulated in order to show the benefits in solving more rigorous models. Available yield data for chemical transformations have been used to model the rest of the petrochemical complex. Multiperiod modeling allows economic parameters, such as demands, to vary with time. The objective function is the maximization of total profit for the entire site, considering economic issues and drivers such as contractual obligations, operating costs and inventory costs subject to constraints on mass balances, bounds on product demands, equipment capacities and conditions, yields and intermediate and final product storage tanks limitations. The most profitable set of product-stock to produce and sell is determined. Final products are intermittently distributed by ships and continuously by trains and trucks while satisfying storage tanks capacities. Discontinuous product delivery is represented by the use of binary variables. Bilinear nonconvex constraints arise from products between inlet flows and product calculated recovery in natural gas plants and from blending and storage of multicomponent streams. The resulting nonconvex large scale mixed integer nonlinear models have been solved with GAMS using an initial point provided by the solution of a linear model (MILP) obtained applying a reformulation–linearisation technique which provides a valid relaxation of the original feasible region. Numerical results provide a realistic insight on supply chain management in the petrochemical complex.
2. Petrochemical complex description Grossmann (2003) who proposed Lagrangean decomposition (spatial and temporal) techniques for the solution of a non-linear programming (NLP) problem that models multisite production planning of production, transportation and sales in a chemical company. Neiro and Pinto (2003) describe a planning model of a refinery complex as a mixed integer nonlinear programming (MINLP) problem, including process unit models and showing significant cost
The petrochemical complex subject of this study has an annual production of 4,800,000 t. Fig. 1 shows its model representation. There are two natural gas processing plants, whose main objective is to extract ethane from natural gas to use it as raw material in the ethylene plants. natural gas plant I, next to the complex, is fed with 24 MM m3 /d of natural gas. Residual gas (mainly methane) is recompressed to pipeline pressure and part of it is taken as feed for the ammonia plant
1308
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Fig. 1. Petrochemical complex model representation.
and the rest is delivered as sales gas. Pure ethane, propane, butane and gasoline are plant products. Natural gas plant II has its cryogenic sector (referred to as demethanizing plant) several kilometers away from the conventional separation train (NGL fractionation plant). The demethanizing plant is fed with 36 MM m3 /d of natural gas. Light gases are separated from the heavy ones (ethane, propane, butanes and gasoline) and injected to the natural gas pipeline. The rich liquefied mixture (5 MM m3 /d of heavy gases) is stored in thermal vessels and pumped to the petrochemical complex where it is fed to containers to equalize the charge, i.e. to damp any pulsation or flow changes that may occur anywhere along the pipeline. The mixture undergoes a conventional distillation train in the NGL fractionation plant to obtain LPG (liquefied petroleum gas: propane, butane and gasoline) and ethane. Ethylene plants process 2300 t/y of pure ethane; ethylene is provided as raw material to polyethylene and VCM plants and the rest is exported. Ethylene plant II is more efficient than ethylene plant I. The ammonia plant production is 2050 t/d, most of which is fed to the urea plant to produce 3250 t/d of urea. Urea, ammonia, polyethylenes and PVC are delivered by ship, train and trucks.
3. Mathematical model The optimization of the supply chain of a petrochemical complex is formulated as a multiperiod MINLP problem resulting from a uniform time discretisation with 1-day time intervals. Binary variables are associated to intermittent product delivery. Two alternative MINLP models have been formulated, with different degree of rigorousness. The first one is formulated with linear models representing the plants, where nonlinearities arise from bilinear equations standing for mixers and splitters. The second one adds nonlinear process models for some plants and operating conditions. The objective is to compare their performance and decide their application depending on the amount of information required. An approximate mixed integer linear model (MILP) has also been formulated applying linearisation techniques to bilinear equations (Quesada & Grossmann, 1995) to obtain a valid initial point required to solve the MINLPs. The objective function is the maximization of total profit, defined as the difference between sales revenue and the total operating cost plus penalties for not meeting demands and inventory targets during the given time horizon. A brief
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316 Table 1 Feed gas composition for natural gas plants Component
Natural gas plant I
Natural gas plant II
CO2 N2 CH4 C2 H6 C3 H8 C4 H10 C4 H10 (iso) C5 H12 C5 H12 (iso) C6 H14
0.0065 0.0144 0.9043 0.0461 0.0176 0.0033 0.0044 0.0015 0.0009 0.0010
0.0212 0.0082 0.9032 0.0381 0.0161 0.0030 0.0049 0.0013 0.0015 0.0025
Table 2 Polyethylenes’ yields Plant
t of product/kmol of C2 H4
LDPE LLDPE HDPE
0.0272 0.0275 0.0273
description of plants and their mathematical models is given below. 3.1. Plants description and models 3.1.1. Natural gas separation plants In both natural gas processing plants, pipeline gas composition has been considered as constant (Tables 1 and 2). The total molar flowrate is a variable and the individual flows have been calculated as (NGI and NGII stand for natural gas plants I and II, respectively): fiti,j = Fiti xii,j
i = NGI, NGII; ∀j; t = 1, . . . , H
(1)
1309
Natural gas plants are fed from several pipelines. The main difference in the composition of the natural gas plants is the carbon dioxide content, as it is shown in Table 1, requiring basic (NGI) and gas subcooled (NGII) process turboexpansion technologies. The cryogenic sector is the core of these plants. In the basic turboexpansion process, inlet gas is cooled by heat exchange with residual gas and demethanizer side and bottom reboilers. The partially condensed gas feed is then sent to a high-pressure separator. The vapor is expanded through a turboexpander to obtain the low temperatures required for high ethane recovery and is then fed to the top of a demethanizer column. The liquid from the high-pressure separator enters the demethanizer at its lowest feed point. Methane and lighter components, the top product, constitute the residual gas, which are re-injected to the pipeline. Carbon dioxide and ethane distribute between top and bottom streams. Heavier hydrocarbons are obtained as bottom product and are further fractionated to obtain pure ethane, propane, butane and gasoline. Fig. 2 shows the cryogenic sector of a basic turboexpansion plant. In a gas subcooled turboexpansion process, a fraction of the vapor from the cold tank is condensed and subcooled by heat exchange with residual gas coming from the demethanizer. The remaining vapor is expanded through the turboexpander and fed to the middle of the column. The subcooled liquid is flashed and fed to the top of the demethanizer column. Linear models have been derived for ethane (ηti,C2 H6 ) and carbon dioxide (ηti,CO2 ) recovery in demethanizer bottom stream, as function of main operating variables in the cryogenic sector: high pressure separation tank temperature and demethanizing column top pressure. These correlations have been obtained based on simulations with a plant rigorous model (Diaz et al., 1995) for conditions not far from actual operating points, with correlation coefficients ranging from
Fig. 2. Cryogenic sector in a basic turboexpansion process.
1310
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
fsptout,s,j = Fsptout,s xspts,j
0.96 to 0.99 ηti,j
=
ai,j Tit
+ bi,j Pit
+ ci,j ,
ai,j , bi,j < 0;
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
(13) (2)
Product flows (foti,j ) are calculated as products between inlet flowrate (fiti,j ) and recovery (ηti,j ), resulting in bilinear
equations
foti,j = fiti,j ηti,j ,
∀j, s, out ∈ Sout ; t = 1, . . . , H
i = NGI, NGII;
j = CO2 , C2 H6 ; t = 1, . . . , H
(3)
In the MILP model, these equations have been replaced by (Quesada & Grossmann, 1995): t,L t t,L t,L t foti,j ≥ fit,L i,j ηi,j + ηi,j fii,j − fii,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
(4)
(5)
(6)
t,U t t,L t,U t foti,j ≤ fit,L i,j ηi,j + ηi,j fii,j − fii,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
(14)
where Sout is the set of outlet streams and component j represents carbon dioxide, ethane, propane, butanes, pentanes and hexane. In the linear model, Eqs. (13) and (14) have been replaced by valid linearisations analogous to Eqs. (4)–(7) (Quesada & Grossmann, 1995). Mass balances in mixers with multicomponent streams are modeled as follows: fmxtin,m,j = fmxtout,m,j ∀j, m; t = 1, . . . , H (15)
Fmxtin,m = Fmxtout,m
∀j, m; t = 1, . . . , H
(16)
in ∈ Si
t,L t t,U t,L t foti,j ≤ fit,U i,j ηi,j + ηi,j fii,j − fii,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
∀j, s; t = 1, . . . , H
in ∈ Si
t,U t t,U t,U t foti,j ≥ fit,U i,j ηi,j + ηi,j fii,j − foi,j ηi,j ,
i = NGI, NGII; j = CO2 , C2 H6 ; t = 1, . . . , H
fsptin,s,j = Fsptin,s xspts,j
(7)
where Sin is the set of inlet streams to mixer m and j represent the same components as for the splitter equations. There are multicomponent storage tanks in natural gas plant II and charge equalizers at the entrance of the fractionation plant. The material balance equation for each intermediate container at each discrete interval is described as the initial volume plus the summation of inflows minus the summation of outflows up to each time interval (Lee, Pinto, Grossman, & Park, 1996): t
t
For propane, butanes, pentanes and hexane, 100% recovery is assumed. Residual gas flowrate (mainly methane and nitrogen) is calculated as the difference between feed gas flowrate and the summation over the heavier components flowrates
vtk,j
foti,j = fiti,j ,
Lee, Pinto, Grossmann, & Park, 1996 modeled crude oil blending operations considering key component balances. In the present work, additional equations to enforce that the composition inside the tank is equal to the outlet stream composition have been added, Eqs. (18)–(21) (Li & Hui, 2002; Reddy, Karimi, & Srinivasan, 2004). All the tanks have a minimum-security inventory level and a maximum volume, Eq. (22) Vkt = vtk,j ∀k; t = 1, . . . , H (18)
i = NGI, NGII;
j = C3 H8 , C4 H10 , C4 H10 (iso), C5 H12 , C5 H12 (iso), C6 H14 ; t = 1, . . . , H Frgti = Fiti −
(8) j=N2 ,CH4
foti,j ,
i = NGI, NGII; t = 1, . . . , H (9)
The stream of liquids to the NGL fractionation plant is composed of carbon dioxide, ethane, propane, butanes, pentanes and hexane. Mass balances in each multicomponent streams splitter s are calculated as follows: fsptout,s,j = fsptin,s,j ∀j, s; t = 1, . . . , H (10) out ∈ Sout
Fsptout,s
=
Fsptin,s
∀s; t = 1, . . . , H
(11)
out ∈ Sout
j
xspts,j = 1 ∀s; t = 1, . . . , H
(12)
=
vok,j
+
fvnin,k,j
−
n=1
fvnout,k,j
∀j, k;
n=1
t = 1, . . . , H
(17)
j
vtk,j = Vkt xvtk,j ∀j, k; t = 1, . . . , H fvtout,k,j ∀k; t = 1, . . . , H Fvtout,k =
(19) (20)
j
fvtout,k,j = Fvtout,k xvtout,k,j VkL ≤ Vkt ≤ VkU
∀j, k; t = 1, . . . , H
∀k; t = 1, . . . , H
(21) (22)
After fractionation, there are storage tanks for natural gas plants products: propane (vC3 ), butane (vC4 ) and gasoline (vC5+ ). Each product is stored in a different tank v (i.e., there is one tank per individual product) therefore Eqs. (18)–(21)
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
are not required and Eq. (17) is replaced by (23) for only one component. These products are intermittently delivered by ship, one ship per product. These tanks have also capacity lower and upper limits Vvt
=
Vvo
t
+
Fvvnin,v
−
n=1
t
Fvvnout,v ,
n=1
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H
(23)
Ship delivery for these products is modeled as shown in Eqs. (20)–(28) (Lee et al., 1996) yatv = 1, v = vC3 , vC4 , vC5+ (24) t
t
t
t
ydtv = 1,
v = vC3 , vC4 , vC5+
(25)
t yatv = TAv ,
v = vC3 , vC4 , vC5+
(26)
t ydtv = TDv ,
v = vC3 , vC4 , vC5+
(27)
TDv − TAv ≥ tvmin v ,
v = vC3 , vC4 , vC5+
(28)
TDv − TAv ≤ tvmax v ,
v = vC3 , vC4 , vC5+
(29)
xwtv =
t n=1
yanv − ydnv ,
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H (30)
t Fvvtout,v ≤ Fvvt,U out,v xwv ,
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H
(31)
v = vC3 , vC4 , vC5+ ; t = 1, . . . , H demvv ≥
t
from the storage tank v to the corresponding ship v. Total product delivery in the ship v cannot exceed its capacity (demvv ), Eq. (33), and penalties have been added to the objective function to account for the shortfall in the ship capacity (deltavv ), Eq. (34). 3.1.2. Ethylene plants Ethane from natural gas plants is fed to both ethylene plants, where ethylene is obtained by thermal cracking. Each plant can be fed from either natural gas plant, fulfilling contract conditions of type TOP (take or pay) with both plants. Fig. 3 shows the flowsheet of one of the ethylene plants. There are eight furnaces in parallel whose feed also includes an important ethane recycle stream. Main products from furnaces are hydrogen, methane, acetylene, ethylene, ethane, propylene, propane, butanes, butylene and pentanes. Cracked gas is compressed and acetylene is converted to ethylene in hydrogenation reactors prior to its entrance to the separation train. The first column in the separation train is the demethanizer: the top product stream is off-gas, which is used as fuel gas in furnaces, and the bottom stream is further separated to obtain pure ethylene, propane, propylene, butane, gasoline and an ethane stream that is partly recycled to the furnaces and partly sent to storage (ethane sphere). The two MINLP models in the present work only differ on the ethylene plant models, as described below. MINLP1 includes simplified linear models and MINLP2, semi-rigorous nonlinear models of the ethylene plants. • Linear model (included in MINLP1) Linear models have been derived for products as function of raw material flowrate, based on simulations with a semi-rigorous plant model (Schulz, Diaz, & Bandoni, 2000) (EPI and EPII stand for ethylene plants I and II, respectively): foti,j = αi,j Fiti + βi,j ,
t Fvvtout,v ≥ Fvvt,L out,v xwv ,
(32)
deltavv = demvv −
t
v = vC3 , vC4 , vC5+ Fvvtout,v ,
(33)
v = vC3 , vC4 , vC5+ (34)
Discontinuous delivery is modeled with binary variables yatv and ydtv to denote the time period t when that ship v starts or finishes loading product its product, respectively. Each ship loads the corresponding product v only once throughout the horizon, Eqs. (24) and (25), and it starts loading at t = TAv and finishes at t = TDv , Eqs. (26) and (27). Loading time is limited, Eqs. (28) and (29); xwtv is a 0–1 continuous variable that is equal to 1 if ship v is loading its product at time t, Eq. (30). Eqs. (31) and (32) are operating constraints on t,U t product transfer rate Fvvtout,v (Fvvt,L out,v ≤ Fvvout,v ≤ Fvvout,v )
i = EPI, EPII;
j = CH4 , C2 H4 , C3 H8 , C3 H6 , C4 H10 , C5 H12 ; t = 1, . . . , H
Fvvtout,v ,
1311
(35)
This model comprises six linear equations to calculate product flowrates per time period, per plant. • Nonlinear model (included in MINLP2) To obtain the profile of the operating variables at the ethylene plants when solving the production chain model, a semi-rigorous model which includes nonlinear correlations for each cracking furnace production (ffoti,h,j ) has been considered, Eq. (36); ethane recycle (Frti ), Eq. (37); furnaces h inlet pressure (Pfith,i ), Eq. (38); and simplified models for the distillation columns in the separation train, Eqs. (39) and (40). These variables are calculated as function of main plant optimization variables: hydrocarbon dilution ratio with steam (Rdti ), ethane conversion in furnace h in period t (Convti,h ), cracked gas compressor suction t ), and other variables as total furnaces pressure (PCGCi t charge, Ffii , ethylene/ethane ratio at the entrance of sepa-
1312
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Fig. 3. Typical ethylene plant (EPI).
ration train Rlti,j . Therefore, this process model comprises 139 equations per ethylene plant per time period (considering eight furnaces in each plant) while the linear model only comprises six equations per plant per time period
as follows:
ffoti,h,j = f (Ffiti , Rdti , Pfith,i , Convti,h ),
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
i = EPI, EPII; ∀j, h; t = 1, . . . , H
(36)
frti = f (ffoti,h,j , PCGC , Rlti , Su,i ), i = EPI, EPII; t = 1, . . . , H
(37)
t 0 = Cpoly + Cpoly
t
Fpinpoly −
t
n=1
Fponpoly ,
n=1
(41)
When inventory levels do not meet given storage targets (IOpoly ) an economic penalty (PENPOLYpoly ) is added to the objective function (Jackson et al., 2003) t PENPOLYtpoly ≥ IOpoly − Cpoly ,
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
(42)
Pfith,i = f (Ffiti , Rdti , Pfoth,i , Convti,h ), i = EPI, EPII; t = 1, . . . , H
(38)
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
t ffttu,i,j = f (PCGCi , Rlti , Su,i ),
i = EPI, EPII; ∀j, u; t = 1, . . . , H
(39)
(43)
Sales cannot exceed the daily forecast demand (demtpoly ) and the difference between the demand and the sales (deltatpoly ), Eq. (45) is penalized in the objective function demtpoly ≥ Fpotpoly ,
t Fbtu,i,j = f (PCGCi , Rlti , Su,i ),
i = EPI, EPII; ∀j, u; t = 1, . . . , H
t − IOpoly , PENPOLYtpoly ≥ Cpoly
(40)
Ethylene plant II is more efficient than ethylene plant I. Therefore ethylene plant II correlations have been adjusted to obtain higher production levels. 3.1.3. Polyethylene plants Ethylene is raw material for polyethylene plants: linear low density polyethylene (LLDPE), low density polyethylene (LDPE) and high density polyethylene (HDPE). Linear inlet–outlet models have been included, based on the information given in Table 2. Polyethylenes bulk storage is modeled
poly = HDPE, LDPE, LLDPE; t = 1, . . . , H
(44)
deltatpoly = demtpoly − Fpotpoly , poly = HDPE, LDPE, LLDPE;
t = 1, . . . , H
(45)
3.1.4. VCM and PVC plants These plants also use ethylene as raw material. Linear models based on the information given in Table 3 have been considered. Product daily demands have been considered.
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316 Table 3 NaOH, VCM, EDC and PVC yields
Table 5 Economical data
Product
t of product/t of C2 H4
NaOH VCM EDC PVC
0.0573 0.0093 0.0087 0.0653
Table 4 Ammonia and urea yields Product
kmol of product/kmol of methane
NH3 (NH2 )2 CO
0.2329 0.9506
3.1.5. Ammonia and urea plants Residual gas (99% methane) from natural gas plant I is used as raw material for the ammonia plant. Methane and air in stoichiometric ratio are used to obtain nitrogen and hydrogen for ammonia synthesis. Carbon dioxide is also obtained and is later used in the urea production. Ammonia is used in urea production and it is also sold as final product. Linear correlations have been derived for both plants production, as function of natural gas consumption based on data from Table 4. Discontinuous delivery of urea and ammonia by ship has been modeled with equations analogous to Eqs. (20)–(28). They are delivered by train and truck according to a forecast daily demand for each product Vvt = Vvo +
t n=1
Fvvnin,v −
t
1313
Fvvnout,v −
n=1
t
Fvtnout,v ,
Natural gas (US$/kmol) Residual gas (US$/kmol) C2 H4 (US$/kmol) C2 H6 (US$/kmol) C3 H8 (US$/kmol) C4 H10 (US$/kmol) C5 H12 (US$/kmol) NH3 (US$/kmol) (NH2 )2 CO (US$/kmol) LDPE (US$/t) HDPE (US$/t) LLDPE (US$/t) NaOH (US$/t) VCM (US$/t) EDC (US$/t) PVC (US$/t) Inventory cost (US$/kmol) Electricity (kWh)
Ethylene plants: Linear correlations have been derived from rigorous plant simulations for required fuel gas for furnaces (both plants) and boilers (only ethylene plant I) as function of ethane flowrate. Remaining plants: Utilities consumption for polyethylenes, VCM, PVC, ammonia and urea plants have been taken from the literature. They are shown in Table 5. Tank inventory costs (ethane, propane, butane, gasoline, ammonia and ethylene) are calculated according to the trapezoidal area, where Civv is the cost coefficient. Urea storage is also handled this way
n=1
v = vNH3 , vurea; t = 1, . . . , H
(46)
There are capacity lower and upper limits for ammonia tank and urea storage room. • Objective function: The objective function is the maximization of net profit along a 20 days horizon, defined as the difference of the incomes from product sales and the costs (operating costs, raw materials costs and tank inventory costs) and the penalties for not satisfying forecast demands, ship capacities and inventory targets profit = sales revenue − raw materials costs −operating costs − Cinv − CIO − Cdem
CInv =
• Operating costs: Natural gas plants: the main operating costs are associated to compression loads. Compressors are driven by gas turbines. In natural gas plant I, compression (33–60 bar) and recompression (18–35 bar) are required. In natural gas plant II, not far from gas wells, no inlet gas compression is required. Linear models for the gas consumption in the compressor turbines have been derived from rigorous simulations as function of raw material flow.
Civv
v
V t + V t−1 k k 2 t
(48)
• Penalty terms: Penalties on deviations from inventory targets for polyethylenes and for not meeting forecast demand (or ships capacity) for the remaining products have been included t θpoly (49) PenPolypoly CIO = t
poly
Cdem = (47)
1.8123 1.63863 11.73 4.1253 7.7625 10.23236 12.70224 2.805 9.9 640 720 600 230 400 240 1080 0.05 0.09454
j
ξj
t
deltatj
(50)
4. Discussion of results Two multiperiod mixed integer nonlinear programming mathematical models have been formulated for the production chain optimization of a petrochemical complex: MINLP1 includes Eqs. (1)–(3), (8)–(35) and (41)–(50) and MINLP2 comprises Eqs. (1)–(3), (8)–(34) and (36)–(50). The main difference between these models is the ethylene plant
1314
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
Fig. 4. Optimal ethylene production profile in ethylene plant II with alternative MINLP models.
model which is simplified and linear in MINLP1 and semirigorous nonlinear in MINLP2. An approximate MILP model has been also formulated replacing bilinear constraints in the MINLPs by valid linearisations (McCormick, 1976), Eqs. (4)–(7), to provide a valid initial point for the MINLPs. Main variables profiles of the MINLPs have been compared, as it is shown in Figs. 4 and 5. Although they show a similar trend, the inclusion of semi-rigorous models for both ethylene plants (MINLP2) renders more uniform production profiles than the obtained with the linear models (MINLP1). In ethylene plant II, total ethylene production linear model has an average error of 6.6% (Fig. 4). Average errors of 8% can be seen in Fig. 5 for the ammonia plant. Linear models for plants have been derived for operating conditions not far from actual ones. Consequently, errors introduced when using linear models are mainly due to process flowrates that are out of range. This could have been improved by replacing linear models by piecewise linear functions. However, as nonlinear models have shown a good performance, further analysis has been made with these ones. Regarding natural gas plant II, propane, butanes and gasoline tanks profiles and the arrival days for the corresponding ships are shown in Fig. 6. After product delivery by ship, the plant production decreases. Natural gas plant I operates at its highest production level throughout the entire time horizon (demethanizer top pressure at its lower bound, 17 bar, and highest ethane recovery). Both ethylene plants require maximum allowed ethane amount from this plant, mainly due to ethane lower price.
Fig. 6. Propane, butanes and gasoline tanks profiles in natural gas plant II.
Figs. 7 and 8 show total furnaces charge and conversion profiles for three of the eight furnaces in ethylene plant II. A slight increase in conversion during the last 4 days is associated to the decrease in total furnaces charge, as it is always seen in actual plant behaviour. As an example of polyethylene plants, Fig. 9 shows optimal LLDPE production profile and inventory level, which most of the time meets the inventory target (360 t, horizontal line). Urea production profile, inventory level and amount daily delivered by trains and trucks are plotted in Fig. 10. Maximum stored amount is 45,150 t and the ship arrival day is 18. Urea amount daily delivered by trains and trucks always meets constant daily demands (900 t). Finally, a different scenario has been analysed with MINLP2 model (Case 2). A 20% increase in HDPE has been specified. Fig. 11 shows that a slight increase in ethylene production throughout the time horizon is enough to fulfil
Fig. 7. Furnaces total feed and ethane recycle profiles in ethylene plant II.
Fig. 5. Optimal urea production profile with alternative MINLP models.
Fig. 8. Furnaces conversion profiles in ethylene plant II.
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
1315
Table 6 Computational details
Constraints Continuous variables Binary variables CPU (s) Solver Objective function (US$)
MINLP1
MINLP2
6932 6123 200 724 DICOPT++ (CONOPT2/OSL) 28624443
16373 18408 200 5002 DICOPT++ (CONOPT2/OSL) 30014774
the increased raw material demand. HDPE forecasted daily demand and HDPE production profile are plotted in Fig. 12, showing that the new demand is satisfied most of the days with an overall increase in HDPE production of 14%.
Fig. 12. Comparison between HDPE production in original Case 1 and Case 2 (20% increase in HDPE demand).
5. Computational details Fig. 9. LLDPE production and inventory profiles. Inventory target at 360 t (horizontal line).
Multiperiod MILP and MINLPs models are coded in GAMS 2.25 (Brooke, Kendrick, & Meeraus, 1992) modelling environment and solved with OSL and DICOPT++ (CONOPT2 and OSL), respectively. MINLP models statistics are shown in Table 6 and both models have required three major iterations to achieve convergence. The approximate MILP model, solved to obtain a valid initial point for MINLP, comprises 13,376 constraints with the same number of continuous and binary variables as MINLP1 model. The problems have been solved in a 1 GHz Pentium III PC, with 256 MB of RAM.
6. Conclusions Fig. 10. Urea production profile, inventory level and amount delivered by trains and trucks.
Fig. 11. Comparison between ethylene production in ethylene plant II in original MINLP2 problem (Case 1) and with a 20% increase in HDPE daily demand (Case 2).
Current petrochemical sites tend to be single integrated facilities to optimize total processing economics for converting natural gas to final industrial products. A model that represents the entire supply chain can provide a better insight on global profit. Two production chain optimization models (MILP1 and MINLP2) of increasing degree of complexity have been formulated. Their solutions provide coordination of responses to demands, production planning, production distribution and inventory level management. The most detailed MINLP solution also determines main plants operating conditions and optimization variables such as conversion in each cracking furnace. Two scenarios corresponding to different product demands have been analysed. An optimization approach in which all the processes are highly integrated has proved to be a useful tool to detect
1316
E.P. Schulz et al. / Computers and Chemical Engineering 29 (2005) 1305–1316
operational and economical improvements. The inclusion of more detailed process models, the exploration of solution methods that can cope with the resulting large-scale MINLP and the inclusion of uncertainty in several process parameters are part of current work.
References Al-Sharrah, G. K., Alatiqi, I., Elkamel, A., & Alper, E. (2001). Planning an integrated petrochemical industry with an environmental objective. Industrial and Engineering Chemical Research, 40, 2103–2111. Brooke, A., Kendrick, D., & Meeraus, A. (1992). GAMS—A user’s guide, release 2.25. Danvers, MA: Boyd and Frase Publ. Co. Bok, J., Grossman, I. E., & Parks, S. (2000). Supply chain optimisation in continuous flexible process networks. Industrial and Engineering Chemical Research, 39, 1279–1290. Diaz, M. S., & Bandoni, A. (1996). A mixed integer optimization strategy for a real large scale chemical plant. Computers and Chemical Engineering, 20(5), 531–545. Diaz, M. S., Serrani, A., de Beistegui, R., & Brignole, E. (1995). An MINLP strategy for the debottlenecking problem in an ethane extraction plant. Computers and Chemical Engineering, 19, 175– 178. Furlonge, H. I., & Young-Hoon, A. (2000). Towards the development of an optimal long-term structure and policy for the development of Trinidad and Tobagoˇıs petrochemical industry. Part I. The methanebased complex. West Indian Journal of Engineering, 22(2). Jackson, J. R., & Grossmann, I. E. (2003). Temporal decomposition scheme for nonlinear multisite production planning and distribution models. Industrial and Engineering Chemical Research, 42, 3045–3055. Jackson, J., Hofmann, J., Wassick, J., & Grossmann, I. E. (2003). A nonlinear multiperiod process optimization model for production planning in multi-plant facilities. In I. E. Grossmann & C. M. McDonald (Eds.).
Fourth International Conference on Foundations of computer Aided Process Operations (FOCAPO), 281–284. Li, W., & Hiu, C.-W. (2002). Scheduling crude oil unloading. storage and processing. Industrial and Engineering Chemical Research, 41, 6723–6734. Lee, H., Pinto, J. M., Grossmann, I. E., & Park, S. (1996). Mixed-integer linear programming model for refinery short-term scheduling of crude oil unloading with inventory management. Industrial and Engineering Chemical Research, 35, 1630–1641. McCormick, G. P. (1976). Computability of global solutions to factorable nonconvex programs. Part I. Convex underestimating problems. Mathematical Programming, 10, 146–175. Neiro, S. M. S., & Pinto, J. M. (2003). Supply chain optimization of petroleum refinery complexes. In I. E. Grossmann & C. M. McDonald (Eds.), Foundations of computer-aided process operations (pp. 59–72). Florida: CACHE AICHE INFORMS. Quesada, I., & Grossmann, I. E. (1995). Global optimization of bilinear process networks with multicomponent flows. Computers and Chemical Engineering, 19(12), 1219–1242. Reddy, P. C. P., Karimi, I. A., & Srinivasan. (2004). A new continuoustime formulation for scheduling crude oil operations. Chemical Engineering Science, 59, 1325–1341. Schulz, E., Diaz, M. S., & Bandoni, A. (2000). Interaction between process plant operation and cracking furnaces maintenance policy in an ethylene plant. Computer Aided and Process Engineering, 8, 487–492. Shah, N. (2004). Process industry supply chains: Advances and challenges. In A. Barbosa-P´ovoa & H. Matos (Eds.). Computer Aided and Process Engineering, 18, 123–138. Turkay, M., & Grossmann, I. E. (1996). Disjunctive programming techniques for the optimization for process systems with discontinuous investment costs—Multiple size regions. Industrial and Engineering Chemical Research, 35, 2611–2623. Turkay, M., Asakura, T., Fujita, K., Hui, C. W., & Natori Y. (1998). Total site optimization of a petrochemical complex. In J. F. Pekny, & G. Blau (Eds.), Foundations of computer-aided process operations, AIChE symposium series no. 320: vol. 94 (pp. 185–189) Production editor, CACHE Publications, B. Carnahan.