Equation-Oriented Approach for Handling the Perturbed-Chain SAFT

Subscriber access provided by - Access paid by the | UCSB Libraries

Process Systems Engineering

An Equation-Oriented Approach for Handling Perturbed-Chain SAFT Equation of State in Simulation and Optimization of Polymerization Processes Jiayuan Kang, Lingyu Zhu, Shenjun Xu, Zhijiang Shao, and Xi Chen Ind. Eng. Chem. Res., Just Accepted Manuscript • Publication Date (Web): 08 Mar 2018 Downloaded from http://pubs.acs.org on March 8, 2018

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

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

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

Industrial & Engineering Chemistry Research

An Equation-Oriented Approach for Handling Perturbed-Chain SAFT Equation of State in Simulation and Optimization of Polymerization Processes Jiayuan Kanga, Lingyu Zhub, Shenjun Xub, Zhijiang Shaoa, Xi Chena∗ a

State Key Laboratory of Industrial Control Technology, College of Control Science and Engineering, Zhejiang University, Hangzhou, Zhejiang, 310027, China b

College of Chemical Engineering, Zhejiang University of Technology, Hangzhou 310014, China

Abstract: The Perturbed-Chain SAFT (PC-SAFT) equation of state (EOS) is a very popular and promising model for fluids. Resolution of the PC-SAFT EOS is generally conducted by subroutine-based calculations that involve logical conditions. This approach can lead to nonsmoothness and convergence issues when used with gradient-based solvers in an equation-oriented (EO) framework. In this paper, we propose a novel EO approach for the PC-SAFT EOS embedded in large flowsheet optimization. First, the mathematical structure of the PC-SAFT EOS is analyzed through a digraph. It reveals that for polymerization systems, the resolution of the PC-SAFT EOS is not only coupled by phase equilibrium, but also by the polymerization kinetics. This feature strongly motivates implementation of both EOS resolution and process optimization in a complete EO framework. A novel EO approach is proposed for selecting the appropriate root of the PC-SAFT EOS by incorporating thermal stability criteria to eliminate the undesirable roots. Numerical examples in simulation and optimization of flash drums are presented to demonstrate that the proposed EO formulation can always select the appropriate root. Finally, the EO formulation is applied for simulation and optimization of an industrial polymerization process. The results show advantages of

∗

Corresponding author. Tel.: 86-571-87953966; Email address: [email protected] 1

ACS Paragon Plus Environment

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

the proposed EO method compared to traditional sequential modular based simulation software in terms of computational robustness and efficiency. Keywords: PC-SAFT EOS; polymerization; process optimization; phase equilibrium; thermodynamics

2


Page 2 of 61

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


1. Introduction An equation of state (EOS) is a property model that relates two or more state functions to describe the thermodynamic behaviors in fluids. Since the pioneering work of Van der Waals [1], a number of EOS have been developed. Cubic equations of state (CEOS) were widely used in the prediction of vapor-liquid equilibria for non-polar components in the 1970s. The mutations of CEOS including Soave-Redlich-Kwong (SRK) equation [2] and Peng-Robinson (PR) EOS [3] improve their predictive performance by modify the alpha functions and mixing rules. Given that the use of CEOS is restricted to systems of simple fluids, more theoretically based EOS have been proposed for complex and macromolecular compounds in recent years. A category of EOS for chain molecules, i.e., the well-known statistical associating fluid theory (SAFT) EOS [4] has been proposed and developed. Among many modifications of the SAFT EOS, perturbed-chain SAFT (PC-SAFT) EOS [5] is among the more popular in recent years. An obvious advantage of PC-SAFT EOS is that it is applicable to both small spherical molecules and chainlike polymers. This model has attracted great attention due to its industrial relevance as a predictive tool. The property models are generally used for process synthesis, design and control. In these cases, the EOS models are used in the so-called “service role” [6]. That is, the required properties are calculated from the EOS when asked. In this framework, the values of the intensive variables (such as temperature, pressure and/or composition) are input to the EOS model, which is treated as a procedure and returns calculated property values. Then, the simulator forms an inner calculation loop that is repeated for every request for property values by an outer calculation loop, which occurs when any of the intensive variable changes. This approach, or the so-called sequential modular (SM) approach, is widely employed in most commercial software packages for process synthesis and design. 3



Resolution of an EOS is the procedure to calculate the required properties once the degrees of freedom of a fluid system have been fixed. Given the high order polynomials in many equations of state, resolution often leads to multiple solutions. How to select the appropriate solution is the key step to avoid inconsistent predictions of phase behavior. Taking the classical CEOS as an example, the cubic equations can be solved either analytically or numerically and in either case, certain guidelines are employed to select the appropriate root. One common guideline is that when three real roots of compressibility factor exist, the largest root is the vapor root and the smallest is the liquid root. This kind of typical if-else condition is embedded in most algorithms for EOS resolution. In contrast to CEOS, the PC-SAFT EOS employs a complicated pressure-explicit function with iterative nature combined with derivatives of state functions. Therefore, PC-SAFT EOS cannot be solved analytically at a specified pressure and temperature, and it is difficult to determine how many roots may exist in theory. In 2005, Yelash et al. [7] mentioned that the PC-SAFT EOS may exhibit “artificial multiple criticality and phase equilibria”. Moreover, Privat et al. [8] pointed out that it in case of pure fluids, the PC-SAFT EOS may exhibit up to five different volume roots. They also highlight that the resolution algorithms should be capable of detecting and managing more than three molar volume roots. Within the framework of SM approach for process simulation and optimization, a number of algorithms for EOS resolution have been proposed. Bausa and Marquardt [9] proposed a method for phase stability test and flash calculation by homotopy continuation. Gani et al. [10] presented a mathematical analysis of classed of property models, including the CEOS and PC-SAFT EOS, and propose solution strategies such as simultaneous or nested iterative loops based on the ordering of these equations. In these methods, an optimizer is coupled with a simulator in which the EOS is treated as a

4


Page 4 of 61

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


black-box procedure. However, because of the interactions between the modules and optimizer, input-output derivatives for the implicit models are not directly available. Instead, finite differences are usually employed, which can lead to inaccuracy and round-off errors. Moreover, in SM mode, each module needs to be solved repeatedly and the intermediate failure in each inner loop is detrimental to the outer loop. In contrast, the equation-oriented (EO) approach allows for simultaneous optimization, property calculation and flowsheet convergence. By employing exact derivatives and advanced mathematical programming algorithms, the EO approach is preferred for complex flowsheets with many nested recycle loops and implicit design specifications. A comprehensive review of SM and EO approach in process optimization is provided in Biegler, Grossmann, and Westerberg [11]. It is suggested that the analytical method will be suitable for the SM approach and the numerical method will be preferred for the EO approach. However, most of the algorithms for EOS resolution utilized analytical procedures to evaluate the equations. These procedure-based calculations may lead to non-smoothness and failure within an EO optimization framework. In 2010, Kamath, Biegler and Grossmann [12] investigated the special mathematical feature of cubic equation and proposed a systematic formulation for CEOS within the EO framework. The benefits of this formulation over the SM method in commercial software were also demonstrated in our previous work [13] through the optimization of an industrial cryogenic air separation unit. In this paper, we focus on an EO framework that combines resolution of PC-SAFT EOS and process optimization. Although algorithms for resolution of PC-SAFT EOS following certain guidelines or heuristics in an external procedure have been applied in commercial software, we argue that these are not suitable for use with an EO approach. Thus, proper constraints are needed so that the PC-SAFT

5



EOS can be expressed in a complete EO form and easily integrated with the rest of the EO process model. An intuitive analysis of the mathematical structure of the PC-SAFT EOS is presented through a digraph. It reveals that for polymerization systems where the molecular weight of a chemical component is determined by reaction, the property calculation is not only coupled with phase equilibrium, but also with the polymerization kinetics. This feature motivates the implementation of both EOS resolution and process optimization in a complete EO framework. PC-SAFT EOS employs a pressure-explicit function and has a complex mathematical formulation (such as the derivatives of the Helmholtz energy). Although efforts have been made in the past to reduce the computational cost and programming efforts for PC-SAFT EOS, the focus was on how to simplify some terms of the model equations. Solms et. al. modified the original PC-SAFT EOS and develop the simplified PC-SAFT EOS by simplifying the radial distribution function and hard-sphere contribution term [14]. Tihic et. al developed a group-contribution simplified PC-SAFT EOS and applied to phase equilibrium calculations of various polymer systems [15]. These modifications simplified part of the expressions in PC-SAFT EOS, but from the point of view of numerical analysis, the model structures remain the same. Thus it is worthy to investigate the proper implementation of PC-SAFT EOS for process simulation and optimization. The paper is organized as follows. In Section 2, the EO framework is presented, with an emphasis on the analysis of the mathematical structure of the PC-SAFT EOS. In Section 3, an EO formulation is presented for the resolution of PC-SAFT EOS. The roots of PC-SAFT EOS are selected by incorporating thermal stability criteria to eliminate the trivial phase equilibrium solutions. This formulation is tested in simulation and optimization of a flash drum. In Section 4, the EO formulation of the PC-SAFT EOS is applied for the simulation and optimization of an industrial polymerization

6


Page 6 of 61

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


process. The model structure of the process optimization presents strong nonlinearity with multilayer coupling among polymerization kinetics, EOS resolution, and flowsheet calculation. The results show advantages of the proposed EO method compared to traditional SM-based commercial software in terms of computational robustness and efficiency.

2. An Equation Oriented Framework for Flowsheet Optimization with PC-SAFT EOS Polyolefins are generally produced in solution, in slurry, or in gas-phase reactors. In solution and slurry reactors, monomer, catalyst and chain transfer agent are mixed in solvent or diluent, where the polymer chain grows. These reactors usually contain both liquid phase and gas phase with gas phase recycle loops. The molecular weight of the polyolefin is controlled by the chain transfer agent (usually hydrogen) dissolved in the liquid phase. Therefore, it is critical to predict the phase behavior of the polymerization reactor for process simulation and optimization. Given that the polymerization system is a mixture of small spherical molecules and chainlike polymers, the PC-SAFT EOS is widely chosen to describe the thermodynamics of the system. 2.1 Mathematical and numerical analysis of PC-SAFT EOS The Helmholtz free energy is the starting point of the PC-SAFT EOS, as all the other properties can be obtained as derivatives of . PC-SAFT EOS presents a dispersion expression for chain

molecules and uses a hard-chain reference fluid. The reduced residual molar Helmholtz energy is then obtained from the following expression: = + +

(1)

where is the mean segment (hard sphere) number in the mixture.

The hard-chain reference contribution, i.e., the first and second term in the right hand side of Eq. (1), 7



Page 8 of 61

represents the Helmholtz energy of an elementary chain by adding a hard-sphere term to a chain term that accounts for the formation of the chain. The third term of the right hand side of Eq. (1) is a dispersion contribution which reflects interactions between hard chains. The molar compressibility factor is derived using the thermodynamic relation:

= 1 + (

) ,

(2)

where is the reduced molar density of the mixture, is the temperature, and is the set of composition of the mixture. Comparing Eqs. (1), (2), and the PVT relation leads to (, ) = 1 + + = !/#

(3)

= +

(4)

where is the molar volume of the mixture, and similarly,

Eq. (2) indicates that the PC-SAFT EOS employs a pressure-explicit function. Due to the complexity of the mathematical formulation (such as the derivatives of the Helmholtz energy), it is necessary to solve this EOS numerically at given pressure and temperature. This feature makes the PC-SAFT EOS quite different from classical CEOS which can be easily solved analytically using the Cardano method [16]. To demonstrate the numerical procedure of solving the EOS, a mathematical analysis of the original PC-SAFT model [5] is presented in an intuitive form in Fig. 1. The equations of the model are given in Appendix A. Fig. 1. The digraph of the variables in PC-SAFT EOS Fig. 1. represents the information flow of the PC-SAFT EOS in a digraph, or directed graph. Draw a directed edge from node a to node b, if and only if, the variable in node a is required explicitly to calculate the variable in node b. Note that each directed edge does not represent a unique equation; they

8


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


just indicate the relationship of the variables. By using a digraph, one can intuitively distinguish the dependence of the variables in a set of nonlinear equations, the assigned input and output variables, and the intermediate variables. Fig. 1. reflects the common approach to solve the PC-SAFT EOS, where an iterative loop (i.e. the directed edges in dashed lines) is designed to calculate the properties at given system pressure and temperature. Despite the algorithm part in dashed edges, it presents three layers in the digraph: the input layer, the intermediate layer, and the output layer. The input layer consists of independent parameters, input variables provided by outer loop, and a guess of the reduced density . The nodes in green in the bottom layer of the digraph are the independent

parameters including the depth of pair potential $ characterizing the dispersion forces, the segment diameter % , the regression coefficients & and '& , where subscript ( is the component index and )

the seven coefficients given in [5]. The nodes in yellow in the bottom layer of the digraph are the input variables including the temperature , the mole fraction of the components * , and the number of

segments per chain in each component . Note that in most of the literature is treated as a

parameter, as it is determined by the inherent attribution of the compound, i.e., the elementary chain formed from elementary segments. In a polymerization flowsheet, however, should be treated as a variable because the length of the polymer chain varies during the reaction and consequently changes the number of segments per chain ( ). It is believed that depends on the molecular weight of a polymer in the form of = + , -!. , ( ∈ {12345+}

(5)

where , is the molecular weight of monomer, -!. is the number average degree of

polymerization, and + is a chain-independent parameter that is regressed individually. 9



Page 10 of 61

The intermediate layer is formed by the white nodes in Fig. 1. This layer is the principle part of the PC-SAFT EOS, where tens of variables and derivatives are defined. As shown in Fig. 1, no internal loop is found in the intermediate layer. This feature indicates complete decoupling in the intermediate layer. All variables can be obtained following the flow of directed edges through forward calculation. In this way, the number density of the mixture 7 and the compressibility factor are calculated.

Given 7 and in the intermediate layer, one can calculate the pressure of the mixture by

where 9 is the Boltzmann constant.

! 8 = 97 × 10 is known and should be the input instead of output in the resolution procedure. Hence, it is suggested by Gross that an iterative loop is required to adjust the reduced molar density until ! 8 = !> [5]. This iterative loop is the so-called algorithm part (i.e. the directed edges in dashed lines in Fig. 1.) in the PC-SAFT EOS. In this loop, a typical tear-variable approach is employed, where the reduced molar density in the input layer is the tear variable and is adjusted in each iteration. Thus, a direct calculation is ensured without solving nonlinear equations. The digraph in Fig. 1. presents the minimum form to solve the PC-SAFT EOS, where only essential variables and equations are involved. Besides the number density and compressibility factor, other

properties such as fugacity coefficient can also be obtained through forward calculation following the sequence of equations [8] in Appendix A.

10


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


2.2 EO framework for flowsheet optimization with PC-SAFT EOS

As shown in Fig. 1., the original PC-SAFT EOS includes an iterative loop to calculate the properties. The EOS is typically implemented under the SM framework for process simulation and optimization. A brief representation of the SM framework is shown in Fig. 2(a) following the instructions on the PC-SAFT implementation documents. Under the SM mode, each operation unit has a self-contained solver. The access to the internal properties is provided by the external property server, which is associated with each port of the unit modular. The property server generally contains of two routines, namely, the equilibrium calculations and property calculations. The equilibrium routine serves as an interface between the unit and properties. An iteration subroutine is usually embedded to ensure the chemical potential of each components to be the same in all existing phase. The iteration subroutine requests fugacity coefficients of all components in all phases from the property model. The property routine contains the PC-SAFT EOS together with the algorithm (see Fig. 1). With different calculation tasks and system specifications (i. e. P-x, P-T or vapor-liquid), the resolution of the PC-SAFT EOS can be conducted by specific subroutines to select the appropriate roots. The calling sequence of the PC-SAFT EOS under the SM framework is investigated by comparing Fig. 1 and 2(a). The input of the PC-SAFT EOS can be obtained from the outer loops: is usually a manipulated variable, * can be provided by phase equilibrium in the outer layer, and is to be

calculated by polymerization kinetics (or further the flowsheet connections for reactor networks) in the outermost layer. From the point of view of model partitioning and the SM approach, the input variables in the PC-SAFT EOS require feedback of two different layers of outer loop. This kind of outer loop dependence will lead to higher computational cost and convergence difficulties. The nature of the polymerization contributes to a coupling of phase equilibrium and polymerization kinetics (or even 11



Page 12 of 61

flowsheet calculations) to solve the EOS. This complexity and coupling feature indicates a strong incentive to implement both EOS resolution and process optimization in a complete EO framework. In the EO framework as shown in Fig. 2(b), the process model is completely expressed in the form of nonlinear equations. Here, the properties are not in the service-role but part of the mathematical model. The model is transcript to numerical equation set and can be simultaneously solved by nonlinear algorithms. The mathematical and numerical analysis of the PC-SAFT EOS reveals the essential equations that should be simultaneously solved. We write the essential equations in an implicit form as !?@AB(!> , , , C, ) = 0

(7)

where C is set of segments number per chain ( )

Fig. 2. SM and EO framework in process simulation Eq. (7) is the EO form of the PC-SAFT EOS. The five (sets of) variables involved in Eq. (7) are the variables in the input and output layer of Fig. 1., i.e., the system pressure !> , the temperature , the

component fraction * , the reduced molar density , and the number of segments per chain in each

component . Eq. (7) is has 2.? + 2 degrees of freedom (DOFs), where .? + 2 variables are

fixed by the outer problem (i.e., , !> , * ), and .? variables are fixed by the component attributes

(i.e., , which may change along the reaction). Therefore, the only variable that must be solved by Eq.

(7) is . The rest of the variables in the EOS are treated as intermediate variables which can be obtained through forward calculation.

In state-of-art EO modelling environment such as AMPL and GAMS, the intermediate variables can be eliminated by presolve and substitution functions. Moreover, the exact derivatives are provided via automatic differentiation. These benefits allow simultaneously solving set of equations with embedded property models. Compared with the nested loops in the SM approach, computational efficiency is

12


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


significantly improved. An EO framework for flowsheet optimization with embedded PC-SAFT EOS is shown in Fig. 3. The EO model consists of mass balance, phase equilibrium, summation equation and energy equilibrium (MESH), together with the reaction kinetics and flowsheet connections. The PC-SAFT EOS is presented as independent equations in the form of Eq. (7), thus avoids the conventional “procedure equations” in an iterative form. Given some DOF of the PC-SAFT EOS, i.e. the variable for polymers are determined by the DPN as described in Eq. (5), the polymerization kinetics and flowsheet connections are required to be simultaneously solved. The explicit equations allow the computation of Jacobian and Hessian of the equation system, which are required by the advanced nonlinear programming (NLP) solvers such as IPOPT [17]. This EO framework has been prototyped in the AMPL environment, which provides exact first and second derivatives via automatic differentiation. IPOPT is used as nonlinear programming solver. A concise EO formulation of a process model is presented as follows: !?@AB(!, , , C, F ) = 0

!?@AB(!, , G, C, H ) = 0 I = I(, , G, , F , H )

J(!, , K, I, , G, C) = L

(8a) (8b) (8c) (8d)

where Eqs. (8a)(8b) represent the essential part of equations in the PC-SAFT EOS that must be solved simultaneously (i.e., solving the reduced molar densities F and H in liquid and vapor phase), Eq. (8c) is the rest of the equations in the PC-SAFT EOS that can be obtained through forward calculation, where K are a set of properties that can be calculated from the reduced molar densities, and J in Eq. (8d) are set of governing equations for the process models with process variables U, including mass

13



balances, energy balances, summations, polymerization kinetics and flowsheet connections. Fig. 3. EO framework for flowsheet optimization with embedded PC-SAFT EOS

3. Simulation and Optimization of a Flash Drum One of the most important applications of the property model is the flash calculation. The flash drum in the polymerization process controls the phase behaviors of the feed flows. Given that the operating conditions and components are usually unable to meet the conditions of the polymerization, only phase equilibrium is considered in the flash drums. The polymer chains in the flash drum will not grow and the variable for polymers are fixed. Therefore, the simulation and optimization of the flash drum is a practical way to have an insight into the PC-SAFT EOS before talking about the large process simulation and optimization. 3.1 Problem statement

In this section, we consider the flash drum for a monomer/diluent separation process. In this process, the polymer product slurry is fed into a separation unit to remove the unreacted monomer and diluent. The monomer and diluent are recovered while the polymer product is dried and pelletized. The manipulated variables in this process are the temperature and pressure of the flash drum to control the phase behavior of fluid. The components in the flash drum are ethylene, hydrogen, hexane, high-density polyethylene (HDPE), titanium tetrachloride (TiCl4, the catalyst) and triethyl aluminum (TEA, the co-catalyst). The three pure-component parameters of the PC-SAFT model (the number of segments per chain in each component . the depth of pair potential $ characterizing the dispersion forces, and the segment diameter % ) are listed in Table 1. The feed flow of the flash drum is shown in Table 2. The DPN of the HDPE is 1093. 14


Page 14 of 61

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


Table 1. Pure component parameters of polyethylene system. Table 2. Feed flow of the flash drum 3.2 PC-SAFT EOS resolution

The simulation of the flash drum can be divided into two procedures, i.e. phase equilibrium and EOS resolution. As analyzed before, the resolution of PC-SAFT EOS predicts the reduced molar density at given system temperature, pressure and component fraction. In this section, an isotherm for ethylene and an isobar for HDPE mixture are presented. Two different algorithms were implemented for PC-SAFT EOS resolution were implemented in MATLAB for a comparison of computation consistency and efficiency. One is the original sequential modular algorithm, implemented according to the instructions provided by the author of PC-SAFT, and the other is the fmincon algorithm applied to the EO PC-SAFT model. Fig 4. presents the pressure-density diagram for pure component of ethylene at T = 250 K. The two algorithms give consistent predictions of the isotherm. A sharp change of density is found at about P = 20 bar. This transition point indicates a phase change from vapor to liquid, where the density of the fluid increase dramatically and then gradually tends to approximate the critical density. Fig. 4. Pressure-density diagram for ethylene at T = 250 K Fig. 5. presents the temperature-density diagram for HDPE mixtures (see Table 2) including ethylene, hydrogen, hexane, HDPE, catalyst and co-catalyst at P = 4 bar. The two resolution algorithms also provide consistent predictions of the isobar. A sharp change of density is found at about T = 353 K. This transition point indicates a phase change from liquid to vapor, where the density of the fluid decrease dramatically and then gradually tends to approximate the critical density. Fig. 5. Temperature-density diagram for HDPE at P = 4 bar 15



Page 16 of 61

Fig. 6. compares the computational consistency and efficiency of the SM and EO algorithm in PC-EOS resolution. The convergence of 40 bar on the ethylene isotherm and 343 K on the HDPE isobar was presented. The initial guess of the reduced molar density in both the two algorithms are set as 0.6. The symbols represent the value of in each iteration of EOS resolution; the triangles and circles

represent the case of ethylene isotherm using SM and EO algorithm, respectively, and the filled and open diamonds represent the case of HDPE isobar using SM and EO algorithm, respectively. For both the two cases, the EO algorithm takes fewer iterations than the SM algorithm in the same convergence tolerance. Fig. 6. Comparison of computational consistence and efficiency of SM and EO algorithm 3.3 Flash calculation using PC-SAFT EOS directly

In the previous section, the phase equilibrium was not considered, thus the isotherm and isobar represent transition regions. These regions indicate the phase change that should be calculated by phase equilibrium equations. In the current research, only vapor-liquid equilibrium (VLE) is considered. Modelling the vanishing phase [18] is not in the scope of this research. The following set of equations is solved: Formulation I:

M = 3 + , ∀( ∈ ? O = ∑ ∈Q 3

3 = * O, ∀( ∈ ? R = ∑ ∈Q

= 4 R, ∀( ∈ ?

4 S H = * S F , ∀( ∈ ?

!?@AB(!, , , C, F ) = 0 16


(9a) (9b) (9c) (9d) (9e) (9f) (9g)

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


!?@AB(!, , G, C, H ) = 0

S F = T (, , C, F ), ∀( ∈ ?

S H = T (, G, C, H ), ∀( ∈ ?

(9h) (9i) (9j)

where M is the overall component flowrate; 3 and represent component flowrates in the liquid and vapor phases, respectively; and C is the set of components; O and R are overall molar flowrates

for both phases; x and y are sets of the compositions * and 4 in the liquid and vapor phases,

respectively; and S is the fugacity coefficient of the i th component where the superscript L and V

represent the liquid and vapor phases, respectively. Eqs.(9a)-(9f) are the general VLE equations including mass balance and phase equilibrium, and Eqs. (9g)-(9j) are the property calculation with embedded PC-SAFT EOS. The DOF of this system is 2.? + 2, where .? variables are fixed by the

overall molar flow rate (i.e. M ), .? variables are fixed by the component attributes (i.e. ), and 2 variables are fixed by the phase equilibrium problem (i.e. and !).

The flash calculations in the SM and EO frameworks are treated differently. In the SM framework, EOS resolution is an individual procedure where the EOS models are used in a service role. The phase equilibrium subroutine requests the property results for a given phase from the EOS model. However, the transition behavior in the region of phase change always indicates two or more property results that satisfy the EOS model. Thus, an important procedure in EOS resolution is choosing the appropriate property result for the given phase. This procedure involves logical conditions such as if-else and min-max, which have derivative discontinuities. In the EO framework, the EOS model is embedded in the phase equilibrium, i.e., two sets of EOS models for each phase are used to calculate the properties. Eqs. (9a)-(9j), denoted as Formulation I is an example of EO formulation for VLE calculation with embedded PC-SAFT EOS. The logical

17



conditions to determine the property results for different phase are thus avoided. It is important to note that specifying a number of variables equal to the DOF in Eqs. (9a)- (9j) does not guarantee uniqueness of the solution. In other words, existence of multiple phase equilibriums at fixed temperature and pressure does not violate the Gibbs’s phase rule. Consequently, specific formulation is required to select the appropriate solution to ensure the physical consistency. Fig. 7(a). presents a case of VLE calculation for the flash drum in section 3.1 using Formulation I. The model was built in the AMPL environment and IPOPT was chosen as the nonlinear solver. The temperature of the flash drum is set as T = 363.15 K and the pressure is set as P = 4 bar. The value of

ranges between 0 for an ideal gas and ≤ 0.74 for the closest packing of segments [5]. Thus it is

reasonable to use an initial guess of H in the interval of [0 0.05] for vapor phase evaluation and F

in the interval of [0.2 0.45] for liquid phase evaluation. The simulation was run 250 times from

different initial guesses of F and H . The convergence tolerance in the IPOPT was set as 10-10. The computation status is shown in Fig. 7(a). The solid dots and the hollow dots represent converged runs and diverged runs with the given initial guess, respectively. The converged solution of F and H are represented by stars. Fig. 7. Flash calculation status using different formulations Perhaps surprisingly, the flash calculation lead to several different solutions depending on the initial guess in Formulation I. All the solutions (star points in Fig. 7(a)) satisfy the constraints of equations in Formulation I, but only one of them is physically consistent. This can be explained by the low sensitivity of the energy parameter of PC-SAFT EOS toward liquid densities of highly nonspherical molecules [19]. This suggests a need for a more specific formulation for flash calculations with long-chain polymers. Another reason lies in the high order polynomials in the PC-SAFT EOS and one

18


Page 18 of 61

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


cannot tell how many roots it will predict for such a complex equation system. Researchers have found that PC-SAFT EOS can predict up to five volume roots for pure components [8], and it is unknown how this conclusion may change for phase equilibrium of mixtures. 3.4 EO formulation to avoid nonphysical phase equilibrium solutions

To select the appropriate root in the flash calculation, it is helpful to recall that the VLE calculation is a Gibbs free energy minimization problem [20]. For most nonideal thermodynamic models, the Gibbs free energy minimization problem is non-convex and may have multiple local solutions. The global minimum is preferred, as it is the most stable energy state. The local minimum does not ensure the smallest Gibbs free energy for a VLE system, and thus tend to be unstable. The VLE calculation using Formulation I is in fact the necessary optimality condition of the original Gibbs free energy minimization problem. Thus formulation I may predict unstable results. Though these local solutions can be identified and compared using multi-start initialization [21], it is not guaranteed to find all the solutions. Alternately, several studies apply rigorous, deterministic global optimization algorithms to VLE calculations [22]. Although global optimization guarantees the correct VLE calculations, it is still difficult to efficiently embed these calculations in an EO flowsheet framework. It is recommended that careful initialization is the best safeguard against trivial solutions [23]. Therefore, this work focuses on computationally efficient, robust and flexible approaches to avoid a specific class of false VLE calculations. The necessary stability condition [8] of a VLE is given as \

] \

_

>0

(10a)

>0

(10b)

F − H > 0 19


(10c)


to avoid the trivial solutions predicted Fig. 7(a). The calculation of

\

is given in Appendix B.

With these conditions, stability analysis of the multiple solutions in Fig. 7(a) is conducted by verifying the derivatives of pressure towards reduced molar density. Fig. 8. shows the VLE stability of the solutions in vapor and liquid phases using Formulation I. It is found that the solutions where F ≤ 0.05 present negative derivatives

\

]

. These solutions correspond to the upper left solution

ridge in Fig. 7(a), where F and H share the same value. This phenomenon indicates a false VLE calculation where the properties are almost the same in both phases and the equilibrium is unstable. Thus, all the solutions in this ridge must be discarded. On the other hand, the solutions where H ≥ 0.007 also present negative derivatives

\

_

as shown in Fig. 8. This suggests that these

solutions should be discarded as well. Fig. 8. Stability analysis of flash calculations using Formulation I By adding the stability conditions Eq. (10) as constraints, a modified VLE formulation is proposed, which is denoted as Formulation II. Formulation II:

Eqs. (9)-(10)

The flash calculation problem in Section 3.3 is solved again using Formulation II. The simulation was run 250 times from different initial guesses of F and H . The computation status is shown in Fig. 7(b). In this diagram, only one ridge of solutions is found. The solutions with indistinguishable phases are eliminated but multiple solutions still exist. The stability constraints are satisfied, so it seems likely that the remaining multiple solutions in Fig. 8. are caused by the low sensitivity of the energy parameter of PC-SAFT EOS toward liquid densities of highly nonspherical molecules. Given that high-molecular-weight compounds, such as polymers, have very low vapor pressure and are non-volatile, another modification of VLE formulation is derived. In

20


Page 20 of 61

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


this formulation, only light components are considered in the phase equilibrium, whereas the PC-SAFT EOS still uses all the components to calculate the properties. This light component VLE formulation is given in the following. M = 3 + , ∀( ∈ ? O = ∑ ∈Q 3

3 = * O, ∀( ∈ ? R = ∑ ∈Q

= 4 R, ∀( ∈ ?

4 S H = * S F , ∀( ∈ ? − b 4 = 0, ∀( ∈ b

!?@AB(!, , , C, F ) = 0

!?@AB(!, , G, C, H ) = 0

S F = T (, , C, F ), ∀( ∈ ? − b

S H = T (, G, C, H ), ∀( ∈ ? − b

(11a) (11b) (11c) (11d) (11e) (11f) (11g) (11h) (11i) (11j) (11k)

where H is the set of non-volatile components (such as HDPE, catalyst and cocatalyst). The phase equilibrium is described by Eq. (11f) instead of Eq. (9f), and the property calculation in Eqs. (11h)-(11i) remains the same. Note that 4 for non-volatile components is fixed to zero by Eq. (11g), so the DOF remains the same as Eq. (9). By adding the stability conditions Eq. (10) as constraints, another modification is derived, which is denoted as Formulation III. Formulation III:

Eqs.(10)-(11)

Formulation III is expected to eliminate the inappropriate solution by combing light component VLE

21



Page 22 of 61

and stability criteria. The flash calculation problem in Section 3.3 is solved again using Formulation III. The simulation was run 250 times from different initial guesses of F and H . The computation status is shown in Fig. 7(c). In this diagram, a unique solution exists. This solution is exactly the flash calculation result conducted by the PC-SAFT packages POLYPCSF embedded in Aspen Plus 7.2. Thus the appropriate root in the flash calculation is selected. The results of flash calculation using different formulations are analyzed. A comparison of the appropriate solution and some abnormal solutions is presented in Table 3. The fifth column is the appropriate result of flash calculation conducted from Formulation III, which is consistent with Aspen Plus 7.2. The abnormal solution presented in the second column is a solution calculated from Formulation I as shown in Fig. 7(a). This solution predicts an unstable H and

\

_

is found to be

negative in Fig. 8. The compressibility factor of vapor phase H is found to be too small compared to

the ideal gas. The mean segment number of the vapor mixture ( H ) is even larger than that of the

liquid mixture ( F ). The calculation result is obviously nonphysical and should be discarded. Another

abnormal solution presented in the third column is a solution calculated from Formulation II as shown in Fig. 7(b). This solution satisfies the stability constraints but predicts a low yield of HDPE in liquid phase, which indicates considerable amount of HDPE go up into the vapor phase. Thus the mean segment number of the vapor mixture ( H ) is almost the same as that of the liquid mixture ( F ). Though the other variables of the VLE calculation is close to the appropriate results in the fifth column, this solution is still nonphysical and should be carefully discarded. The abnormal solution presented in the fourth column is a solution calculated from Formulation III without stability constraints, i.e., the light component VLE (Eqs. 11) only. Though Eqs. 11 restrict the mole fraction of heavy component such as HDPE in vapor phase to zero, it predicts almost the same value of F and H . Thus most of 22


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


the properties in the vapor and liquid phase, i.e. S H , S F , H and F are almost the same. This unstable equilibrium is nonphysical and should be discarded. Hence, it is necessary to introduce the stability constraints in the light component VLE formulation. Table 3. Comparison of nominal solution and abnormal solutions 3.5 Flash drum optimization case study

This case study deals with steady state optimization of a flash drum. The flash drum is used to remove the unreacted monomer from the HDPE slurry product. The objective is to minimize the amount of monomer in the liquid phase. The manipulated variables in this case are the temperature and pressure of the flash drum to control the phase behavior of fluid. The two variables play very important roles in the resolution of PC-SAFT EOS. Thus it is better to test the reliability of the modified EO formulation is tested through optimization than the simulation case. The optimization problem is presented as follows: min *Qfgh O ,\

s. t. 335 ≤ T ≤ 370 3≤P≤5

,gn\o ≤ p

Eqs. (10)-(11). (Formulation III)

(12)

where c is the solid fraction constraint for HDPE. We expect that Formulation III can select the appropriate solution for flash calculations during any iteration of the optimization or at the optimal solution. Three cases differing in the solid fraction constraint for HDPE are presented as shown in Table 4. The same initial guess is used for all the three cases, i.e., an infeasible point of T = 313.15 K, P = 2 atm, F = 0.3, and H = 0.001. All three cases could be reliably solved using the NLP solver IPOPT with CPU times of less than 0.1 s. The optimization results are validated by Aspen Plus 7.2 to test the computational consistence. As shown in 23



Table 4., the optimal objectives and operating conditions of the EO formulation consist with the results of Aspen Plus 7.2. This indicates that the proposed EO formulation is robust since the solver can find the appropriate optimal solutions from an infeasible starting point. Table 4. Three cases of flash optimization

4. Simulation and Optimization of an Industrial Polymerization Process In this section, the EO formulation of PC-SAFT EOS is validated in a slurry polymerization. In the slurry reactors, monomer, catalyst and chain transfer agent are mixed in diluent. These reactors present both liquid phase and gas phase with recycle gas loops. The molecular weight of the polyolefin is controlled by the concentration of chain transfer agent (usually hydrogen) dissolved in the liquid phase. Therefore, it is critical to predict the phase behavior of the polymerization reactor for process simulation and optimization. Different from the flash drum case where the property variables are fixed, the polymer chains

grow and for polymers are to be calculated by the polymerization kinetics (or further the flowsheet

connections in reactor networks). As analyzed in Section 2, the nature of the polymerization process contributes to a coupling of phase equilibrium and polymerization kinetics. This interaction leads to higher computational complexity in SM mode rather than a simple superposition of equations or units. Implementation of EOS resolution and process optimization in an EO framework is motivated. The EO approach was implemented in the AMPL environment and IPOPT was used as the nonlinear solver. Formulation III was adopted for flash calculation. The advantage of EO approach is demonstrated through the polymerization process simulation and optimization. 4.1 Problem statement

Fig. 9. illustrates the flowsheet of an industrial HDPE slurry process with two CSTRs in series. The 24


Page 24 of 61

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


five-site Ziegler–Natta catalyst system and the other feed streams, including ethylene, hydrogen, and n-hexane, are fed continuously into the CSTRs. The product of the first reactor is fed into the second reactor. The vapor streams leaving the reactors are recycled to the feed streams through coolers, flash drums, and compressors to achieve high monomer conversion. The final product that leaves the second reactor is dried and pelletized. Fig. 9. Flowsheet of HDPE slurry process Table 5. Kinetic mechanism of homopolymerization. The kinetic mechanism is summarized in Table 5, where Cp denotes a potential active site of the catalyst, A is the co-catalyst, M is the monomer, H2 is hydrogen, P0 is an active site, Cd is a deactivated site, Dn is a dead polymer of chain length n, Pn is a living polymer chain of length n, j represents the index of active sites, and kaA, ki, kp, ktM, ktH, ktA, kt, and kd are the kinetic constants of different reactions. The pseudo kinetic rate constant of transfer and deactivation is defined as qn ()) = 9rs ())[t]F + 9rg ())[bf ]F + 9ru ())[A] + 9r ()) + 9 ()),

(13)

where [∙] denotes the concentration of the components and subscript L denotes the liquid phase. For living chains of length w = 1, the net reaction rate is

+\x(&) = 9 ())[t]F [!= ())] + 9rs ())[t]F y = ()) − 9 ())[t]F [!z ())] − qn ())[!z ())].

(14)

+\{(&) = 9 ())[t]F ([!|z ())] − [! ())]) − qn ())[! ())].

(15)

+nx (&) = qn ())[!z ())] − 9rs ())[t]F y= ()),

(16)

For living chains of length w ≥ 2, the net reaction rate is

Similarly, the net reaction rate for dead chains of length w = 1 and w ≥ 2 can be derived as

+n{ (&) = qn ())[! ())].

The mth moment of the living and dead polymers at the jth active site are, respectively,

25


(17)


∞

Page 26 of 61

∞

Y m ( j ) = ∑ n m [ Pn ( j )] , X m ( j ) = ∑ n m [ Dn ( j )] n =1

(18)

n =1

The net reaction rates of the zeroth-order moments of living and dead polymer chains are given by +} ~ (&) = 9 ())[t]F [!= ())] + 9rs ())[t]F y = ()) − qn ())y= ()), + ~ (&) = qn ())y = ()) − 9rs ())[t]F [!z ())].

(19) (20)

Similarly, for the first- and second-order moments of living and dead polymer chains, the net reaction rates are provided by +} x (&) = 9 ())[t]F [!= ())] + 9 ())[t]F y = ()) + 9rs ())[t]F y = ()) − qn ())yz ()), + x (&) = qn ())yz ()) − 9rs ())[t]F [!z ())],

+} (&) = 9 ())[t]F [!= ())] + 9rs ())[t]F y = ()) + 9 ())[t]F 2yz ()) + y = ()) − qn ())y f ()), + (&) = qn ())y f ()) − 9rs ())[t]F [!z ())].

(21) (22) (23) (24)

The number average molecular weight is calculated by ∑ (} x (&) x (&))

tw = ∑x (} ~ (&) ~ (&)) x

(25)

where Ns is the number of active sites. The number average degree of polymerization is DPN = tw/,

(26)

where mw is the molecular weight of monomer. The weight average molecular weight is calculated by ∑ (} (&) (&))

t, = ∑x

x x x(} (&) (&))

(27)

Obviously, Mn and Mw are determined by the concentration of moments. The calculation of moments requires a series of equations that describe the thermodynamics, phase equilibrium, mass balance and flowsheet connections. The detailed process model can be found in previous work [24]. The information flow of the process model is shown in Fig. 10. It presents strong nonlinearity with

26


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


multilayer coupling among polymerization kinetics, EOS resolution, and flowsheet calculation. From the point of view of model partitioning and sequential modular approach, thermodynamics calculation that playing the “service role” requires feedback from two different layers of outer loop (see the arrows in red). The simulator requires higher convergence accuracy for the inner calculation loop than the outer loop. However, in addition to the phase equilibrium in the outer loop, the inner loop for thermodynamic calculation also requires the polymer property calculation in the outermost layer to obtain the property variable for polymers. This kind of multi-layer dependence will lead to more computational cost and convergence difficulties. Moreover, the lack of exact derivative information may lead to divergence in each layer for derivative-based algorithms in optimization problems. The loops are solved repeatedly and any intermediate failure in each inner loop can be detrimental to the outer loop. Hence, EO approach is more suitable for this kind of complex flowsheets with multi-layer dependence and implicit design specifications. Figure. 10. Model structure of HDPE slurry process 4.2 Simulation of the polymerization process

To demonstrate the advantages of the EO approach for flowsheet simulation of the HDPE process shown in Fig. 10., several cases were studied based on industrial operating conditions. The operating conditions of two different polymer grades can be found in our previous work [25]. The EO and SM approaches are both tested for the flowsheet simulations. The simulation for flowsheet with large recycle loops is difficult for the SM approach and careful initialization is required. A common strategy is to tear the reflux stream and adjust the tear variables to initialize the flowsheet simulation. Two simulation cases based on the two polymer grades, namely Grade G and O were performed. Table 6 presents the simulation status of the two cases solved by the EO and SM approaches, respectively. 27



The first case is the simulation of Grade G. The initialization of the EO and SM approach were set to a point with low constraint violation of 8.66. The EO approach succeeded in the simulation in 7 iterations costing 0.023 CPU s. But the SM approach failed to converge with 200 iterations (loops) costing 234 CPU s. However, when the failed points are used as the initial guess of the EO approach, the simulation succeeded in 8 iterations costing 0.035 CPU s. The second case is the simulation of Grade O. The initial guess of the EO and SM approach were set to the EO converged point at Grade G with constraint violation of 1.77e2. The EO simulation succeeded in 6 iterations costing 0.019 CPU s. But the SM approach failed to converge after reaching the maximum 200 iterations (loops) costing 491 CPU s. However, when the failed points are used as the initial guess of the EO approach, the simulation succeeded in 14 iterations costing 0.055 CPU s. Table 6. Flowsheet simulation status The key process variables of the simulation results of Grade G were compared as shown in Table. 7. The SM approach in the third case predicted significantly large fugacity coefficients for light components in the liquid phase, so almost all the light components go into the vapor phase. The mole fraction of heavy component, such as HDPE, was as large as 98%. This phenomenon is inconsistent with the process mechanism. It seems that the SM approach finds a solution in the single phase region for light components. Note that the reflux stream does not affect the mass balance of the flowsheet, the liquid outlet mass flow of the reactor should be consistent with the feed flows. The inconsistent result of liquid outlet flow of the reactor BF in Table 7 indicates an imbalance of the tear streams. In other words, the failure of the SM approach in flowsheet simulation with reflux streams may be explained. The reflux streams are torn in the SM approach, and the inner thermodynamic calculation is repeated for every request for property value by an outer calculation loop. The PC-SAFT EOS is resolved when

28


Page 28 of 61

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


any of the intensive variables change. The intermediate phase equilibrium procedure affects the polymerization kinetics and thus changes the average chain length of the polymers (i.e. the DPN). Hence the key variable in the PC-SAFT EOS is changed and predicts a different property value. The intermediate procedures may traverse infeasible thermodynamic regions in the iterative loops, and the SM simulator fails to continue the outer loop iteration. Thus the tear streams cannot be balanced. Table 7. Flowsheet simulation case study 4.3 Optimization of the polymerization process

To achieve different end-use polymer specifications, a polymerization process generally produces products with various grades through the same equipment. Industrially, the designed operating conditions for a certain grade are often derived by experience, which may not result in optimal economic performance. In this section, the optimization maximizes the production rate of HDPE constrained by the process and polymer quality specifications. The optimization problem is formulated as follows: max Bgn\o

s. t. t, ≤ t, ≤ t, F

twF ≤ tw ≤ tw F ≤ z ≤

g(t,, tw, z) = 0

(28)

where Bgn\o is the flow rate of the HDPE product, and z represents the decision variables, including

temperature and pressure of the reactor, the feed flow rates of hydrogen and ethylene. g(t,, tw, z) =

0 represents the complete EO model for this polymerization process, including the PC-SAFT EOS with Formulation III and the process model.

Three cases differing in the Mw and Mn specifications are presented as shown in Table. 8. The operating conditions of Grade G are used as the initial guess for all the three cases. All three cases 29



could be reliably solved using IPOPT. The optimization problems were solved in about 50 iterations with less than 0.5 CPU s. This indicates that the EO framework is robust since the solver can find optimal solutions having three different polymer quality specifications from the same starting point. We also attempted to test these cases in the SM approach using Aspen Plus, but unfortunately all the cases failed. In the SM mode, SQP approach is adopted for optimization after the inner simulation converged. Any intermediate failure of convergence in the each inner loop is detrimental to the optimization. In these cases, the inner simulations of the reactor diverged after 200 iterations using Wegstein method and the optimizer returned bad results after 3-6 loops. The EO formulation for PC-SAFT EOS exhibits its benefits in process optimization and motivates more potential applications in process modelling. Table. 8. Three cases of flowsheet optimization

5. Conclusions We present an EO framework that combines resolution of PC-SAFT EOS and process optimization. First the mathematical structure of the PC-SAFT EOS is analyzed through digraph of variables. The digraph reveals that for polymerization systems, the resolution of the PC-SAFT EOS is not only coupled with the phase equilibrium, but also interacts the polymerization kinetics. In order to improve the computational efficiency, we implement both EOS resolution and process optimization in a complete EO framework. Then a modified EO formulation for PC-SAFT EOS is proposed to select the appropriate roots in the flash calculation. It is shown that the proposed EO formulation always select the appropriate root through simulation and optimization of a flash drum. The results are validated by the conventional resolution approach in commercial software. Finally, the modified EO formulation for PC-SAFT EOS is applied in the simulation and optimization of an industrial polymerization process. The results demonstrate advantages of the proposed EO method compared to traditional SM-based 30


Page 30 of 61

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


commercial software in terms of computational robustness and efficiency. Appendix A: Presentation of the original PC-SAFT model

Essential equations of the PC-SAFT EOS that must be simultaneously solved: Temperature-dependent segment diameter: () = %( 1 − 0.12 exp −

3$(

9

, ( ∈ [1, wp ]

Mixing rules: % & =

1 % + %& , {(, )} ∈ [1, w ]f 2

$ & $ $& = ∙ (1 − 9 & ), {(, )} ∈ [1, w ]f 9 9 9 Mean segment number:

= *

z

Definition of coefficients:

f $% < = * * (

& &

z &z

$ & < )% 9 &

$ f $ f % < = * * ( & )f % <

& & 9 &

= = + ' = '= +

z &z

−2 −1 − 1 z + f , ( ∈ [0, 6] −1 − 1 −2 'z + 'f , ( ∈ [0, 6]

Number density of molecules:

6 7 = ( * < )|z

z

Definition of coefficients, note that the reduced density is defined by = <

= 7 * , w ∈ [0, 3] 6

z

31



Definition of power series: ¡

z =

= ¡

f = '

=

Definition of coefficient: ?z = (1 +

8 − 2f 20 − 27f + 12< − 2h |z + (1 − ) ) (1 − )h ((1 − )(2 − ))f

Radial distribution function of the hard-sphere fluid:

£ & =

& =

& , {(, )} ∈ [1, w ]f + &

1 3f 2ff f + & + , {(, )} ∈ [1, w ]f

& 1 − < (1 − < )f (1 − < )

Equation-Oriented Approach for Handling the Perturbed-Chain SAFT

Recommend Documents