Comparison of Methods for Computing Crude Distillation Product

Oct 20, 2015 - fixed-cut swing method. 1. INTRODUCTION. Crude distillation units (CDUs) are the first processing unit in oil refineries. CDUs separate...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/IECR

Comparison of Methods for Computing Crude Distillation Product Properties in Production Planning and Scheduling Gang Fu and Vladimir Mahalec* Department of Chemical Engineering, McMaster University, 1280 Main Street West, Hamilton, Ontario L8S 4L8, Canada ABSTRACT: Production planning and scheduling optimize feedstocks usage in refineries by using crude simplified distillation models which relate feed true boiling point (TBP) curves to product TBP curves and from these calculate product yields and properties. We compare product yield and properties calculated by the swing-cut methods (fixed-cut, weight/volume transfer, and light/heavy) with results from the pseudocuts TBP-based method. The latter uses product yields and TBP curves computed via the hybrid distillation model. Swing-cut methods use assumptions which lead to lower accuracy of predicting yields vs hybrid model yields, resulting in errors in computed bulk product properties. If swing-cut methods yields are replaced by the correct yields from hybrid model, all four methods have approximately the same accuracy. Therefore, for a mixture of crudes, by using accurate yields from the hybrid model, product bulk properties can be computed by blending rules as simple as those used by the fixed-cut swing method.

1. INTRODUCTION Crude distillation units (CDUs) are the first processing unit in oil refineries. CDUs separate crude oil into intermediate fractions according to boiling ranges. These fractions are processed by the downstream units to upgrade their qualities and are finally blended into the final products. CDU is a complex distillation process, producing several fractions and having many degrees of freedom which can be used to optimize the operation. It is an energy intensive process; on average a CDU utilizes around 20% of a refinery’s total energy consumption. Optimization of CDU operation is a very important component of production planning, of production scheduling, and of real-time operations. Numerous versions of simplified CDU models have been proposed with the aim of having a model simple enough to be used for planning and scheduling and also able to predict accurately the properties of the products. In this work, we describe product property calculation methods which rely on accurate product true boiling point distillation curve predictions, such as those from a hybrid model of a crude unit. Figure 1 shows an example of a CDU (AspenTech, 2006)1 consisting of a preflash tower (which removes light components from the feed), atmospheric distillation (which operates at atmospheric pressure and separates the bulk of the crude into several products), and vacuum distillation (which operates under vacuum to separate the heavy end of the crude into several products). This CDU is used as a case study for comparing different property prediction methods. If a CDU can separate sharply each of the crude fractions, then each product from the CDU will have the yield corresponding to the width of the cut and its TBP curve will overlap its section of the crude TBP curve. In reality, product distillation curves differ significantly from their respective sections of the crude TBP curve. The back-end a product TBP curve is above the crude TBP curve and the front end of the product TBP is lower than the crude TBP curve. The middle point of a product TBP curve, as a rule, is not located on the crude TBP curve. Figure 2 shows crude TBP and product © XXXX American Chemical Society

distillation curves for the CDU in Figure 1. The entire TBP curve is divided into nonoverlapping sections (“cuts”). Other crude properties, e.g. % sulfur or gravity or viscosity, also vary from one temperature range to another temperature range (from one cut to another), as shown in Figure 3 and Figure 4. One should note that the “S” like curves of TBP, gravity, and sulfur are of similar shape. That is because all these properties are determined by distribution of hydrocarbons with different boiling points. Intuitively we can conclude that an accurate prediction of specific gravity and sulfur content of products can be computed via the product TBP curve. A complete CDU model has three parts. They are (a) crude oil assay data modeling in a form required by the CDU model, (b) the CDU model, which predicts product yields and TBP curves, and (c) a model for computation of additional product properties (e.g., sulfur, specific gravity, cetane index, etc.). Crude oil assays divide crude intp distinct cuts (e.g., heavy naphtha, kerosene, etc.), assuming perfect separation between adjacent cuts. Each cut is characterized by its initial and final TBP values, its yield, and an average value for each of the properties provided by the assay (e.g., specific gravity, % sulfur, octane number, etc.). These assay data are often used to construct property curves (e.g., % sulfur vs TBP or % sulfur vs LV% vaporized) or are used directly, depending on the needs of the CDU model. With respect to the prediction of product yields and TBP curves, there are two widely used assumptions in many published works on simplified crude distillation modeling or in refinery planning software. These assumptions are the basis for simple models which predict product TBP curves, part (b) of the model. They assume that (i) the back end and front end points of the adjacent products are equidistant from the crude Received: August 7, 2015 Revised: September 29, 2015 Accepted: October 19, 2015

A

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 1. Example crude distillation unit.

Figure 2. Feed and product TBP curves for a crude distillation unit. Figure 4. Sulfur content of crude oil and CDU products.

distillation curve and (ii) the midpoint of a product distillation curve lies on the crude distillation curve. However, inspection of actual CDU product TBP curves (from a rigorous simulation or from a plant) shows that the back end of the lighter cut and the front end of the adjacent heavier cut are not equidistant from the crude TBP curve. In addition, midpoints of product TBP curves do not lie on the crude TBP distillation curve. The hybrid model for accurate prediction of product yields TBP curves is described only briefly in this paper, since it is a subject of a separate paper. Interpreting crude assay model and TBPcut-based property prediction models will be described in detail in this paper. Our goal is to examine different methods for computing product yields and bulk product properties, evaluate their accuracies, and determine which method is most suited (i.e., introduces the least amount of nonlinear constraints) for use in multiperiod planning models. Specifically, we evaluate bulk properties calculation methods based on TBP-cuts, traditional swing-cut, light and heavy swing-cut, and volume (mass) transfer methods. Accurate

Figure 3. Specific gravity of crude oil and CDU products.

B

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Each production mode had a different set of deltas. However, this model still does not deal well with processing different types of crude oil feed, since deltas also depend on the crude properties, not just on the CDU operation. Another possibility is to use delta-based models based on operating conditions. Pinto et al. (2000)3 and Neiro and Pinto (2004)4 proposed the use of nonlinear models to derive delta models for flows and quality values. The second type of delta-based models are swing-cut methods. They require that one estimate for each product the swing-cut size and TBP range. From this information and crude properties in the corresponding TBP range, a method-specific procedure is used to calculate product properties. Zhang et al. (2001)5 used the traditional swing-cut model, which assumed that the size of the swing-cut is fixed and that the properties of the swing part are constant across the cut. Li et al. (2010)6 used the flexible swing-cut method, which calculated the swing-cut based on different operating modes (max naphtha, max diesel, or max gasoline) using a “weight transfer ratio” method. The properties of swing-cut and fractions were determined by regression models based on crude properties. The regression coefficients needed to be updated if the properties of the feed changed. Guerra et al. (2011)7,8 also used a volumetric transfer ratio-based swing-cut method. In order to correct product yields, they introduced a bias term which minimized deviations from a rigorous simulation. Menezes et al. (2013)9 improved the swing-cut method by separating the swing-cut into two parts (light and heavy part). The properties of these two parts were calculated individually using interpolated qualities relative to their light and heavy swing-cut quantities. Several prior works have used the full model approach to development of simplified CDU models. Caballero and Grossmann (1999)10 proposed an aggregated model for the striping and rectifying sections. Gadalla and Smith (2003)11 extended this aggregated model to use a cascade column to represent a complex distillation column. Based on the cascade column structures, the nonlinear model of CDU has been proposed for planning by Alattas et al. (2011).12 They represented a distillation unit as a sequence of flash processes using the fractionation index to determine the separation in each section. All the models above rely on the assumption that the distance from the crude TBP curve to the back end of the light fraction is equal to the distance from the crude TBP curve to the front end of the neighbor heavier cut. Mahalec and Sanchez (2013)13 introduced a hybrid model of an atmospheric distillation tower; the model does not rely on this middle distance assumption. The model was designed for real time applications, so the temperature profile along the tower was estimated from several available tray temperature measurements. Then these data were used to calculate internal reflux ratios at various points along the tower. The model was demonstrated to predict product TBP with less than 1% error when compared to the rigorous tray-to-tray distillation model. The model is not applicable for planning and scheduling, since the tray temperatures are not available for the future operations. Fu and Mahalec (2015) developed a modified hybrid model suitable for planning, for scheduling, or for RTO. The model predicts product TBP curves from the feed properties, product yields, and operating conditions (furnace outlet temperatures, pumparound duties, stripping steam flows). A more detailed review of efforts to develop reduced order crude distillation models was presented by Ochoa-Estopier et

product TBP curves are computed via the hybrid model of the CDU. This model does not rely on the equidistance assumption, and the product TBP curves correspond to those that are observed in actual CDUs. We illustrate how product and crude TBP curves and property distribution curves can be used to compute the bulk properties (e.g., % sulfur) of the product streams. TBP curves provide 9 points (1%, 5%, 10%, 30%, 50%, 90%, 95%, and 99%) for each product. Results computed by the TBP-cut method combined with yields from the hybrid model are compared with those from a rigorous trayto-tray model and with the results from various swing-cut methods. In Section 2 of this paper we present a brief review of the prior work on crude distillation units related to the products properties prediction for applications in planning, scheduling, and RTO. Section 3 describes a method to estimate pseudocomponents and properties for each pseudocomponent starting from crude assay data. This makes it possible to compute properties of crude mixtures as required for the CDU distillation model. Section 4 describes briefly the hybrid CDUs model for TBP prediction. Section 5 describes computation of other stream properties (e.g., specific gravity and sulfur). Comparisons of product properties predictions via TBP computed from the hybrid model with rigorous tray-to-tray model predictions and with predictions from swing-cut methods are given in Section 6. Section 7 presents the conclusions.

2. LITERATURE REVIEW The rigorous distillation tower model provides accurate and robust models for predicting CDU product yields and properties. However, rigorous distillation models have many equations and are highly nonlinear and therefore are not suitable for planning and scheduling purposes. In order to find an optimal solution in reasonable solution times for planning and scheduling models, crude distillation units have traditionally been represented by various forms of linearizations or by simplified nonlinear representations. On the other hand, real time optimization (RTO) applications use a tray-to-tray rigorous model for production operation optimization. This separate modeling path lead to inconsistent decisions between planning, scheduling, and RTO. Increased accuracy of CDU models for planning and scheduling is highly desirable, since it can enable better quality of decisions with respect to crude purchase and crude processing. If such models are also accurate enough for RTO, that would help eliminate the discrepancies between the three decision processes. One approach to simplified CDU modeling is to use incremental changes (“delta-based” models) to product TBP curves deviations from some base case operation. Another approach (we shall call it “full model”) is to compute product TBP curves via a model which uses operating conditions, crude TBP data, and sometimes approximations of the CDU structure. There are two types of delta-based models. One of them is to modify the back end and the front end of the product by adding or subtracting some delta differences. Since the CDU unit can operate under different operating modes, the deviations from the base case are not constant but depend on the operating mode. Product yields and the cutpoints are calculated by using the middle point of the modified front and back ends of TBP curves of adjacent products. Brooks et al. (1999)2 improved this kind of model by defining different production modes. C

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research al. (2014).14 They also developed a very accurate neural network-based model of a crude distillation unit and compared its results to a rigorous simulation. Once product yields and TBP curves are estimated, one can proceed to compute other product properties. Recent work by Guerra et al. (2011)7,8 fitted crude properties (e.g., sulfur, cetane index, etc.) by a fourth order or lower polynomial. The bulk properties of the products were then computed by volumetric or mass averaging contributions to the product property from each of the crudes. Kelly et al. (2014)15 used monotonic interpolation methods for interpolation to avoid Runge’s phenomenon. They compared two different monotonic interpolations (PCHIP and Kruger) with simple linear interpolation in TBP curve prediction. Four case studies were given applying monotonic interpolation for blending property prediction and cutpoint optimization. Menezes e al. (2013)9 introduced light/heavy swing-cut and combined it with TBP microcuts to compute bulk product properties; the example presented in their paper showed that the predictions from both methods were similar. The TBP microcuts method is very similar to the methods that have been used in rigorous simulators since the late 1980s.

Table 2. Crude Oil 1 Assay Data TBP

100 to 800 F 800 to 1200 1200 to 1400 1400 to 1640

25 50 100 120

6.8 10 30 50 62 70

130 180 418 650 800 903

76 90

1000 1255

TBP

Table 1. Pseudo-component Definition Rule in This Paper Increment

LV%

name Methane Ethane Propane Isobutane N-Butane 2-Methylbutane N-Pentane Water

API curve LV (%)

LV%

sulfur curve LV (%)

wt %

0.1 0.15 0.9 0.4 1.6 1.2

5 10 15 20 30 40

90 68 59.7 52 42 35

2 5 10 20 30 40

0 0.01 0.013 0.05 1.15 1.62

1.7 0

45 50 60 70 80 bulk

32 28.5 23 18 13.5 31.4

45 50 60 70 80 bulk

1.9 2.15 2.54 3 3.7 2.3

Table 3. Crude Oil 2 Assay Data

3. COMPUTATION OF CRUDE OIL PROPERTIES CURVES AND CRUDE MIXING FROM CRUDE ASSAY DATA Crude oil typically consists of a large number of compounds, and its chemical composition is not known. The petroleum refining community has adopted crude characterization in a form of crude assays. This paper uses the used pseudocomponent method, in which the crude oil is cut into pseudocomponents based on boiling ranges. The method is used by the rigorous simulation software (HYSIS, AspenPlus, Pro/II). We have adopted the pseudocomponent definition based on the rules shown in Table 1.

Boiling-point range

light end

temp (F)

LV%

temp (F)

6.5 10 20 30 40 50

120 200 300 400 470 550

60 70 80 90 95 98 100

650 750 850 1100 1300 1475 1670

light end

Methane Ethane Propane Isobutane N-Butane 2-Methylbutane N-Pentane water

API

sulfur

LV (%)

LV (%)

API

LV (%)

wt %

0.2 0.5 0.5 1 1 0.5

2 5 10 20 30 40

150 95 65 45 40 38

2 5 10 20 30 40

0 0.01 0.015 0.056 1.3 1.7

2.5 0.1

50 60 70 80 90 95 98 bulk

33 30 25 20 15 10 5 34.8

45 50 60 70 80

2 2.3 2.7 3.2 3.8

bulk

2.5

squares optimization as shown by eq 2. The results of extrapolation for the Crude oil 1 TBP curve are shown in Figure 5 and Table 4 compared with AspenPlus results. We can see that the four-parameters beta function can give results

In order to illustrate the method, two crude oils are used in this paper. The assay data of Crude oil 1 and Crude oil 2 are shown in Table 2 and Table 3, respectively. Generation of property curves for a crude oil. The procedure to generate a crude oil TBP curve and associated properties curves from the assay data is presented here by using Crude oil 1 as an example: 1. Extrapolate the TBP curve based on the crude assay data. It is very common that crude assays do not address the high boiling point range. So we need to extrapolate the incomplete TBP curve to cover the boiling point range. Sanchez et al. (2007)16 reviewed several different probability distribution functions to the fit distillation curves of petroleum fractions. They concluded that the cumulative beta function with 4 parameters ensures a good extrapolation. Therefore, in this work the beta cumulative density function and least-squares are used to perform extrapolation of crude oil TBP. The cumulative beta density function is given by eq 1. The parameters used in this equation are calculated by least−

Figure 5. TBP curve beta function (four parameters) extrapolation results for crude oil 1. D

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Table 4. Beta Extrapolation Result Compared with AspenPlus LV (%)

AspenPlus

beta extrapolate

0 1 5 10 30 50 70 90 95 99 100

−75.63 −13.01 94.22 180.13 418.26 650.04 903.60 1256.54 1410.12 1548.51 1561.44

1.75 29.75 104.81 177.27 418.78 650.01 907.00 1256.42 1393.42 1578.60 1686.45

which are differ from those from AspenPlus only at the very ends of the boiling range. f (x , α , β , A , B ) x≤B Γ(α + β) ⎛ x − A ⎞α − 1 1 ⎜ ⎟ = A B − A Γ(α)(Γ(β)) ⎝ B − A ⎠ ⎛ B − x ⎞ β−1 ⎟ ×⎜ dx ⎝B − A⎠

Figure 7. Specific gravity for each pseudocomponent for crude oil 1.



x volume ; α , β , A , B , are four parameters

(1)

n

min ∑ (LV %pred , k − LV %assay , k )2 k=0

k = 0...n known crude assay data

(2)

2. Use linear interpolation to compute the volumetric percentage of each pseudocomponent. Once we define the boiling point range for each pseudocomponent, the volumetric percentage of each pseudocomponent can be calculated (Figure 6). 3. Compute other crude properties for each pseudocomponent by linear interpolation of the corresponding crude oil property curve. The middle of the pseudocomponent TBP range is used as the boiling point of this pseudocomponent. Resulting property curves for specific gravity and sulfur are shown in Figure 7 and Figure 8, respectively.

Figure 8. Sulfur content for each pseudocomponent for crude oil 1.

After computing the property distribution curve vs TBP, the bulk properties of the crude oil can be estimated as a volumetric or a weight-based average along the curve. For volumetricbased properties (such as specific gravity), eq 3 is used. For mass-based properties (such as sulfur), eq 4 is used. The bulk properties prediction results for Crude oil 1 and Crude oil 2 are shown in Table 5. ps = 0

BulkPropvol =



volps*Propps

n

ps = pesudocomponent ps = 0

BulkPropmass =

∑n

(3)

volps*SGps*propps

ps = 0

∑n

ps = pesudocomponent

volps*SGps (4)

Crude oil mixing. The procedures for crude oil mix modeling 1. Calculate volumetric amount of each light component and pseudocomponent in the mixture. Volumetric amounts are

Figure 6. Volumetric percentage for each pseudocomponent for crude oil 1. E

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

The graphic illustration of these two computational steps is given in Figure 9. Taking kerosene as an example, the middle

Table 5. Crude Properties Computation Compared with AspenPlus and Assay Data crude oil 1

crude oil 2

mixture (crude oil 1:2 = 0.2:0.8)

SG (g/cm3)

S (wt %)

SG (g/cm3)

S (%)

SG (g/cm3)

S (wt %)

0.864 0.869 0.869

2.328 2.264 2.300

0.847 0.851 0.851

2.545 2.482 2.500

0.851 0.854 NA

2.503 2.438 NA

this work AspenPlus assay data

added for each component proportionally to the fraction of the crude in the mixture, eq 5. Following that, the TBP curve of the crude mixture is calculated by accumulating the volumes of each light component and of each pseudocomponent. crd = 0

volps =



volps , crd*ratiocrd

n

ps = pesudocomponent ; crd = 1, 2 (crude oil 1, crude oil 2)

Figure 9. Two steps TBP prediction via the hybrid model.

(5)

2. Calculate volume-based properties (such as specific gravity) for each (pseudo) component by using eq 6.

section of the distillation curve is defined by the atmospheric tower feed distillation curve and the yield information on the fractions. Then a straight line is fitted through 30, 50, and 70 vol %. Partial least-squares (PLS) models which relate the operating variables to product distillation curves were used to predict the vertical deviation between front and back sections of the curve and this straight line. The majority of the equations in the model are linear. An exception is the equation which computes the reflux from the [reflux/(reflux + distillate)] variable in the atmospheric tower model. Even though the model is small (about 160 equations to describe preflash, atmospheric, and vacuum towers), the model predicts very accurately the product TBP curves with 0.5% RMSE over a wide range of conditions. Sample TBP predictions for two crude mixtures and their comparison with the results from a rigorous model in AspenPlus are shown in Table 6. In order to develop a hybrid model for use in planning and scheduling, an engineer would have to collect the plant data or generate such data from a rigorous simulation (crude mix being processed, operating conditions, product yields, product distillation curves, other properties) and use partial leastsquares (PLS) to construct the model for prediction of product TBP curves. Once these curves are known, methods described in this work can be used to compute other product properties.

crd = 0

propps =



ratiovol , crd *propps , crd

n

ps = pesudocomponents ; crd = 1, 2 (crude oil 1, crude oil 2)

(6)

3. Calculate weight-based properties (such as sulfur) by using eq 7. crd = 0

propps =

∑n

ratiovol , crd *SGps , crd*propps , crd crd = 0

∑n

ratiovol , crd *SGps , crd

ps = pesudocomponents ; crd = 1, 2 (crude oil 1, crude oil 2)

(7)

Once the property curve is obtained, the bulk properties for the crude mixture can be calculated using eq 3 and eq 4.

4. HYBRID MODEL OF CDU FOR PREDICTION OF PRODUCT TBP CURVES Here we provide a brief overview of the hybrid model of a CDU. A detailed description of this hybrid model for prediction of product TBP curves will be published elsewhere. The main idea behind the hybrid model of the CDU is to predict first the straight line which passes through the middle bulk of the product TBP curve and then compute deviations from the line at the front and at the back of the product TBP range (Mahalec and Sanchez, 2012).13 The line though the middle section is based on TBP points corresponding to 30 vol %, 50 vol %, and 70 vol % of a product. Extensive rigorous simulations have shown that the middle section is related only to the product yields (cuts) and the corresponding points on the feed TBP distillation curve. The line through the middle of the product TBP curve is determined by the design of the distillation tower. The operating conditions only have an effect on the front and back sections of the products TBP curve. The second step is to predict the front and back section (TBP1%, 5%, 10%, 90%, 95%, and 99%) vertical deviation from the middle straight line.

5. TBP-CUT-BASED PRODUCT PROPERTIES PREDICTION METHOD Product property curves are distributed along LV % of the product according to the property value that each pseudocomponent has. For monotonically increasing properties (specific gravity, viscosity) the product property curve has an “S” shape similar to the product TBP curve. This is because distillation is a physical process, which only changes the distribution of pseudocomponents in products. The properties for each pseudocomponent do not change during processing in CDUs. Once the product TBP curve is known, it can be used to calculate other properties (such as specific gravity or sulfur) at any point along the vol % axis by linear interpolation. The property calculation model can be one of three types: F

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Table 6. Hybrid Model Product TBP Predictions from Hybrid Model and from AspenPlus Crude oil 1:Crude oil 2 = 0.3:0.7 Mixture of 2 crude oils PF

Aspen Plus

Hybrid model

Err (%)

Aspen Plus

Hybrid model

Err (%)

F F F F F F F F F F F F F F

467.1 380.0 311.0 1393.0 680.0 360.0 520.0 640.0 351.7 461.6 538.7 882.1 1183.4 775.5

463.6 376.7 312.2 1391.9 676.4 360.7 519.2 638.1 350.0 458.9 536.7 875.3 1183.4 773.3

−0.77 −0.86 0.38 −0.08 −0.53 0.18 −0.15 −0.29 −0.49 −0.58 −0.37 −0.77 0.00 −0.28

472.0 380.0 328.2 1428.0 680.0 360.0 520.0 640.0 349.7 463.2 540.9 912.5 1187.5 803.4

472.6 377.0 327.7 1429.5 674.8 359.7 519.4 637.6 346.6 460.0 537.2 912.5 1190.1 802.7

0.12 −0.79 −0.14 0.10 −0.76 −0.07 −0.11 −0.38 −0.88 −0.69 −0.68 0.00 0.21 −0.09

PF-COT Naphtha TBP95 AP-Feed TBP05 AP-Feed TBP95 AP-COT AP-HNaphtha TBP95 AP-Kero TBP95 AP-Diesel TBP95 AP-Kero TBP05 AP-Diesel TBP05 AP-AGO TBP05 LVGO 95 HVGO 95 HVGO 05

AP

VP

8

volumetric-based (like specific gravity), mass-based (like sulfur content), and blending-index-based (like pour point index). In this paper, computation of volumetric-based and mass-based properties (specific gravity and sulfur content) will be discussed. The following steps illustrate computation of product properties: 1. For each point of this product TBP curve, calculate the volumetric-based properties curve for the product (specific gravity and sulfur) via linear interpolation of the points on the corresponding property curve of the crude feed (crude mixture in our example) via eq 8 and eq 9:

Sulfuri = [ ∑ (SGi , j + 1 + SGi , j)*(LVj + 1 − LVj)/2 j=1

*(sulfuri , j + 1 + sulfuri , j)/2]/SGi*98 i = 0, 1, 2, 3 (nph , kero , dsl , AGO) j = 1, ..., 9 (1%, 5%, 10%, 30%, 50%, 70%, 90%, 95%, 99%)

i = 0, 1, 2, 3 (nph , kero , dsl , AGO) J = 1, ..., 9 (1%, 5%, 10%, 30%, 50%, 70%, 90%, 95%, 99%)

(8)

Sulfur i , j = linear interpolation (TBPi , j ,Sulfurmixcrude) i = 0, 1, 2, 3 (nph , kero , dsl , AGO) j = 1, ..., 9 (1%, 5%, 10%, 30%, 50%, 70%, 90%, 95%, 99%)

(9)

2. Calculate the bulk volumetric-based property of the product i from the cumulative property curve (e.g., specific gravity) as described by eq 10:

Table 7. Specifications for Four Different Modes of Operation

8

∑ j = 1 ((SGi , j + 1 + SGi , j)*(LVj + 1 − LVj))/2 98

standard

i = 0, 1, 2, 3 (nph , kero , dsl , AGO)

max naphtha

j = 1, ..., 9 (1%, 5%, 10%, 30%, 50%, 70%, 90%, 95%, 99%)

(11)

6. CASE STUDIES The crude unit model described in the AspenPlus (2006)1 manual on petroleum process modeling is used in these case studies. The TBP-based property prediction method, which uses the hybrid model results, is compared with the fixed swingcut method, the weight transfer ratio method, light/heavy swing-cut, and the rigors results generated by AspenPlus. Three sets of test have been carried out for four different modes of operation (standard operation, maximum naphtha, maximum kerosene, and maximum diesel). The results are presented in the English system of units, since they are prevalent in refinery operations in North America. Test #1: Verification of hybrid model prediction accuracy for TBP curves. TBP 5% point and TBP 95% points of each products were specified in an AspenPlus tray-totray model to obtain the volumetric flow rate of each product. Details of the specifications are shown in Table 7. This is different from Li et al. (2005),6 who specified only the TBP 99% for each product in AspenPlus. The initial boiling points

SGi , j = linear interpolation (TBPi , j , SGmixcrude)

SGi =

Crude oi 1:Crude oil 2 = 0.7:0.3

unit

max kerosene

(10)

3. Calculate the bulk mass-based property of the product i from the property curve (e.g., sulfur) as given by eq 11:

max diesel

G

TBP05 TBP95 TBP05 TBP95 TBP05 TBP95 TBP05 TBP95

unit

nph

kero

dsl

AGO

F F F F F F F F

NA 360 NA 380 NA 340 NA 360

350 520 360 520 330 540 350 520

460 640 460 640 470 640 460 660

540 NA 540 NA 540 NA 560 NA

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research (IBPs), final boiling points (FBPs), and flow rates were calculated by AspenPlus (Table 8 and Table 9, respectively).

products. For instance, under standard operation, the kerosene 99% TBP point is 550.0, while the diesel 1% TBP is 406.0. The middle point between them is 478.0, as shown in Table 11. The

Table 8. IBPs and FBPs for Each Fraction standard max naphtha max kero max diesel

TBP01 TBP99 TBP01 TBP99 TBP01 TBP99 TBP01 TBP99

unit

nph

kero

dsl

AGO

F F F F F F F F

43.1 380.3 60.6 400.9 28.4 361.1 43.2 380.2

309.6 550.0 320.7 550.0 280.6 564.4 309.7 550.3

406.0 664.6 406.3 664.6 411.5 662.9 403.0 685.1

448.3

Table 11. Middle TBP Point between Adjacent Products

standard

448.3 449.6

max naphtha

467.2 max kero

Table 9. Products Flow Rate for Four Operating Modes standard max naphtha max kero max diesel max min standard max naphtha max kero max diesel max min

unit

nph

kero

dsl

AGO

bbl/day bbl/day bbl/day bbl/day bbl/day bbl/day LV% LV% LV% LV% LV% LV%

5618 6704 4318 5620 6704 4318 7.1 8.4 5.4 7.1 8.4 5.4

16934 15729 20784 16813 20784 15729 21.3 19.8 26.1 21.1 26.1 19.8

10707 10849 7605 13109 13109 7605 13.5 13.6 9.6 16.5 16.5 9.6

12458 12436 13010 10175 13010 10175 15.7 15.6 16.3 12.8 16.3 12.8

max diesel

middle point overlap middle point overlap middle point overlap middle point overlap

unit

cut point (nph-kero)

cut point (kero-dsl)

cut point (dsl-AGO)

F

344.9

478.0

556.4

F F

70.6 360.8

144.0 478.2

216.2 556.5

F F

80.2 320.8

143.7 488.0

216.4 556.2

F F

80.5 345.0

152.9 476.6

213.3 576.1

F

70.5

147.3

217.9

equidistance assumption means that the cutpoint between kerosene and diesel corresponds to the 478 F point on the crude TBP curve, which is 24.6 LV%, as shown in Table 12. Table 12. Cumulative Product Cutpoints Based on Equidistance Assumption vs Actual

standard

These product flow rates, stripping steam flows, pumparound duties, and reflux/(reflux + distillate) were inputs to the hybrid model, which computed its own predictions of the product TBP curves. Table 10 compares the results from the hybrid model with the rigorous AspenPlus model. TBP points computed by the hybrid model are less than 0.5% different from the results of the AspenPlus rigorous CDU model. The only exception is the 1 LV% point for naphtha. Test #2: Evaluation of error introduced by equidistance assumption. Previously published simplified CDU models rely on the assumption that the crude TBP is in the middle of the distance between FBP and IBP of the TBP curves of the adjacent products. Four different modes of operation (standard, max naphtha, max kerosene, and max diesel) are used to test this assumption. Results from the rigorous model (Table 8) were used to calculate the midpoint of adjacent

max naptha

max kerosene

max diesel

Equidistant midpoint Actual (AspenPlus) difference Equidistant midpoint Actual (AspenPlus) difference Equidistant midpoint Actual (AspenPlus) difference Equidistant midpoint Actual (AspenPlus) difference

unit

nphkero

kerodsl

dslAGO

LV% LV% % LV% LV% % LV% LV% % LV% LV% %

7.6 7.1 0.6 8.7 8.4 0.3 5.9 5.4 0.5 7.6 7.1 0.6

24.6 28.3 −3.7 24.7 28.2 −3.5 26.1 31.5 −5.4 24.4 28.2 −3.7

35.3 41.8 −6.5 35.3 41.8 −6.5 35.3 41.1 −5.8 37.8 44.7 −6.9

Rigorous CDU simulation shows that the cutpoint between kerosene and diesel is actually at 28.3 LV%; that is, the error is 3.7% of the crude volume, which correspond to 17% error in the kerosene yield. Complete comparison for all four modes of operation is presented in Table 12. It can be seen that errors in

Table 10. TBP Prediction by Hybrid Model vs AspenPlus

AspenPlus

Hybrid

error

TBP01 TBP05 TBP95 TBP99 TBP01 TBP05 TBP95 TBP99 TBP01 TBP05 TBP95 TBP99

unit

nph

kero

dsl

AGO

residue

F F F F F F F F % % % %

43.1 122.8 360.0 380.3 44.8 123.5 361.1 382.5 −3.71 −0.49 −0.31 −0.58

309.6 350.0 520.0 550.0 311.0 349.3 519.3 549.0 −0.44 0.19 0.14 0.19

406.0 460.0 640.0 664.6 403.3 457.9 639.2 663.7 0.67 0.47 0.12 0.13

448.3 540.0 798.8 850.8 445.9 538.7 800.5 850.9 0.54 0.24 −0.21 −0.02

620.4 692.8 1509.3 1587.1 617.1 690.3 1513.7 1605.9 0.55 0.37 −0.29 −1.17

H

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research product yields are significant for kerosene, diesel, and AGO. This confirms conclusions by Guerra et al. (2011a,b)7,8 who employed yield offset in their model in order to match rigorous simulation. They did not explore how the offset varies with crude type. Test #3: Evaluation of methods to compute product bulk properties. There are four test cases in this test, as shown in Table 13. The typical size of the swing-cut as

Table 15. Specific Gravity Predictions vs AspenPlus

Table 13. Methodology Used in Four Cases case

swing-cut size

1. fixed cut 2. VTR

fixed VTR/WTR

fixed empirical correlation

3. L/H cut

fixed

4. TBP cuts

NA

Mixing of light/heavy parts of each cut TBP-based

properties

source Zhang (2001)5 Li (2005),6 Guerra (2011)7,8 Menezes (2013)9

specific gravity

unit

nph

kero

dsl

AGO

AspenPlus fixed cut VTR L/H cut TBP cuts fixed cut error VTR error L/H error TBP cut error

g/cm3 g/cm3 g/cm3 g/cm3 g/cm3 % % % %

0.763 0.769 0.766 0.769 0.754 0.820 0.443 0.750 −1.154

0.827 0.824 0.824 0.824 0.823 −0.369 −0.394 −0.412 −0.497

0.858 0.848 0.848 0.848 0.858 −1.153 −1.170 −1.147 −0.033

0.889 0.879 0.879 0.879 0.894 −1.154 −1.112 −1.105 0.515

Table 16. Sulfur Predictions vs AspenPlus

This paper

introduced by Zhang (2001)5 is used in this paper. The swingcut between naphtha and kerosene is 5%, that between kerosene and diesel is 7%, and that between diesel and AGO is 9%. In this case study, volumetric transfer ratio (VTR) is used because the volumetric flow rates are used. In case 1 (“fixed cut”), the size and the properties of the swing-cut are fixed. In case 2 (“VTR”), the size of the swing-cut is calculated by VTR. The properties are determined by high order polynomial correlations. In case 3 (“L/H cut”), VTR is used due to the volumetric flow rate. In case 3, the size of the swing-cut is fixed, but the properties of this swing-cut are separated into the light and heavy parts. The properties of the light and heavy parts are calculated as a linear relationship between the properties at their adjacent hypothetical interfaces. For “fixed cut” and “L/H cut” methods properties are computed by blending the swing-cut properties with the bulk of the corresponding product. In case 4 (this work, “TBP” method), the size of the swing-cut does not need to be specified; properties are calculated as described in Section 5. Table 14 shows yield, specific gravity, and % sulfur for the middle section and the cuts corresponding to each product. These data were used to compute the product properties shown in Table 15 and Table 16 (specific gravity and % sulfur, respectively). Also shown in these tables are the results from the rigorous model in AspenPlus and from the TBP method described in this paper. For data used in this work, for monotonically increasing properties (e.g., specific gravity), error relative to AspenPlus is approximately the same for all cases

sulfur content

unit

AspenPlus fixed cut VTR L/H cut TBP cuts fixed cut error VTR error L/H error TBP cut error

wt wt wt wt wt % % % %

% % % % %

nph

kero

dsl

AGO

0.17 0.16 0.05 0.14 0.20 −9.67 −70.58 −19.75 18.82

1.49 1.40 1.40 1.39 1.46 −5.80 −5.74 −6.31 −1.54

2.23 1.99 1.99 2.00 2.21 −10.65 −10.79 −10.38 −0.92

2.83 2.64 2.64 2.65 2.83 −6.56 −6.63 −6.29 0.01

(four methods of computing cuts and product properties); see Table 15. Properties which do not increase monotonically with (LV%, TBP) points (e.g., sulfur) are predicted much more accurately by using the TBP method combined with the accurate TBP predictions from the hybrid model, as illustrated by Table 16. The distribution of specific gravity is shown in Figure 10 and Figure 12, while sulfur distribution is shown in Figure 11 and Figure 13. Test #4: What is the major source of error for swingcut methods? Since there is a fairly large error in computing product yields via swing-cut methods, this test is to evaluate its impact on the accuracy of the product properties predictions. Specifically, we wanted to determine if use of the yields from the hybrid model would improve the accuracy of property computation by the swing-cut methods. In order to compare with the previous results, in the hybrid we set the stripper steam to be the same as in the AspenPlus model. Product yields in the hybrid model were manipulated to minimize deviations in TBP 95% between the hybrid model and AspenPlus. Table 17 shows the TBP 95% values, while Table 18 shows the corresponding product flows in both AspenPlus and the hybrid model. Instead of using product cumulative cuts (cumulative yields) from the

Table 14. Product Properties Computed by Swing-Cut Methods fixed cut naphtha sw1-l sw1-h kerosene sw2-l sw2-h diesel sw3-l sw3-h AGO

VTR

L/H cut

size (vol %)

yield

SG

sulfur

yield

SG

sulfur

size

yield

SG

sulfur

3.0 4.6 0.4 14.0 2.6 4.4 5.0 1.3 7.7 14.4

7.6

0.74 0.79 0.79 0.82 0.84 0.84 0.85 0.86 0.86 0.89

0.03 0.23 l 1.36 1.78 1.78 2.07 2.38 2.38 2.78

7.6

0.77

0.05

7.6

17.0

0.82

1.40

10.7

0.85

1.99

22.1

0.88

2.64

3.0 4.6 0.4 14.0 2.6 4.4 5.0 1.3 7.7 14.4

0.74 0.79 0.80 0.82 0.84 0.84 0.85 0.86 0.87 0.89

0.03 0.21 0.64 1.36 1.68 1.84 2.07 2.23 2.40 2.78

17.0

10.7

22.1

I

17.0

10.7

22.1

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 13. Sulfur content of product streams based on the yields from the hybrid model.

Table 17. Hybrid Model TBP Constraints Matching with AspenPlus

Figure 10. Specific gravity of product streams based on yields using the equidistance assumption.

TBP A+ hybrid error

F F %

nph TBP95

kero TBP95

dsl TBP95

kero TBP05

dsl TBP05

AGO TBP05

360.0 360.5 0.14

520.0 520.5 0.09

640.0 639.7 −0.04

350.0 349.3 −0.21

460.0 459.6 −0.08

540.0 540.5 0.10

Table 18. Hybrid Model Yields vs Aspen Plus when TBP Constraints Are Matched flow rate Aspen plus Hybrid Aspen plus Hybrid

bbl/day bbl/day % %

Nph

Kero

Dsl

AGO

5618 5599 7.06 7.03

16934 17145 21.27 21.54

10707 10561 13.45 13.27

12458 12873 15.65 16.17

the original swing-cut methods (Table 15) are already around 1%. Sulfur results show much more improvement (Table 21); note that the swing-cut methods, if they are provided with the correct yields, are almost as accurate as the TBP-based method.

Figure 11. Sulfur content of product streams based on yields using the equidistance assumption.

7. CONCLUSIONS This work compares four methods for computation of product bulk properties in simplified crude distillation unit models: (i) conventional swing-cut, (ii) light/heavy swing-cut based on TBP-cuts, (iii) volume (weight) transfer, and (iv) accurate product TBP-based method based on TBP-cuts. The first three methods employ variations of rules for transferring selected parts of the crude feed to individual products in order to compute the product yields and bulk properties. The fourth method uses accurate product TBP curves and yields computed via a hybrid model of the CDU; bulk properties are computed by blending TBP-cuts which are present in a given product. All methods were compared to results from rigorous CDU simulation in AspenPlus. Swing-cut and volume transfer methods have difficulties computing correct product yields, which leads to significant errors in predicted product bulk properties. The hybrid model/ TBP-cut method benefits from accurate yield predictions and consequently has higher accuracy of computing bulk properties. When accurate yields from the hybrid model are used in calculations of swing-cuts and volume transfer methods, the errors in results are much smaller and on par with those from hybrid model/TBP-cut methods. An exception is prediction of sulfur level for naphtha product, which is an indication that for

Figure 12. Specific gravity of product streams based on the yields from the hybrid model.

swing-cut methods, we substituted the yields from the hybrid model, and based on these yields, the corresponding (swing) cuts were calculated as shown in Table 19. For instance, for the VTR method, the yield of kerosene changed from 17.0% to 21.5%. This eliminates the error in the cumulative cuts and makes it possible that the property calculation methods [part c of the crude unit model] use the crude property data from the correct region of the crude assay. Specific gravity results are shown in Table 20. There is an improvement in accuracy, even though the prediction errors in J

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Table 19. Product Properties Computed by Swing-Cut Methods When Yields Are Provided from the Hybrid Model fixed cut naphtha sw1-l sw1-h kerosene sw2-l sw2-h diesel sw3-l sw3-h AGO

VTR

yield (vol %)

SG

sulfur

yield (vol %)

SG

sulfur

size

yield (vol %)

SG

sulfur

3.0 4.0 1.0 14.0 6.6 0.4 5.0 7.8 1.2 15.0

7.0

0.74 0.79 0.79 0.82 0.84 0.84 0.85 0.86 0.86 0.89

0.03 0.23 0.23 1.36 1.78 1.78 2.07 2.38 2.38 2.79

7.0

0.76

0.05

7.0

21.6

0.83

1.46

13.2

0.86

2.25

16.2

0.89

2.76

3.0 4.0 1.0 14.0 6.6 0.4 5.0 7.8 1.2 15.0

0.74 0.79 0.80 0.82 0.84 0.85 0.85 0.86 0.87 0.89

0.03 0.18 0.58 1.36 1.77 1.93 2.07 2.36 2.51 2.79

21.6

13.2

16.2

specific gravity

unit

nph

kero

dsl

AGO

AspenPlus fixed cut VTR L/H cut TBP cuts fixed cut error VTR error L/H error TBP cut error

g/cm3 g/cm3 g/cm3 g/cm3 g/cm3 % % % %

0.763 0.765 0.764 0.766 0.754 0.241 0.097 0.400 −1.154

0.827 0.836 0.826 0.826 0.823 1.123 −0.064 −0.082 −0.497

0.858 0.847 0.858 0.858 0.858 −1.323 −0.045 −0.011 −0.033

0.889 0.915 0.885 0.886 0.894 2.931 −0.421 −0.331 0.515

Table 21. Sulfur Prediction When Using Yields from the Hybrid Model sulfur content

unit

AspenPlus fixed cut VTR L/H cut TBP cuts fixed cut error VTR error L/H error TBP cut error

wt wt wt wt wt % % % %

% % % % %

nph

kero

dsl

AGO

0.17 0.15 0.05 0.12 0.20 −13.61 −71.75 −32.53 18.82

1.49 1.44 1.46 1.45 1.47 −3.16 −1.65 −2.38 −1.54

2.23 2.25 2.25 2.24 2.22 0.73 0.93 0.42 −0.92

2.83 2.76 2.76 2.77 2.84 −2.32 −2.51 −2.00 0.01



21.6

13.2

16.2

NOMENCLATURE IBP initial boiling point FBP final boiling point TBP ture boiling point TBPxx true boiling point of xx% percent distillated (liquidvolume-based) SG specific gravity wt % weight percentage LV% liquid volume percentage VTR volumetric transfer ratio WTR weight transfer ratio ps pseudo component crd crude oil npg heavy naphtha kero kerosene dsl diesel AGO atmospheric gasoil Prop property BulkProp bulk property sw-1 swing-cut 1 (swing-cut between heavy naphtha and kerosene) sw-2 swing-cut 2 (swing-cut between kerosene and diesel) sw-3 swing-cut 3 (swing-cut between diesel and ago) sw1-l light part swing-cut 1 sw1-h heavy part swing-cut 1 sw2-l light part swing-cut 2 sw2-h heavy part swing-cut 2 sw3-l light part swing-cut 3 sw3-h heavy part swing-cut 3 L/H improved swing-cut method which separates swingcut into light and heavy i hydrocarbon (feed or products) j point on TBP curve (01,05,10,30,50,70,90,95,99)

Table 20. Specific Gravity Prediction When Using Yields from the Hybrid Model

highly nonlinear crude property curves one may need to employ the hybrid model/TBP-cuts method. The results of this work show that in most cases the least amount of nonlinearity will be introduced in the planning model of CDU if the conventional swing-cut method based on the yields from the hybrid model is employed to compute product bulk properties.



L/H cut

size (vol %)



REFERENCES

(1) Getting started: modeling petroleum processes; Aspen Technology Inc.: Cambridge, MA, 2006. (2) Brooks, R. W.; Van Walsem, F. D.; Drury, J. Choosing Cut-points to Optimize Product Yields. Hydrocarbon Processing 1999, 78 (11), 53−60. (3) Pinto, J. M.; Joly, M.; Moro, L. F. L. Planning and Scheduling Models for Refinery Operations. Comput. Chem. Eng. 2000, 24 (9− 10), 2259−2276. (4) Neiro, S. M. S.; Pinto, J. M. A General Modeling Framework for the Operational Planning the Petroleum Supply Chain. Comput. Chem. Eng. 2004, 28, 871−896. (5) Zhang, J.; Zhu, X.; Towler, G. A Level-by-Level Debottlenecking Approach in Refinery Operation. Ind. Eng. Chem. Res. 2001, 40, 1528− 1540.

AUTHOR INFORMATION

Corresponding Author

*Tel.: +1 905 525 9140 ext. 26386. E-mail address: mahalec@ mcmaster.ca. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS Support by Ontario Research Foundation and McMaster Advanced Control Consortium is gratefully acknowledged. K

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research (6) Li, W.; Hui, C.-W.; Li, A. Integrating CDU, FCC and product blending models into refinery planning. Comput. Chem. Eng. 2005, 29 (9), 2010−2028. (7) Guerra, O. J.; Le Roux, A. C. Improvements in Petroleum Refinery Planning: 1. Formulation of Process models. Ind. Eng. Chem. Res. 2011, 50, 13403−13418. (8) Guerra, O. J.; Le Roux, A. C. Improvements in Petroleum Refinery Planning: 2. Case studies. Ind. Eng. Chem. Res. 2011, 50, 13419−13426. (9) Menezes, B. C.; Kelly, J. D.; Grossmann, I. E. Improved Swing Cut Modeling for Planning and Scheduling of Oil-Refinery Process. Ind. Eng. Chem. Res. 2013, 52, 18324−18333. (10) Caballero, J. A.; Grossmann, I. E. Aggregated Models for Integrated Distillation Systems. Ind. Eng. Chem. Res. 1999, 38, 2330− 2344. (11) Gadalla, M.; Jobson, M.; Smith, R. Shortcut models for retrofit design of distillation columns. Trans I ChemE, Part A 2003, 81, 971− 986. (12) Alattas, A. M.; Grossmann, I. E.; Palou-Rivera, I. Integration of Nonlinear Crude Distillation Unit Models in Refinery Planning Optimization. Ind. Eng. Chem. Res. 2011, 50, 6860−6870. (13) Mahalec, V.; Sanchez, Y. Inferential Monitoring and Optimization of Crude Separation Units via Hybrid Models. Comput. Chem. Eng. 2012, 45, 15−26. (14) Ochoa-Estopier, L. M.; Jobson, M.; Smith, R. The use of reduced models for design and optimization of heat-integrated crude oil distillation systems. Energy 2014, 75, 5−13. (15) Kelly, J. D.; Menezes, B. C.; Grossmann, I. E. Distillation blending and cutpoint temperature optimization using monotonic interpolation. Ind. Eng. Chem. Res. 2014, 53, 15146−15156. (16) Sanchez, S.; Ancheyta, J.; McCaffrey, W. C. Comparison of Probability Distribution Functions for Fitting Distillation Curves of Petroleum. Energy Fuels 2007, 21, 2955−2963.

L

DOI: 10.1021/acs.iecr.5b02919 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX