Comparative Evaluation of Rigorous Thermodynamic Models for the

Dec 4, 2015 - *E-mail: [email protected]. Tel. ... Taking a closer look SRK-MHV2 shows slightly better correlative capabilities with lowe...
1 downloads 9 Views 2MB Size
Article pubs.acs.org/IECR

Comparative Evaluation of Rigorous Thermodynamic Models for the Description of the Hydroformylation of 1‑Dodecene in a Thermomorphic Solvent System Victor Alejandro Merchan* and Günter Wozny Chair of Process Dynamics and Operations, Technische Universität Berlin, 10623 Berlin, Germany S Supporting Information *

ABSTRACT: In this contribution implications of using different rigorous thermodynamic models for the description of the liquid−liquid equilibrium (LLE) appearing in an hydroformylation process of 1-dodecene in a decane/N,N-dimethylformamide thermomorphic solvent system are discussed. Besides a fully predictive GE model (UNIFAC Dortmund), a representative of cubic equations of state with mixing rules based on GE-models (SRK-MHV2 with nonrandom two liquid (NRTL)) and a representative of equations of state based on statistical associating fluid theory (heterosegmented PCP-SAFT) are evaluated. The influence of the chosen thermodynamic model on the simulated process-wide performance is quantified, and the simulation results are compared against published experimental data from a miniplant with full recycle. Special attention is given to equations of state approaches with fitted binary interaction parameters, since they allow a good correlation of available experimental data and are applicable to process-wide simulations under consideration of gaseous compounds. A short analysis regarding the computational effort related to the evaluation of the equation of state models is additionally provided. The obtained results show that both selected equation of state models, in conjunction with the parameters obtained in this contribution, have a good agreement with available phase equilibrium and miniplant experimental data. Hence, both approaches appear to be suitable for process design and optimization tasks. Taking a closer look SRK-MHV2 shows slightly better correlative capabilities with lower computational effort, although requiring a larger amount of binary interaction parameters. However, while phase equilibrium compositions are almost equally well described by both approaches, particular care should be taken when calculating volumetric properties. For the latter properties the heterosegmented PCP-SAFT approach provides better results.

1. INTRODUCTION Within the field of homogeneous catalysis, and in particular for groups working on tunable solvent systems, the hydroformylation of long chain olefins is a highly interesting model reaction. On the basis of the simultaneous goals of high product selectivity and efficient catalyst recycle, several process approaches have been developed, among others, the hydroformylation in thermomorphic multicomponent solvent (TMS) systems,1−3 micellar solvent systems (MSS),4,5 ionic liquids (IL),6,7 CO2-expanded liquids (CXL), and gas-expanded liquids (GXL).8,9 Toward the development of novel processes with innovative solvent systems there is a great demand for the application of advanced model-based methods. These offer a systematic way of finding process improvements while reducing the effort for complex and expensive experiments. Since the impact of such methods is strongly linked to the existence and quality of high fidelity process models, a major condition for the success of model-based methods for processes related to tunable solvent systems is the availability of proper thermodynamic models. The computational costs to properly describe the phase behavior with a selected thermodynamic model should be reasonable. In contrast to other tunable solvent systems, for example MSS,10 the TMS system consisting of N,N-dimethylformamide (DMF) and decane, as used among others in refs 3, 11, and 12 for the hydroformylation of 1-dodecene, represents a reacting system for which a consistent thermodynamic description of © XXXX American Chemical Society

phase behavior using standard engineering modeling approaches can be expected. Though the modeling and simulation of such gas−liquid, size-asymmetric reacting systems containing highly polar and nonpolar components at a relatively high pressure level of 2 MPa still represent a major challenge, it is evident that the aforementioned hydroformylation TMS system represents a promising system for the application of modelbased methods. Consequently some modeling efforts have been performed. Up to now these have either focused on the description of phase behavior with advanced engineering equations of state (EoS) models (see ref 13 for gas−liquid equilibrium (GLE) and ref 14 for liquid−liquid equilibrium (LLE)), on reactor related issues such as optimal identification of kinetic model structure and parameters,11,15 and on optimal reaction routes.15,16 While the methodology discussed by Hentschel et al. in ref 16 clearly shows the wide field of applications of model-based methods toward optimal process-wide decisions, it has not been discussed in detail whether the used models and the obtained results reflect the real process behavior accordingly. The simple fact that no temperature dependence is considered for the LLE of a TMS system illustrates this problem clearly. Before the optimal results obtained with Received: September 7, 2015 Revised: December 3, 2015 Accepted: December 4, 2015

A

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

Article

Industrial & Engineering Chemistry Research

and hence catalyst rich phase resulting from the phase separation is recycled to the reactor, while the organic phase (product phase) can be further purified in subsequent downstream steps. Table 1 lists the main components that appear in the considered process. Following the suggestion from Markert et

approximated LLE calculations are analyzed, it may be helpful to find out how the quality of the LLE description may affect simulation results concerning an overall hydroformylation process. Constrained by the available LLE data, in this contribution we study the implications of using particular LLE modeling approaches on the steady state behavior of a hydroformylation process as described in ref 12. Besides technical aspects such as process-wide aldehyde yield at different process conditions (recycle rates, decanter temperature, reactor outlet concentration) computational aspects regarding equation-oriented solution are discussed. This contribution is organized as follows: Section 2 gives a short introduction of the hydroformylation process, an overview of available experimental equilibrium data and current modeling approaches, and discusses the model assumptions that are applied within the process-wide analysis and validation. Because the current description of LLE with aldehyde systems does not show a satisfactory accordance with experimental data, section 3 shows alternative modeling approaches that allow a better accordance. After introducing the model implementation strategy based on algorithmic differentiation, alternative parametrizations for the heterosegmented PCP-SAFT model and an EoS/GE approach17 are discussed. Section 4 presents the results of this contribution. Besides the results of the examined parametrization cases, results of simulation studies are presented that clearly show the influence of the chosen thermodynamic model on the process-wide simulation results. At this point also a validation against experimental miniplant data is presented and computational issues related to the equation of state models are discussed. Finally section 5 summarizes the main results and draws conclusions regarding the applicability/advantages of the different thermodynamic models for the particular hydroformylation process.

Table 1. Considered Components Appearing in the Hydroformylation of 1-Dodecene in Decane/DMF TMS System component number

component name

component short name*

1 2 3 4 5 6 7 8 9 10 11 12

1-dodecene n-tridecanal iso-dodecene iso-tridecanal dodecane N,N-dimethylformamide n-decane hydrogen carbon monoxide active catalyst cat precursor Rh(acac) (CO)2 BIPHEPHOS

nC12en nC13al iC12en iC13al nC12an DMF nC10an H2 CO cat prec

*

Though these short names appear as subscripts or superscripts they are not listed separately in the nomenclature.

al.19 and according to the reaction kinetics derived by Kiedorf et al.11 all appearing isomeric alkenes and aldehydes were respectively lumped into single pseudocomponents designated as iso-dodecene and iso-aldehyde. Dodecane results from the hydration of 1-dodecene, which is also taken into account in the reaction network. Though the catalyst and ligand concentration have a decisive influence on the reaction performance, the components 10−12 are neglected in the thermodynamic description of the system. On the basis of the very low amounts of catalyst complex in the mixture, it is assumed that these components can be neglected without effecting the phase equilibrium calculations. Their molar fractions differ from the main components molar fractions by 3 to 4 orders of magnitude. 2.1. Overview of Available Experimental Data and Models. 2.1.1. Liquid−Liquid-Equilibrium/Gas−Liquid-Equilibrium. Table 2 gives an overview about the published LLE

2. CONSIDERED PROCESS AND AVAILABLE INFORMATION The hydroformylation process with a BIPHEPHOS modified rhodium catalyst, which is mainly soluble in polar solvents, is summarized in Figure 1. This process has been presented in

Table 2. Considered LLE Experimental Data

Figure 1. Simplified flowsheet of hydroformylation process with stream numeration for process modeling. Dotted line makes clear that no decanter gas output is taken into account in the mathematical model.

detail in.12,18 In addition to the reactants 1-dodecene and syngas, the apolar solvent n-decane and the polar solvent DMF are required in order to conduct the process according to the TMS approach: The multicomponent liquid system, which is completed by the resulting products, behaves in such a way that the reacting liquid phase is homogeneous at reaction temperature (around 90 °C), while showing a phase split at lower temperatures (normally between 5 and 25 °C). The polar

mixture type

components

ref

temp[K]

points givena

binary binary ternary

DMF/decane DMF/1-dodecene DMF/decane/ dodecene

3b 3 3

283.15−347.65 288.15−333.15 298.15

15 6 6

3

333.15 343.15 283.15

7 10 8

20

288.15 298.15 298.15

6 9 12

ternary

ternary

DMF/decane/ n-dodecanal

DMF/decane/ n-tridecanal

a Points given per phase. bOnly measurements made with the analytical method were used.

B

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

Article

Industrial & Engineering Chemistry Research

Figure 2. Liquid−liquid equilibrium experimental data and modeling results of n-dodecanal/DMF/decane system. (Left) At 283.15 K; (center) at 288.15 K; (right) at 298.15 K. Symbols represent either experimental data from ref 3 (black circles) or LLE calculations with the respective models where the feed lies in the center of the experimental tie lines (see Appendix B.1). Lines are binodal curves calculated using the respective model.

Table 3. Inventory of Available and Missing Experimental Data Regarding Single Binary Pairs

a

component

i/j

1

1-dodecene n-tridecanala iso-dodecene iso-tridecanal dodecane DMF decane H2 CO

1 2 3 4 5 6 7 8 9

× × × × √ ○ √ √

2

3

4

5

6

7

8

9

×

× ×

× × ×

× × × ×

√ ○ × × √

○ ○ × × × √

√ √ × × √ √ √

√ √ × × √ √ √ √

× × × ○ ○ √ √

Available n-dodecanal is considered as n-tridecanal. at all.



× × × × × ×

× × × × ×

Binary GLE or LLE data available.

√ × √ √

√ √ √



√ √

√ ×

Ternary LLE data available. No information available

approach or using an alternative model that gives a better fitting while allowing application for multicomponent systems. A detailed study on gas solubility in the considered TMS system is given in a series of papers by Vogelpohl et al.13,24 Besides experimental work, the gas solubility was modeled using the standard PCP-SAFT EoS. It was shown that single, temperature independent binary parameters gained from binary experiments (one BIP per gas/liquid binary pair) are sufficient to correctly predict the temperature-dependent GLE of multicomponent systems, consisting of up to four liquid components and two gases, with good accuracy. With the exception of measurements of CO solubility in a n/ iso-tridecanal mixture, n-dodecanal was treated as the product aldehyde. Similar to the LLE case the GLE of CO, H2 and syngas in isomeric liquid components have not been separately studied. 2.1.2. Missing Information. As pointed out in section 2.1.1, there is the general problem that experimental data for systems containing n-tridecanal are very scarce. Fortunately the few available experimental LLE and GLE data show that ndodecanal and n-tridecanal have, even quantitatively speaking, a highly similar phase equilibrium behavior with the remaining compounds. Hence n-dodecanal binary measurements can be considered equal to the n-tridecanal ones. On the basis of this and taking into account the main compounds described in Table 1, Table 3 summarizes the available data from which BIPs are available or can be derived. From Table 3 it can be seen that for the compounds appearing in relatively high amounts in the liquid phase (decane, DMF, n-tridecanal, 1-dodecene), almost all binary

measurements taken into account in this contribution. All experimental data come from the work of Schäfer and Sadowski,3,14,20 who not only measured the LLE but also performed modeling work with the PCP-SAFT EoS21,22 model. While all considered binary systems and the first ternary system of Table 2 are very well described using the standard homosegmented PCP-SAFT approach with a relatively low amount of estimated binary interaction parameters, the description of systems containing the long chain aldehydes required the application of the copolymer PC-SAFT concept,23 also known as the heterosegmented PC-SAFT approach, to obtain a better accordance between experimental and simulation results. The resulting equation of state model is specified as heterosegmented PCP-SAFT (hs-PCP-SAFT). Though the modeling approach based on hs-PCP-SAFT introduced in3 brought a significant improvement in correlating the measurement data in comparison to the homosegmented approach, a detailed comparison of measurements and simulation data (see Figure 2) shows that the heterosegmented approach along with the proposed binary interaction parameters (BIPs) may not be a proper model selection for a rigorous plant-wide model. The major influence of the aldehyde molar fraction on the phase behavior is not described properly: The simulated mixing gaps are in general larger than those in the experiments, so that an acceptable fitting is only given at low aldehyde molar fractions. Besides that, the obtained information on the temperature dependence is of poor quality. The relative large deviations between experiments and simulation in systems containing aldehyde suggest the need of either finding a better parametrization for the hs-PCP-SAFT C

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

Article

Industrial & Engineering Chemistry Research

conceivable approach for treating the binary interactions with isocompounds would be using the same BIPs as already used for their linear counterparts. Further on, taking into account that the treated linear compounds and their corresponding, branched counterparts have very similar volumetric properties and taking into account that the considered mixing rules (see section 3.1) do not exhibit the Michelsen−Kistenmacher syndrome,26 an intuitive simplification would be to lump the corresponding n- and isocompounds into a single general alkene and a single general aldehyde. Neglecting dodecane and syngas, the decanter separation could be simplified to four instead of six compounds. While lumping n- and isocompounds prevents the possibility of quantifying the alkene isomerization and hydroformylation of alkane isomers in a precise way, the overall equilibrium compositions are not directly affected by this. In fact, the obtained molar fractions of the general alkene and aldehyde can be thought to consist of arbitrary combinations of a linear and a branched compound. Dodecanal/Tridecanal Issue. According to the modeling goal it is important to describe the temperature dependency of the LLE taking place in the decanter, hence in a system containing n-tridecanal. However, there is no temperaturedependent ternary LLE data for those systems, but for systems with n-dodecanal. Since the few available experimental data points of ternary n-tridecanal systems (at one single temperature) show a similar phase equilibrium behavior to the corresponding n-dodecanal system, ternary n-tridecanal systems will be considered to have the same (temperature-dependent) phase behavior as systems with n-dodecanal. When modeling LLE with equations of state additional care needs to be taken, since phase equilibrium calculations not only depend on BIPs but also on pure compound parameter values. In this particular case, however, simulation studies show that there are only negligible differences in multicomponent LLE calculations with either n-dodecanal or n-tridecanal using their respective pure component parameters with BIPs obtained from n-dodecanal ternary systems. Because of this and in order to further simplify the LLE evaluations, for the simulation studies in section 4.3 n-dodecanal is considered directly. Though replacing n-tridecanal with n-dodecanal would be a major modeling mistake in a first principle model, (the overall mass balance would be violated), for our modeling goal, it is an acceptable assumption. Laying the focus on the quantification of the influence of the decanter phase separation conditions on the whole process, as done in this contribution, it is irrelevant whether n-dodecanal or n-tridecanal is taken into account, since both are considered to have the same phase behavior. 2.3. Simplified Model of Reactor−Separator System. A standard flash formulation as given in Appendix B is required to realize the first modeling goal, that is, the phase equilibrium calculation of ternary d systems that contain n-dodecanal. Regarding the coupled reactor/separator system (Figure 1) a model formulation consisting of submodels for the mixer, the reactor, and the decanter (all models without energy balance) is studied according to the numbering of the stream given in Figure 1. Since for our objectives neither the heat exchanger nor the recycle pump need to be considered, the equalities of respective flows and compositions (Fs=6 = Fs=8, xs=6,c = xs=8,c, Fs=10 = Fs=11 and xs=10,c = xs=11,c) follow implicitly. To clarify some of the assumptions on the model formulation, the mathematical models for the reactor and decanter are shown next. Analogously the model of the mixer,

system information is given, only the interaction 1-dodecene/ntridecanal is missing. However, for some other pairs, in particular those containing iso-components, no information is given at all. Following a rigorous approach, however, these influences need to be evaluated as well, either through an estimation or a predictive approach. Regarding the complex multicomponent phase equilibria at high pressure (2 MPa) that take place in the decanter, no experimental data is available. The thermodynamic model can only be indirectly validated with the experimental data from the miniplant with full recycle as described in ref 18. Such a comparison is provided in section 4.3. 2.2. Modeling Approach/Modeling Goal. At a first stage the modeling goal is to obtain a correlation of the available ndodecanal−DMF−decane LLE data that is better than the state of the art.14 This should be done either by using further parametrizations of the hs-PCP-SAFT EoS or alternative, rigorous thermodynamic models (see section 3.1). Further on, holding on to rigorous thermodynamic models and being constrained by the missing information discussed in section 2.1.2 and further assumptions discussed in section 2.2.1, a mathematical model for the coupled multicomponent reaction/separation process with recycle (Figure 1) is built. The model allows not only a general analysis of the process behavior, but also to some extent a comparison with the experimental data presented in ref 18. It is of particular interest to quantify the influence of the solvent feed distribution, feedrecycle ratios, and separation temperature in the decanter using different rigorous thermodynamic models. The following simulation cases (modeling goals) should be supported by the model: • Find reactor performance/alkene conversion XR (see eq 10), or process performance/process-wide aldehyde yield Yprocess C13al (see eq 11) for fixed decanter feed. This way it would be possible to study the process behavior with reactor conditions as obtained from batch experiments.11 • Find XR or Yprocess C13al that guarantees stationary operation under fixed process feed conditions, fixed decanter conditions, and fixed recycle flow rate. 2.2.1. Model Assumptions. Besides neglecting the catalyst complex as justified at the end of section 2, the following assumptions are considered: Neglecting Gas Phase in Decanter Separation. Though it is already well-known that gas solubility and gas distribution have a major role on the reactor performance and product distribution,11 for the discussed modeling goals and simulation studies a detailed examination of these factors is not necessary. Simulation studies with rigorous models not only show that the molar fractions of gas in the liquid phases have very low values, but also make clear that a syngas feed has a negligible influence on the resulting liquid phases and product distribution leaving the decanter. Neglecting Dodecane. Neglecting the hydrogenation reaction and hence the presence of dodecane is best justified by the experimental fact that only very small dodecane amounts are produced. In case of slightly higher amounts neglecting dodecane amounts could additionally be justified by the fact that binary information is only available for mixtures with DMF.25 Simplifying Isocompounds. Taking into account that the isocompounds already represent lumped compounds and that there is no information on binary interactions at all, a D

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

Article

Industrial & Engineering Chemistry Research

yield of linear aldehyde and accordingly xs=9,iC13al ≪ xs=9,nC13al. Hence, the yield calculated according to eq 11 can be compared with the yield toward linear aldehydes that was obtained experimentally. A comparable statement is not possible for the conversion of dodecene, since there is not such a large difference in the absolute values of the real molar fractions xs=9,nC12en and xs=9,iC12en. The model equations were implemented in the web-based modeling environment MOSAIC.27 For simulation cases considering the UNIFAC−Dortmund and the SRK-MHV2 models code generation with external function calls for gPROMS28 was used.29 The thermodynamic packages required for function calls were generated with Aspen Plus V8.2 as CAPE-OPEN property packages. For the simulation cases regarding hs-PCP-SAFT and for some SRK-MHV2 evaluations code generation for MATLAB was used. In this case the fugacity coefficient calls were embedded as external functions (see section 3.2).

which also needs to be considered, consists of component balances for each compound and summation equations. For a better readability the compounds are organized in sets, LIQ = {C12en,C13al,DMF,C10an} for liquid compounds and GAS = {H2,CO} for the gases. Additionally COM = LIQ ∪ GAS. C12en and C13al, represent a general olefin and aldehyde with the pure compound properties of 1-dodecene and ndodecanal, respectively. Reactor with Full Gas Separation. 5

7

∑ Fs ,c·xs ,c − ∑ Fs ,c·xs ,c + νc·r = 0,

c ∈ COM

s=4

s=6

(1)

xs = 7, c = 0,

c ∈ LIQ

(2)

xs = 6, c = 0,

c ∈ GAS

(3)



xs , c = 1,

s = 4, 5, 6, 7 (4)

c ∈ COM

3. PHASE EQUILIBRIUM MODELING Though the modeling goal only requires a proper model for the description of the LLE in the decanter, keeping in mind the option of a later plantwide simulation under consideration of all compounds, it is important to consider models that not only allow a good description of the LLE (e.g., a GE-model), but also allow the consideration of supercritical components. For the description of such complex multiphase equilibria, models based on association theories,30 such as the cubic plus association equation of State (CPA EoS)31 or the perturbed chain statistical associating fluid theory (PC-SAFT)21 should perform well.13,14 Another alternative is given by the use of socalled EoS/GE models, which combine a cubic EoS model with a mixing rule based on a model for the excess Gibbs energy.17 Though these models in comparison with advanced association models normally require a larger number of binary interaction parameters, they also provide in many cases excellent correlative capabilities. An advantage is given by the fact that this latter type of model is currently available in most commercial process simulators. 3.1. Considered Models. 3.1.1. PCP-SAFT EoS. Being a modification of the statistical associating fluid theory (SAFT) by Chapman et al.,32 PC-SAFT, developed by Gross and Sadowski,21 is a particularly successful equation of state model, which has shown impressive capabilities in modeling thermophysical properties of pure components and phase equilibria of multicomponent systems. PC-SAFT provides an expression for the residual Helmholtz energy Ares which consists of several additive contributions accounting for different intermolecular forces. Besides the hard chain contribution (hc) that accounts for repulsive interactions, relevant attractive forces such as London dispersion (disp), dipole−dipole interaction (dipole) and association (assoc) are usually considered (eq 12). The use of different expressions for the dipole−dipole contribution Adipole leads to different denominations of the PC-SAFT approach. Using the dipole interaction term proposed by Gross and Vrabec22 leads to the so-called perturbed chain polar SAFT (PCP-SAFT) EoS.

Decanter Unit. 10

Fs = 6·xs = 6, c −

∑ Fs·xs ,c = 0,

c ∈ COM (5)

s=9

φs = 9, c ·xs = 9, c = φs = 10, c ·xs = 10, c ,

c ∈ LIQ

(6)

xs = 9,H2 = 0

(7)

xs = 10,CO = 0

(8)



xs , c = 1,

s = 9, 10 (9)

cinCOM

Reactor Conversion/Process Aldehyde Yield. XR =

Fs = 4 ·xs = 4,C12en − Fs = 6·xs = 6,C12en

process YC13al =

Fs = 4 ·xs = 4,C12en

(10)

Fs = 9·xs = 9,C13al Fs = 1·xs = 1,C12en

(11)

Following the model assumptions listed in section 2.2.1 regarding the neglected hydrogenation reaction and the lumped alkene and aldehyde compounds there is only the need to take into account one reaction rate r. The fugacity coefficients φs,c are evaluated according to eq 17. It should be remarked that for the simulation cases listed in section 2.2 the assumption of full syngas separation through stream 7 (see eqs 2 and 3) can lead to a simpler model formulation. Since in all considered cases the reaction rate r results directly from the component balance of aldehyde (C13al) around the reactor and thus determines the gas consumption in the reactor, it can be shown that the syngas feed has no influence on the reactor or decanter performance. Hence it is viable to completely neglect syngas in our modeling approach. In fact, the model implementation used for the simulation cases discussed in section 4.3 neglects syngas. Regarding the reactor and process performance indicators XR process and Yprocess C13al it should be remarked that the YC13al is particularly valuable for the validation against experimental data. For the process under consideration (in particular for the ligand under consideration) experimental work11,15,18 shows that the yield of isomeric aldehydes is negligibly small in comparison with the

Ares = Ahc + Adisp + Adipole + Aassoc

(12)

In refs 3 and 13 it was found that for the description of the considered hydroformylation mixtures (see Table 1) dipole− dipole interactions need to be considered, but no association forces. On this basis, the resulting expression for the reduced E

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

Article

Industrial & Engineering Chemistry Research residual Helmholtz energy ãres = Ares/(N·k·T)) is given by eq 13. hc disp dipole a ̃res(ρmix , xmix ⃗ , T) = ã + ã + ã

within this contribution, we imply NRTL as the underlying activity coefficient model. In addition to the already referred hs-PCP-SAFT and SRKMHV2 models, for some LLE evaluations the UNIFAC− Dortmund (UNIFAC-DO)36 model was used. Among the UNIFAC variants available through Aspen Plus V8.2, UNIFAC−Dortmund fits the available experimental data the best. Flash calculations with UNIFAC−Dortmund are done with Aspen Plus V8.2 or using a CAPE-OPEN property package generated by it. 3.2. Model Implementation and Evaluation Methods. For the comparison of the correlative capabilities of the different models, single flash calculations at fixed inlet concentrations, pressure, and temperature, also called pTflash calculations, are evaluated. Both models, SRK-MHV2 and hs-PCP-SAFT, were coded in MATLAB release 2013a (The MathWorks, Inc.) as functions for the evaluation of the reduced residual Helmholtz energy according to eq 13 and eq 25. The partial derivatives of ãres with respect to ρmix or vmix and with respect to xmix,i, which are required for the evaluation of the compressibility factor Zmix (eq 16) and the fugacity coefficient φmix,i (eq 17), respectively, are obtained by the automatic differentiation tool ADiMat.37 For the sake of a simpler notation, the subscript “mix” introduced in eq 13 is neglected in the following equations.

(13)

The notation in ãres(ρmix,xm ⃗ ix,T) indicates the computational implementation of ãres as a function of the total number density of molecules of the mixture ρmix, the component molar fractions of the mixture xm ⃗ ix and the temperature T (see section 3.2). To evaluate the expressions given in eq 13 according to the hsPCP-SAFT approach, for each segment α of an apolar component i, three pure component parameters are needed, the number of segments per chain mi,α, the segment diameter σi,α, and the segment dispersion energy ϵi,α/k. In the case of polar components, additionally to the three previously mentioned parameters, the dipole moment μi,α is required. For the calculation of mixtures, mixing rules are required. As usual the Lorentz−Berthelot mixing rules are applied, whereas in this case they refer to the interaction of a segment α of component i with a segment β of component j (eqs 14 and 15). It should be remarked that with exception of the long chain aldehydes, which consists of two segments, all other compounds are considered to be homosegmented. Further only heterosegmented components with one polar segment are considered. σi , α , j , β =

1 ·(σi , α + σj , β) 2

ϵi , α , j , β = (1 − ki , α , j , β) · ϵi , α /k·ϵj , β /k

⎛ ∂a ̃res ⎞ ⎛ ∂a ̃res ⎞ ⎟ Z = 1 + ρ·⎜ ⎟ = 1 − v·⎜ ⎝ ∂v ⎠ ⎝ ∂ρ ⎠

(14) (15)

⎛ ∂a ̃res ⎞ ln φi = a ̃res + ⎜ ⎟− ⎝ ∂xi ⎠

The model equations building the hs-PCP-SAFT approach, as considered in this work, can be found in refs 33 and 14. Minor typing errors appearing in refs 33 and 14 can be located by proofing their consistency with the polar term published in ref 22. Sauer et al.34 provide alternative formulations of the heterosegmented dipole contribution within the context of group contribution methods. In this contribution, one of these additional approaches is considered as well. 3.1.2. SRK-MHV2. An alternative, thermodynamically consistent method able of correlating all available experimental data while taking into account possible pressure effects is given by so-called EoS/GE models.17 In this contribution, an approximate zero reference pressure model, the MHV2 mixing rules,35 comes into consideration. Well known drawbacks of such zero reference pressure models, like low performance in description of size asymmetric systems containing gases with hydrocarbons of different lengths, are often related to the fact that mixing rules very often use group contribution methods as its underlying activity coefficient model. Such problems also appeared in preliminary studies when EoS/GE models with different mixing rules and different UNIFAC variants were evaluated. In those cases a large deviation between model and measurements not only appeared in mixtures of gases with long hydrocarbons but also in ternary systems with n-dodecanal, DMF, and decane. To avoid such deviations and to get a good correlation of all available experimental data, in this study the SRK EoS with MHV2 mixing rules and nonrandom two-liquid (NRTL) as underlying activity coefficient model is used. With this model it is also possible to obtain a satisfactory fit of the GLE measurements available. The considered model equations describing the SRK-MHV2 model are given in Appendix A. Please note that whenever we refer to the SRK-MHV2 model

NC

(16)

⎛ ∂a res ⎞ ̃ ⎟ ⎟ + Z − 1 − ln Z ∂ x ⎝ j ⎠

∑ xj·⎜⎜ j=1

(17)

For the LLE pT-flash calculations a standard equilibrium stage model as formulated in Appendix B was implemented. The LLE calculations were initialized with the direct substitution method38 and switched to a simultaneous Newton based method after a few iterations. The latter solution is carried out with MATLAB’s “fsolve” using the trust region dogleg algorithm with default settings. 3.3. Parameter Estimation. While for the PCP-SAFT approaches it is only necessary to find a better description of the n-dodecanal−DMF−decane system, in the case of the SRKMHV2 model it is necessary to estimate parameters that allow a proper correlation of all available measurements listed in Table 2. 3.3.1. hs-PCP-SAFT. A possible way of obtaining a better description of the ternary n-dodecanal−DMF−decane system is given by simultaneous parameter estimation against ternary LLE-data and pure compound properties of n-dodecanal. The chosen structure of the heterosegmented n-dodecanal molecule plays a major role in defining the degrees of freedom of the parameter estimation problem. What can be estimated is not only the pure compound parameters of the n-dodecanal segments but also binary interaction parameters between the n-dodecanal segments, DMF, and decane. It can be argued that simultaneously estimating pure compound properties and LLE data can lead to arbitrary parameter values, since a relatively high number of free parameters need to be estimated. This is certainly true, but following the main goal of correlating the available LLE data in a better way and under consideration of missing binary data for most binary subsystems, it appears to be an acceptable option. F

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

Article

Industrial & Engineering Chemistry Research

sets were successfully checked against stability analysis. For this latter point, the formulation proposed in ref 41 was considered, but only local optimization methods were used. 3.3.2. SRK-MHV2. In addition to the ternary system with ndodecanal, DMF, and decane, for the description of the required LLEs with the SRK-MHV2 model, it is necessary to estimate BIPs for the remaining binary and ternary systems listed in Table 2. As common for EoS/GE models the BIPs to be estimated correspond to the ones of the activity coefficient model (see eqs 36 and 37). As shown in Table 3 the BIPs for the DMF/decane and 1-dodecene/DMF pairs can be directly estimated from binary measurements. Additional interaction parameters can be determined from the ternary measurements. From the ternary 1-dodecene-DMF-decane system BIPs for dodecene/decane binary pair result. From the ternary ndodecanal−DMF−decane system BIPs can be obtained for the n-dodecanal/decane and n-dodecanal/decane pairs. According to eqs 36 and 37 there is a relatively large number of parameters that can be fitted. For the most interesting ndodecanal−DMF−decane subsystem in section 4.2 several cases considering a different amount of parameters are compared. The BIPs were estimated using the Aspen Plus data regression system42 with the Maximum-Likelihood objective function. The strategy used for the parameter estimation is discussed in more detail in the Supporting Information. Additionally, parameter estimation studies with MATLAB using eq 20 as objective functions were carried out. The results obtained with Aspen were confirmed.

Furthermore, by allowing the appearance of additional parameters, the correlative strength of hs-PCP-SAFT, or its model structure, can be evaluated. The already mentioned problem of unrealistic parametric values can be partially handled by using parameter estimation algorithms which allow the definition of variable bounds and setting reasonable variable bounds. A safer approach would be the application of methods discussed in ref 39. The parameter estimation problem is solved through the minimization of a least-squares objective function OF, which consists of two different parts (eq 18), one for the pure compound properties OFpure (eq 19) and one for LLE OFLLE (eq 20). The equilibrium compositions xorg,cal , xpol,cal and the n,i n,i pure compound properties saturated liquid density ρSat,cal and L,n vapor pressure pLV,cal, are obtained from the solution of the pTflash equation system provided in Appendix B.1 and Appendix B.2, respectively. OF = OFpure + OFLLE NPM

OF

pure

=

∑ n=1

(18)

2 ⎛⎛ LV,exp ⎛ ρ Sat,exp − ρ Sat,cal ⎞2 ⎞ − pnLV,cal ⎞ ⎜⎜ pn L,n ⎟⎟ ⎟ + ⎜ L,n LV,exp Sat,exp ⎟ ⎜⎜⎜ ⎜ ⎟ ⎟⎟ pn ρL , n ⎠ ⎝ ⎠⎠ ⎝⎝

(19)

⎛⎛ org,exp ⎞2 − xnorg,cal xn , i ,i ⎜ ⎜ ⎟ ∑ ⎜⎜ org,exp ⎟ x n , i ⎝ ⎠ i=1 ⎝

NLLE NC

OF

LLE

=

∑ n=1

⎛ x pol,exp − x pol,cal ⎞2 ⎞ n,i n,i ⎟⎟ + ⎜⎜ ⎟⎟ pol,exp x ⎝ ⎠⎠ n,i

4. RESULTS AND DISCUSSION 4.1. Alternative PCP-SAFT Parametrization for Ternary Systems Containing Long Chain Aldehydes. Table 5 lists the resulting pure compound parameters along with their respective average relative deviations (ARD) and maximum relative deviation (RD). While the homosegmented approach clearly provides the best fit, it should be remarked that further approaches show reasonable fit with at least acceptable deviations. Also the parameter values of most cases lie in an acceptable range. Only the results for the parameter ϵ/k in case E present unusual values, situated on the parameter bounds defined in the optimization problem. Table 6 lists the corresponding BIPs. Their respective correlative capabilities are evaluated in Figure 3, which shows the average relative deviation of the molar fractions of each component in each phase on a percentage basis (ARD[%]). To allow a comparison with the correlations obtained for further binary and ternary systems (see section 4.2) the ARD of the molar fractions are averaged over all components and all phases at each single temperature according to eq 21. The results of this latter evaluation are given in Table 7.

(20)

In this study the potential of additional parametrizations is studied. Besides the homosegmented approach (parameterization case A) and the heterosegmented approach by Schäfer et al.3 (parameterization case B) three additional parametrization cases are considered. Table 4 summarizes all cases. Table 4. Considered Parametrization Cases Evaluating PCPSAFT and Its Derivatives case

characteristics/estimated parameters

A

nC12al as homosegmented compound. m, σ and ϵ/k from PLV and ρSat L data for fixed μ = 2.88 D nC12al as heterosegmented compound. All parameter values as in ref 14 nC12al as heterosegmented compound with decane as tail. m, σ, ϵ/k and μ of head, knC12al−tail,nC12al−head and knC12al−head,DMF from PLV, ρSat L , and LLE data nC12al as heterosegmented compound. m, σ, and ϵ/k from tail and head, μ from head and 5 BIPs from PLV, ρSat L , and LLE data same as case D but using 1st alternative model for ãdipole presented in ref 34

B C D E

While the LLE measurement data are directly taken from ref 14, the pure component data pLV and ρSat L are obtained as evaluated results from the NIST ThermoData Engine version 7.140 (available in Aspen Plus V8.2) in a temperature range between 320 and 570 K. The parameter estimation problem is solved with MATLAB’s “lsqnonlin” using the trust region reflective algorithm. In all cases the parameter values published in14 were used as initial guess for optimization. Additionally, several further initial guesses were tried, which led to the same solutions. The results obtained by simulation studies using the obtained parameter

ARDTorg,pol = NM NC

⎛ |x org,exp − x org,cal| n,i n,i

∑ ∑ ⎜⎜ n=1 i=1

100 · 2·NM ·NC



xnorg,exp ,i

+

|xnpol,exp − xnpol,cal |⎞ ,i ,i ⎟ ⎟ xnpol,exp ⎠ ,i

(21)

All results show that the hs-PCP-SAFT approach from ref 3 (case B) clearly outperforms the homosegmented case A, and again that the cases C, D, and E outperform case B, guiding to a much better accordance between experiments and simulation. In particular the results of cases D and E, the ones with the G

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

Article

Industrial & Engineering Chemistry Research

Table 5. Pure Compound Parameters of n-Dodecanal from Different Parametrization Cases from Table 4 and Their Influence on Pure Component Properties Estimation vapor pressure case

m

σ

ϵ/k

μ

ARD

max RD

ARD

max RD

see Table4

[−]

[Å]

[K]

[D]

[%]

[%]

[%]

[%]

6.3108

3.5401

252.73

2.88

0.36

1.08

0.10

0.40

4.6627 1.5599

3.8384 3.0601

243.87 220.56

4.72

5.19

8.70

9.51

2.88

4.6627 1.6906

3.8384 3.4410

243.87 175.78

1.62

2.59

16.12

17.71

3.42

4.5183 1.1854

3.9087 2.9032

236.51 181.36

1.81

6.03

11.86

14.35

3.23

3.1010 2.4346

4.5477 2.7745

350.00 50.00

0.93

3.32

16.12

17.72

3.97

a

A B −tail −head C −tail −head D −tail −head E −tail −head a

sat. liq. density

New Fit Due to Major Deviations When Using Parameters From ref 3.

Table 6. Binary Interaction Parameters Resulting from Different Parametrization Cases from Table 4 case

binary interaction parameters

see Table4

DMF

A

DMF n-dodecanal DMF n-dodecanal DMF n-dodecanal DMF n-dodecanal n-dodecanal DMF n-dodecanal n-dodecanal

B C D

E

0.01 head

0.0075

head

−0.561 + 0.0017·T

head tail

0.006 −0.002

head tail

0.075 −0.001

decane 0.1159−0.000315·T 0.015 0.1159−0.000315·T −0.032 0.1159−0.000315·T −0.107 0.1159−0.000315·T −0.075 −0.074 0.1159−0.000315·T −0.075 −0.075

n-dodecanal tail

0.1159−0.000315·T −0.031 0.1159−0.000315·T −0.107 −0.002 0.044 −0.001 0.075

Figure 3. ARD analysis of molar fractions of each component in each phase based on different parametrization cases from Table 4: (left) T = 283.15 K; (center) T = 288.5 K; (right) T = 298.15 K.

most parameters that are fitted, show the best agreement. A closer look shows that the fitting improvement is particularly noticeable for the compositions of the organic phase, while the polar phase still shows relatively large deviations. Additionally it can be seen that the fitting is worse at lower temperatures. This particularly occurs because for 283.15 K there are relatively many measurements near the critical point that do not get fitted properly. It should be remarked that these results are not meant to be interpreted as better parameter values for n-dodecanal and its binary interactions. They rather represent parameter values

Table 7. ARD Analysis of Different PCP-SAFT Parametrizations at Different Temperatures (See Equation 21) ARDorg,pol [%] T

case see Table 4

T = 283.15 K

T = 288.15 K

T = 298.15 K

A B C D E

30.17 24.36 23.08 16.80 14.20

24.02 18.71 18.15 12.44 12.41

21.30 18.06 14.62 10.18 9.59 H

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

Article

Industrial & Engineering Chemistry Research

Table 8. Modeling Results for DMF−Decane, 1-Dodecene−DMF and 1-Dodecene−DMF−Decane and Number of Required BIPs ARDorg,pol T

temp

a

mixture

[K]

PCP-SAFT

SRK-MHV2

UNIFAC-Do

DMF/decane 1-dodecene/DMF 1-dodecene/DMF/decane

283.15−347.65 288.15−333.15 298.15 333.15 343.15

9.15a 7.74a 12.67b 2.89 12.03

2.67b 2.28b 4.70c 5.45 4.19

10.55 14.43 5.81 23.07 26.84

Total number ofrequired BIPs: 2. bTotal number ofrequired BIPs: 4. cTotal number ofrequired BIPs: 11.

obtained by applying state of the art local minimization techniques after augmenting the search space according to Table 4. Though some BIPs show relatively large absolute values, compared to other BIPs used in modeling with PCPSAFT, they are still in an acceptable range for purely fitting purposes. From another perspective, these relatively large values of ki,α,j,β could be seen as an indicator that the model structure is not as good as required for a better parameter identification. A particular interesting case is given by the results of the parametrization case E. Though it leads to the best fit, the fact that the obtained parameter values are in most cases in the parameter bounds (which represent rather unrealistic values) could lead to the conclusion that the model structure is not proper. However, if the variable bounds are set close to the solution of case D the results of the parametrization case E and those of case D regarding the parameter values and the resulting ARD are very similar. For further evaluations the results of the parametrization D will be considered. 4.2. LLE Results of Different Models and Parametrizations. While the former section focused on the description of the n-aldehyde containing system with PCPSAFT and its derivatives, this section focuses on the comparison of the models discussed in section 3. On the basis of the available experimental data (Table 2) and the estimation procedure described in the Supporting Information, NRTL BIPs for the description of the experimental data with the SRK-MHV2 model were obtained. The BIPs are listed in Table 1 of the Supporting Information. First, the binary and ternary systems are compared that have been successfully modeled by PCP-SAFT.3 After that, the particularly important ternary system containing the long chain aldehyde is discussed. 4.2.1. Binary Systems and Dodecene/DMF/Decane System. Table 8 shows a comparison of the average relative devation of the calculated molar fractions of different models against the experimental data as defined by eq 21. In addition to the ARD, the number of estimated BIPs is provided. While the SRK-MHV2 model provides a better fit than the PCP-SAFT EoS, it should be considered that the latter approach requires less BIPs parameters still providing an acceptable accuracy. The difference in the number of parameters is mainly due to the fact that the SRK-MHV2 approach requires BIPs for the subsystem dodecene−decane (obtained from ternary LLE data), while PCP-SAFT predicts the ternary system successfully based only on BIPs of the subsystems DMF−decane and dodecene−DMF. From Table 8 it can also be seen that the UNIFAC−Dortmund approach provides acceptable results for the binary systems over the considered temperature range, but no good temperature dependence for the ternary system. It can be concluded that

both the SRK-MHV2 and the PCP-SAFT model describe these systems with acceptable accuracy. 4.2.2. Ternary n-Dodecanal/DMF/Decane System. In section 4.1 the influence of the number of BIPs on the quality of the LLE correlation was analyzed for the hs-PCP-SAFT approach. An analogous procedure can be done with the SRKMHV2 model. Considering the BIPs for DMF−decane as already known and fixed, according to eqs 36 and 37 there is a total of up to 20 BIPs that can be considered to describe the binary interactions between n-dodecanal−DMF and n-dodecanal−decane. Following the approach described in the Supporting Information three different converging parametrization cases are found, which differ from each other by the amount of parameters to be estimated. Table 9 lists the Table 9. Considered Parametrizations using SRK-MHV2 with NRTL for the Description of Ternary n-Dodecanal (1), DMF (2), Decane (3) System. BIPs for (1) and (2) Were Previously Obtained from Binary LLE Data case

estimated parameters

A1 B1 C1

b1,2, b1,3, b2,1, b3,1 a1,2, a1,3, a3,1, b1,2, b1,3, b2,1, b3,1 a1,2, b1,2, b1,3, b2,1, b3,1, c1,3, f1,2, f1,3, f 2,1, f 3,1

parameters that were estimated in each case. To simplify matters, components are numbered as follows, n-dodecanal (1), DMF (2), and decane (3). The obtained parameter values are listed in Table 2 of the Supporting Information. To compare the results of these parametrizations against the ones obtained by PCP-SAFT, Figure 4 shows the ARD of the molar fraction of each component in each phase (over all temperatures) against the best PCP-SAFT parametrization

Figure 4. ARD analysis of molar fractions of each component in each phase weighted over all temperatures. Comparison of best PCP-SAFT parametrization (D) against NRTL (A1), (B1), and (C1) (see Table 9). I

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

Article

Industrial & Engineering Chemistry Research

Figure 5. Comparison of correlative capabilities of considered parametrizations for the hs-PCP-SAFT (a−c) and the SRK-MHV2 (d−f) EoS models. Diagrams a, b, and c correspond to parametrization cases B, C, and D (see Table 4), diagrams d, e, and f correspond to cases A1, B1, and C1 (see Table 9).

• SRK-MHV2 (kC12en,C13al = 0): using estimated parameters from Table 1 of the Supporting Information but neglecting binary interaction between 1-dodecene and ndodecanal • SRK-MHV2 (kC12en,C13al ≠ 0): using estimated parameters from Table 1 of the Supporting Information; BIPs between 1-dodecene and n-dodecanal were estimated with UNIFAC Dortmund • hs-PCP-SAFT (PAR D): with pure component parameters and BIPs from parametrization case D (see Tables 5 and 6) 4.3.1. Comparison of Similar Reactor Product Flows. To study the influence of the decanter feed conditions on the process behavior three different (but to a certain extent comparable) decanter feed concentrations listed in Table 10 are

(case D). Figure 5 additionally gives an overview of the correlative capabilities of the most relevant considered parametrizations while giving an insight into the temperature dependence. As shown the large deviations by the parametrization cases B and C are particularly evident. Regarding the molar fractions related to the organic phase, all parametrizations shown in Figure 4 present a comparable accuracy. On the other hand it is clear that the SRK-MHV2 allows a much better correlation of the molar fractions in the polar phase. While the ARD shows relatively comparable values for all parametrization cases with SRK-MHV2 (Figure 4), the results for case A1 (Figure 5) could lead to the assumption that the temperature dependence can be neglected. However, the results of the cases B1 and C1 emphasize the temperature dependence available in the experimental data. Since there is no major difference in the results of cases B1 and C1 in the subsequent simulation studies, the BIPs obtained by the parametrization case B1 will be considered. 4.3. Process-Wide Simulation Studies. The simulation case studies are determined by the modeling goals discussed in section 2.2, the modeling assumptions from section 2.2.1, and the model formulation from section 2.3. The evaluated thermodynamic models, with their corresponding parametrizations are listed next. It needs to be considered that no binary information for the 1-dodecene/n-aldehyde pair is available: • UNIFAC-Do: using default functional group distribution and functional group parameters from Aspen Plus V8.2 • hs-PCP-SAFT: with parameters from ref 14; equivalent to case B from Table 4

Table 10. Considered Decanter Feed Conditions feed

wC12en

wC13al

wDMF

wC10an

1 2 3

0.01 0.03 0.05

0.05 0.15 0.25

0.47 0.41 0.35

0.47 0.41 0.35

considered. The compositions are chosen so as to study the influence of the amount of produced aldehyde at similar DMF/ decane ratios. Note that by fixing the decanter’s feed conditions, outlet pressure, and temperature, not only the decanter but also the considered simplified model of the reactor/decanter system (see eq 2.3) become fully specified. After solving the decanter model for the aforementioned J

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

Article

Industrial & Engineering Chemistry Research

Figure 6. Resulting process-wide operation parameters subject to different decanter feed conditions (Table 10) and decanter separation temperatures.

the description of the process behavior, since it leads to almost identical results as in the case neglecting those interactions. Only in the operation region characterized by high aldehyde compositions and relatively high temperatures do the differences become noticeable. Regarding the temperature dependence of the evaluated process parameters, it can be seen that the hs-PCP-SAFT approach using parameters from ref 3 presents in most cases similar trends as the UNIFAC−Dortmund calculations. On the other hand both SRK-MHV2 approaches and the hs-PCPSAFT with new fitted parameters show comparable trends. Though differences between the models appear, for most considered cases, in particular in regions of low temperatures between 270 and 285 K, temperature changes do not seem to play a major role. This changes however at higher temperatures where differences become more relevant. In particular, the simulation case with the highest aldehyde feed (Feed 3) process and aldehyde presents improving process yields YC13al separation efficiencies εC13al with higher temperatures. These results, however, should not lead to the conclusion that higher separation temperatures lead in general to better process operation. The observed strong changes are rather related to the fact that the considered operation region, characterized through relatively high temperatures and aldehyde concen-

specifications, the component balance for the aldehyde around the reactor (eq 1 for c = C13al) can be solved to obtain the corresponding reaction rate r required to achieve the desired operation with full recycle. Having obtained the value of r all further unknown variables follow from mass balances around the mixer and reactor. Additionally, it should be noted that because of the considered model formulation all results are size independent and therefore independent of the exact value of the decanter feed flow. In addition to the overall process yield Yprocess C13al defined in eq 11, two additional performance indicators are considered to evaluate the process behavior. These are the product-recycle ratio Fprod/Frec = Fs=9/Fs=10 and the aldehyde separation efficiency εC13al (eq 22) which defines the fraction of produced aldehydes that reaches the product phase. εC13al =

Fs = 9·xs = 9,C13al Fs = 8·xs = 8, C13al

(22)

Figure 6 summarizes the influence of the decanter feed conditions from Table 10 and the decanter separation temperature on the aforementioned process parameters. Several facts follow from Figure 6. First, the results clearly show that considering BIPs for the dodecene/aldehyde binary system as estimated by UNIFAC−Dortmund almost has no influence on K

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

Article

Industrial & Engineering Chemistry Research

coupled reactor/decanter system. In this study the recycle flow is fixed. Though different values of the recycle flow lead to different corresponding reactor conversions XR as shown in Figure 8, the application of a basic model analysis, explained

trations, is near to the critical region that divides the mixing gap from the homogeneous regime. Such operation regions are usually characterized by strong phase distribution changes with temperature, as indicated in Figure 6c. 4.3.2. Validation against Available Miniplant Data. For evaluation purposes, the simplified model described in section 2.3 in conjunction with the thermodynamic models listed at the beginning of section 4.3 is validated against published experimental data from ref 18. The process information provided in ref 18 is by no means extensive and basically limited to time profiles of plantwide conversions and yields that reach a stable, nearly steady state behavior. Though in some of the following analyses values of the reactor conversion XR are presented in order to characterize the process behavior, it should be noted that an experimental value of XR cannot be obtained from the published data. This would require a more detailed measurement setup that provides at least accurate and reliable flow and composition data of the recycle stream. The basic specifications from ref 18 required for the model validation are listed in Table 11. For the sake of simplicity, these conditions are translated into corresponding total inlet flows and compositions as given in Table 12.

Figure 8. Reactor conversion dependent on a thermodynamic model at conditions from Table 11.

next, demonstrates that the process-wide aldehyde yield Yprocess C13al is independent of the recycle flow rate. Due to the conservation of mass, for given total flow rate Ffeed and feed molar fractions feed feed of the solvents xfeed DMF and xC10an it follows that xs=9,DMF = xDMF and xs=9,C10an = xfeed . Knowing the temperature, pressure, and C10an two molar fractions of a two-phase system in equilibrium, in this case xs=9,DMF and xs=9,C10an, all further equilibrium molar fractions, among other xs=9,C12en and xs=9,C13al, result from the Gibbs’ phase rule. As shown in eq 23, these latter molar fractions are decisive for the evaluation of the Yprocess C13al , which is shown to be independent of the recycle flow. xs = 9,C13al process YC13al = xs = 9,C12en + xs = 9,C13al (23)

Table 11. Process conditions18 V̇ feed,org wfeed C10an wfeed nC12en Tdecanter a

60 mL h−1

V̇ Make‑Up

4.5 mL h−1

0.724 0.276 278.15 K

wMake‑Up DMF

0.9927a

pprocess

20 bar

For the sake of simplicity wDMF is considered to be 1.

Table 12. Process Feed Conditions from Table 11 as Total Inlet Molar/Mass Flows with Corresponding Molar/Mass Fractions. Density Was Calculated with PCP-SAFT Ffeed xfeed nC12en xfeed DMF xfeed C10an

0.351 mol h−1

Ṁ feed

47.942 g h−1

0.204 0.162 0.634

wfeed nC12en wfeed DMF wfeed C10an

0.252 0.087 0.661

The results of Figure 7 clearly show that the modeling approaches derived in this contribution agree much better with the experimental yield than the modeling approach proposed in ref 3. In particular the SRK-MHV2-based approaches show a better agreement with the experiments. Note that there are no results based on UNIFAC−Dortmund. With this latter model no physically feasible solution could be found for the considered feed conditions. 4.3.3. Influence of Density Calculation. A central step of the model validation discussed in section 4.3.2 is given by the conversion of volumetric flow data in molar or mass specific flow. For this step a value of the average feed density is required. Because pure compound densities calculated with PCP-SAFT show the best agreement with available experimental data, in the model validation discussed above the average density was calculated with PCP-SAFT. Table 13 lists the corresponding process conditions when the average density is calculated with the SRK-MHV2 model. Clearly the model

Following the model assumptions given in section 2.2.1 the best way of validating the model is to compare the simulated process-wide aldehyde yield Yprocess C13al against the corresponding experimental value. Figure 7 shows these results based on the conditions given in Table 11 using the flows according to Table 12. Having fixed the values of the feed flows of all components and the temperature and pressure of the decanter, an additional degree of freedom needs to be fixed in order to simulate the

Table 13. Process Conditions Provided by Reference 18 as Total Inlet Molar/Mass Flows with Corresponding Molar/ Mass Fractions. Average Density Calculated with SRKMHV2 Ffeed xfeed nC12en xfeed DMF xfeed C10an

Figure 7. Comparison of experimentally obtained aldehyde yield against simulations using different thermodynamic models. L

0.281 mol h−1

Ṁ feed

38.929 g h−1

0.211 0.134 0.655

wfeed nC12en wfeed DMF wfeed C10an

0.256 0.071 0.673

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

Article

Industrial & Engineering Chemistry Research

yields at lower solvent compositions, and high yields at higher ones. At lower decane and DMF molar fractions, both models show rather comparable yield values, while in regions of higher DMF and decane molar fractions, the hs-PCP-SAFT approach predicts higher yields. Figure 10 also shows how relatively small changes in solvent composition can have a major impact on the process performance. As an example, using the SRK-MHV2 models, for molar fractions xs=9,DMF = 0.15 and xs=9,C10an = 0.65 a yield of Yprocess C13al = 0.62 results. Increasing the DMF molar fraction by 10% and the decane molar fraction by 4% leads to a yield of Yprocess C13al = 0.83, that is, to an increase of yield of 33.9%. Similarly decreasing the solvent molar fractions by the same percent leads to a yield of Yprocess C13al = 0.47, that is, to a yield decrease of 24%. Using the hs-PCP-SAFT approach (PAR D) for the same reference conditions and changes leads to yields of process Yprocess C13al = 0.90 and YC13al = 0.46, with changes of +48% and −26%, respectively. These results clearly show the importance of using the right density calculation in model validation and additionally show a simplified method to quantify the requirements of solvent distribution for a given process performance. Regarding the considered multicomponent phase equilibria, both, the SRK-MHV2 and hs-PCP-SAFT in conjunction with parameters presented in this work, provide similar, feasible results in the relevant operation regions that have been studied experimentally. Having validated the model against experimental data, an important aspect that should be evaluated previous to the use of such models for advanced model-based methods is related to the computational effort. In particular, it is of great interest to quantify how the selection of a certain rigorous thermodynamic model affects the overall computational times. There are several obvious factors that can make the computational comparison of the considered thermodynamic models particularly awkward, for example, the efficiency of the implemented code and derivatives, the fact whether analytic solutions exist, among others. To allow for a fair comparison, both models are implemented as explained in section 3.2 following the published expressions for the reduced residual Helmholtz energy ãres. For the evaluation of computational effort, the CPU times required for the solution of single converging Newton steps are measured. The considered Newton steps refer to the solution of the pT-flash model as formulated in Appendix B. For different equilibrium calculations concerning the hydroformylation process described in Section 2, Table 14 shows the

selection has a major influence on the solution of the density root problem. The influence on the process behavior is plotted in Figure 9. Since at these feed conditions a solution with UNIFAC−Dortmund could be found, the results with this model are considered as well.

Figure 9. Reactor conversion dependent on thermodynamic model at feed conditions from Table 13.

A comparison of Figures 8 and 9 shows major differences in the simulated process behavior. Comparable differences are also noticeable in the process-wide aldehyde yield Yprocess C13al , which has now values around 0.54 for phase equilibria calculated with SRK-MHV2. Here a drawback of the SRK-MHV2 approach can be noticed. Though the phase equilibria is very well correlated, care should be taken in the evaluation of volumetric properties. The results of different approaches show the major influence of density calculations on the component stream flows and hereby on the molar fractions going into the process. As previously discussed in section 4.3.2 these molar fractions play a major role on the process-wide aldehyde yield Yprocess C13al . Because different models for the density root problem lead to different feed compositions and these feed compositions directly affect the process-wide performance, we analyze this performance as a function of the overall feed composition, not only to quantify the uncertainty related to the density root problem but also to capture the influence of the solvent composition on the operation. Figure 10 shows the results of such an analysis considering two different thermodynamic models, SRK-MHV2 (kdoce,trid = 0) and hs-PCP-SAFT (PAR D), and a relatively narrow domain of solvent compositions. Regarding the influence of the solvent feed composition, both models show relatively simple and similar trends: Low

Figure 10. Process-wide aldehyde yield as a function of total feed composition: (left) using SRK-MHV2 model; (right) using hs-PCP-SAFT. M

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

Article

Industrial & Engineering Chemistry Research

acceptable results. Given the predictive nature of the UNIFAC−Dortmund approach, the results provided by it are satisfactory, especially in the region of low process temperatures. For these systems the PCP-SAFT approach with a relatively low number of BIPs also provides a good agreement with the measurements and reproduces the temperature dependency better than UNIFAC−Dortmund. The SRKMHV2 approach provides the best fit, however, with the use of a relatively large number of BIPs. For the ternary systems containing aldehydes the situation is different. Here neither UNIFAC−Dortmund, nor the currently published best parametrization using hs-PCP-SAFT20 describe the experimental data in a satisfactory way, as required for process-wide simulations. Therefore, parametrization studies were carried out with the hs-PCP-SAFT and the SRK-MHV2 approach. An expansion of the parameter space proposed by Schäfer et al.3 leads to parameters that show a better accordance with experimental LLE data. For the SRK-MHV2 approach several parametrization cases were considered in order to find a good fit with a low number of variables. Though both models show a good accordance with experimental data, also regarding the temperature dependence of the LLE, in regions with high aldehyde molar fractions, close to the critical point, SRK-MHV2 provides a better fit. A comparison of the computational effort clearly shows that the evaluation effort required by equilibrium calculations using hs-PCP-SAFT can be significantly larger than that for the corresponding calculation using SRK-MHV2. The ratio between the cost of evaluating the hs-PCP-SAFT model and the cost of evaluating the SRKMHV2 grows with the number of compounds, the number of polar compounds, and with the appearance of heterosegmented components (see Table 14). To evaluate the performance of the thermodynamic models against experimental miniplant data, and constrained by incomplete information, for example, regarding isomeric components, a simplified model has been developed and evaluated in several simulation cases. Simulation studies varying important parameters, such as the decanter output concentrations and the decanter separation temperatures, give a first impression of the influence of the thermodynamic LLE modeling approach on the process-wide behavior. In general, results show that UNIFAC−Dortmund and hs-PCP-SAFT as parametrized in ref 14 show opposite temperature dependency trends as hs-PCP-SAFT and SRK-MHV2 with parameters estimated in this contribution. Interestingly these simulation studies indirectly show operation regions in which the model selection does not have a major influence on the predicted process behavior. These operation regions are characterized by low operation temperatures and low aldehyde molar fractions. However, in feasible operation regions with higher aldehyde concentrations major differences between the approaches appear. This is particularly valid for the comparison between experimental miniplant and process simulation results, which clearly shows how the selected model for the LLE description has a major role on the quantification of the process performance. While the SRK-MHV2 and the hs-PCP-SAFT approach with parameters estimated in this contribution show a very good agreement with the experimental miniplant data, UNIFAC−Dortmund and hs-PCP-SAFT with parameters from ref 3 results highly differ from these. On the basis of the obtained results, both the hs-PCP-SAFT and the SRK-MHV2 with parameters estimated in this

Table 14. Comparison of Computational Effort Required for Evaluation of Multicomponent Equilibrium Calculations. The Considered Times Refer to the Evaluation of a Newton Step number of BIPs used

decane/CO DMF/H2 decane/DMF/ 1-dodecene/DMF/ 1-dodecene/DMF/decane n-dodecanal/DMF/ decane DMF/decane/H2/CO 1-dodecene/decane/DMF/ndodecanal

hs-PCP-SAFT

SRK-MHV2

ths − PCP − SAFT tSRK − MHV2

1 1 2 2 4 7

4 6 4 4 14 14

2.352 3.619 3.793 3.980 5.615 9.594

6 9

22 24

9.447 17.020

ratio between the cost of evaluating phase equilibrium calculations (pT-flash) with hs-PCP-SAFT and the cost of evaluating the same equilibrium calculations with the SRKMHV2 model. It should be further remarked that hs-PCPSAFT is implemented in such a way that it covers the requirements of the corresponding mixtures to which it is applied. For example, in all cases, where no aldehyde appears, the hs-PCP-SAFT reduces to the homosegmented PCP-SAFT approach,22 if additionally no polar component is part of the mixture (as in the case of CO−decane solubility) the equation system corresponds to the most basic PC-SAFT.21 Regarding the computational effort Table 14 shows clearly that the computational effort of using hs-PCP-SAFT in comparison to SRK-MHV2 grows with the number of components, number of polar components and with the appearance of heterosegmented components. From a purely correlative perspective, for the considered simulation studies it is not necessary to use advanced approaches such as the EoS/GE models or hs-PCP-SAFT. It would have been enough to use a properly parametrized GE model with properly parametrized parameters to obtain similar results with an even lower computational effort. However, since a full rigorous process simulation should also consider the supercritical components, we concentrate on models that can handle the complete multicomponent system in a thermodynamically consistent way. The approach of using rigorous thermodynamic models also follows the trend in process engineering toward the use of “complex” design models in advanced model-based applications.

5. SUMMARY/CONCLUSIONS On the basis of available LLE data for binary and ternary systems consisting of DMF, 1-dodecene, decane, and ndodecanal (n-tridecanal) and on published experimental data from a miniplant for the hydroformylation of 1-dodecene in a corresponding TMS system, the performance of a predictive GE-model (UNIFAC−Dortmund), an EoS/GE model (SRKMHV2 with NRTL) and a SAFT-based EoS (hs-PCP-SAFT) were evaluated. Particular attention was put on the EoS based approaches, since they are directly applicable to a process-wide simulation, including the consideration of supercritical components. Regarding the binary and ternary systems excluding the aldehyde (Table 8), all examined models provide at least N

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

Article

Industrial & Engineering Chemistry Research

For supercritical compounds, characterized by Tr,i ≥ 1, the Boston−Mathias Extrapolation, eq 31, is considered

contribution seem to be suitable for process design and optimization tasks. On the one hand, taking a closer look the SRK-MHV2 model shows slightly better correlative capabilities with a lower computational effort, however, with a larger number of binary interaction parameters. On the other hand, while phase equilibrium compositions are almost equally well described by both approaches, particular care should be taken with the calculation of volumetric properties. For the latter the hs-PCP-SAFT approach provides better results.

Ωipure = [exp(ci· (1 − (Tr , i)di ))]2

where di = 1 +

ci = 1 −

gE = R·T (25)

Gj , i = exp( −αj , i·τj , i)

∑ xi·bipure

For the description of the mixing rule of amix a new variable αmix = amix/(bmix·R·T) is introduced. To calculate amix the larger root of the following quadratic equation needs to be solved with respect to αmix (MHV2 mixing rule)

i=1 NC

=

gE + R·T

τj , i = aj , i +

) + q2 ·((αmix ) −

∑ xi·(αi i=1

i=1

(27)

where q1 = −0.4783 and q2 = −0.0047. The pure component parameter bpure is calculated as i bipure = 0.08664·

= [1 +

(37)

For the comparison of the hs-PCP-SAFT against the SRKMHV2 a standard flash model was solved. For the sake of simplicity only the system with two phases (A and B) in equilibrium is formulated, where ΦA is the phase fraction of phase A. The following equations build the standard two-phase flash model: Component balances (i = 1, ···, NC)

(28)

(R ·Tic)2 pure ·Ωi pic

mipure · (1

(36)

B.1. Multicomponent LLE

R ·Tic pic

(29)

where Ωpure here represents the so-called Soave alpha function, i which is a function of the acentric factor ωi and the reduced temperate Tr,i = T/Tci . Following the default setting of the SRK-MHV2 implementation in Aspen Plus,42 for undercritical compounds, characterized by Tr,i < 1, the Schwartzentruber−Renon−Watanasiri alpha function, eq 30), is considered. Ωipure

+ ej , i ·log(T ) + f j , i ·T

B. GENERAL FLASH MODEL

The pure component parameter apure is calculated by i aipure = 0.42747·

T

where τ j,i is in general asymmetric τj,i ≠ τi,j and so are the parameters aj,i, bj,i, ej,i, and f j,i. These parameters are estimated to fit LLE data. αj,i and hence cj,i and dj,i are symmetrical parameters. While these parameters can also be estimated, recommended values for different types of mixtures exist.

))

⎛ bmix ⎞ pure ⎟ ⎝ bi ⎠

bj , i

αj , i = cj , i + dj , i·(T − 273.15 K)

pure 2

∑ xi·log⎜

(35)

and eq 37, respectively.

NC

∑ xi·αi

(34)

Using the definition of interaction parameters as considered in Aspen Plus, τj,i and αj,i are evaluated with eq 36

(26)

i=1

q1·(αmix −

⎛ ∑NC τ ·G ·x ⎞ j=1 j,i j,i j ⎟ ∑ ⎜⎜xi· NC ∑k = 1 Gk , i ·xk ⎟⎠ i=1 ⎝ NC

where

NC

2

(33)

In addition to eqs 28 to 30 defining and requiring pure component properties, in eq 27 a value for gE is additionally required. We consider the NRTL activity coefficient model

For bmix a linear mixing rule is considered

pure

1 di

NRTL as Underlying Activity Coefficient Model

a ̃res(vmix,xmix , T ) =

NC

(32)

mpure is calculated as a function of the acentric factor: mpure = i i 0.48508 + 1.55171·ωi − 0.15613·(ωi)2.

or by its respective expression for the reduced residual Helmholtz function ãres38

bmix =

mipure − (P1i + P 2i + P 3i ) 2

and

A. THE SRK-MHV2 MODEL The SRK equation of state for a mixture is given by either its pressure explicit representation amix R·T p= − vmix − bmix vmix ·(vmix + bmix ) (24)

a res ̂ R·T ⎛ ⎛ bmix ⎞ amix b ⎞ = − log⎜1 − ·log⎜1 + mix ⎟ ⎟− vmix ⎠ bmix ·R ·T vmix ⎠ ⎝ ⎝

(31)

xiin = Φ A ·xiA + (1 − Φ A ) ·xiB

Isofugacity conditions (i = 1, ···, NC) xiA ·φi A = xiB·φi B

(39)

and summation equation NC

∑ (xiA − xiB) = 0

0.5

− (Tr , i) )

i=1

2

(38)

2

(40)

The fugacity coefficients φAi and φBi are evaluated according to eq 17. Following the need to find the values of density ρ (or

− (1 − Tr , i) ·(P1i + P 2i ·Tr , i + P 3i ·(Tr , i) )]

(30) O

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

Industrial & Engineering Chemistry Research



molar volume v) and compressibility factor Z for each phase A and B at a given pressure p and temperature T, for each phase a density root problem consisting of eqs 41 and 42 needs to be solved. Here only the formulation for the PCP-SAFT EoS (w.r.t. ρ) is given. An analogous formulation for the molar volume v results from eq 16. For the sake of a simpler notation, the superscript indicating the phase is neglected in the following equations. Density root problem ⎛ ∂a ̃res ⎞ Z = 1 + ρ·⎜ ⎟ ⎝ ∂ρ ⎠

(41)

⎛ 10 Å ⎞3 p = Z ·k·T ·ρ ·⎜10 ⎟ ⎝ m⎠

(42)

ARD = Average relative deviation [%] ARD = 100·

Within the pure component parameter estimation discussed in section 4.1, it is necessary to calculate the vapor pressure pLV and saturated liquid density ρSat,M (in kg/m3) at a given L temperature T. These properties can be obtained by a reduced variant of the equation system shown above, consisting of eqs 41 and 42 for the vapor and liquid phase (j = V,L) and an isofugacity condition for the pure component (eq 43).

max RD = 100·max

(43)

After solving the equations system, at a given temperature can be obtained directly from the compressibility factor according to eqs 44 and 45, where the varialble M stands for the molar mass of the considered component.



ρLSat, M = M /vLSat

(45)

ASSOCIATED CONTENT

yiexp

Greek Letters

S Supporting Information *

α = EoS parameter, α = a/(b · R · T) αj,i = nonrandomnes parameter for binary j/i interactions ϵ = depth of pair potential [J] ϵ/k = PC-SAFT pure component parameter [K] ε = separation efficiency [−] μ = dipole moment [D] ν = stoichiometric coefficient [−] Ω = Soave alpha function ω = acentric factor [−] Φ = phase fraction [−] ρ = total number density of molecules [Å−3] ρ = mass density [kg/m3] σ = segment diameter [Å] τj,i = Temperature dependent energy interaction parameter between j-i pair of molecules [−]

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.iecr.5b03328. Estimation of NRTL binary interaction parameters (PDF)



yical − yiexp

MHV2 = modified Huron/Vidal mixing rules of second order MSS = micellar solvent system NRTL = nonrandom two-liquid local composition model PC-SAFT = Perturbed Chain Statistical Associating Fluid Theory PCP-SAFT = Perturbed Chain Polar-Statistical Associating Fluid Theory pT-flash = equilibrium stage with fixed pressure and temperature RD = relative deviation [%] SAFT = Statistical Assosication Fluid Theory SRK-MHV2 = SRK EoS with MHV2 mixing rules SRK = Soave−Redlich−Kwong TMS = thermomorphic multicomponent solvent UNIFAC−DO = UNIFAC Dortmund model

ρSat,M L

(44)

NM y cal − y exp 1 · ∑ i exp i NM n = 1 yi

BIP = binary interaction parameter CPA = Cubic-Plus-Association CXL = CO2-expanded liquids DMF = N, N-dimethylformamide EoS = equation of state EoS/GE = EoS model with mixing rules based on a GE model GLE = gas−liquid equilibrium GXL = gas-expanded liquids hs-PCP-SAFT = heterosegmented PCP-SAFT IL = ionic liquids LLE = liquid−liquid equilibrium max RD = maximum relative deviation [%]

B.2. Pure Component

vLSat = ZL·R ·T /p LV

NOMENCLATURE

Abbreviations

In many solution strategies, the density root problem is solved within an internal loop when solving eqs 38 to 40.

φL = φV

Article

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel.: +49 30 314 29945. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work is part of the Collaborative Research Centre ”Integrated Chemical Processes in Liquid Multiphase Systems” coordinated by the Technische Universität Berlin. Financial support by the Deutsche Forschungsgemeinschaft (DFG) is gratefully acknowledged (TRR 63). The authors thank Dr.-Ing. Tilman Barz for enriching discussions on parameter estimation.

Indices

α = segment number α β = segment number β c = component number c i = component number i j = component number j n = measurement number n

P

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

Article

Industrial & Engineering Chemistry Research s = stream number s Latin Letters

A = Helmholtz free energy [J] a = EoS energy parameter in physical term [J m3/kmol2] ã = reduced Helmholtz free energy [−] â = molar Helmholtz free energy [J/mol] aj,i = asymmetric BIP in eq 36 b = EoS covolume parameter [m3/mol] bj,i = asymmetric BIP in eq 36 ci = pure component parameter in eq 31 cj,i = symmetric BIP in eq 37 di = pure component parameter in eq 31 dj,i = symmetric BIP in eq 37 ej,i = asymmetric BIP in eq 36 F = molar flow [kmol/s] f j,i = asymmetric BIP in eq 36 Gj,i = coefficient as defined in eq 35 gE = specific excess Gibbs energy [J/mol] GE = excess free enthalpy [J] k = Boltzmann constant k = 1.38066 × 10−23 J/K ki,α,j,β = BIP between segment α of component i and segment β of component j M = molar mass [kg/kmol] Ṁ = mass flow [kg/s] m = number of segments per chain [−] m = pure component parameter in SRK EoS, mi = f(ωi) N = total number of molecules [−] NC = number of components [−] NLLE = number of measured tie lines NM = number of measurements [−] NPM = number of experimental pure data OF = value of objective function P1 = polar parameter in eq 30 P2 = polar parameter in eq 30 P3 = polar parameter in eq 30 pLV = vapor pressure [Pa] p = pressure [Pa] q1 = mixing rule constant, see eq 27 q2 = mixing rule constant, see eq 27 R = universal gas constant [8.314 J/(mol K)] r = reaction rate [kmol/s] T = temperature [K] t = time [s] V̇ = volume flow [m3/s] v = specific molar volume [m3/mol] w = mass fraction [kg/kg] X = conversion [−] x = molar fraction [mol/mol] x⃗ = vector of molar fractions Y = yield [−] y = general variable y Z = compressibility factor [−]



feed = related to the feed hc = contribution of hard chain system L = liquid phase low = related to lower bound M = related to mass Make-Up = related to Make-Up stream mix = mixture org = organic phase pol = polar phase process = process-wide prod = related to product stream pure = pure compound property R = related to the reactor r = reduced term (relative to a critical value) rec = related to recycle stream res = residual term (deviation from ideal gas property) Sat = saturated property T = at a fixed Temperature T up = related to upper bound V = vapor phase

REFERENCES

(1) Behr, A.; Fängewisch, C. Temperature-Dependent Multicomponent Solvent Systems - An Alternative Concept for Recycling Homogeneous Catalysts. Chem. Eng. Technol. 2002, 25, 143−147. (2) Shaharun, M. S.; Dutta, B. K.; Mukhtar, H.; Maitra, S. Hydroformylation of 1-octene using rhodium-phosphite catalyst in a thermomorphic solvent system. Chem. Eng. Sci. 2010, 65, 273−281. (3) Schäfer, E.; Brunsch, Y.; Sadowski, G.; Behr, A. Hydroformylation of 1 - Dodecene in the Thermomorphic Solvent System Dimethylformamide/Decane. Phase Behavior - Reaction Performance - Catalyst Recycling. Ind. Eng. Chem. Res. 2012, 51, 10296−10306. (4) Haumann, M.; Koch, H.; Hugo, P.; Schomäcker, R. Hydroformylation of 1-dodecene using Rh-TPPTS in a microemulsion. Appl. Catal., A 2002, 225, 239−249. (5) Müller, M.; Kasaka, Y.; Müller, D.; Schomäcker, R.; Wozny, G. Process design for the separation of three liquid phases for a continuous hydroformylation process in a miniplant scale. Ind. Eng. Chem. Res. 2013, 52, 7259−7264. (6) Riisager, A.; Eriksen, K. M.; Wasserscheid, P.; Fehrmann, R. Propene and 1-octene hydroformylation with silica-supported, ionic liquid-phase (SILP) Rh-phosphine catalysts in continuous fixed-bed mode. Catal. Lett. 2003, 90, 149−153. (7) Sharma, A.; Lebigue, C. J.; Deshpande, R. M.; Kelkar, A. a.; Delmas, H. Hydroformylation of 1-octene using [Bmim][PF6]-decane biphasic media and rhodium complex catalyst: Thermodynamic properties and kinetic study. Ind. Eng. Chem. Res. 2010, 49, 10698− 10706. (8) Koch, D.; Leitner, W. Rhodium-catalyzed hydroformylation in supercritical carbon dioxide. J. Am. Chem. Soc. 1998, 120, 13398− 13404. (9) Jin, H.; Subramaniam, B. Homogeneous catalytic hydroformylation of 1-octene in CO 2-expanded solvent media. Chem. Eng. Sci. 2004, 59, 4887−4893. (10) Schrader, P.; Culaguin-Chicaroux, A.; Enders, S. The water +nonionic surfactant C12E8 + 1-dodecene phase behavior across a wide temperature range. Chem. Eng. Sci. 2012, 93, 131. (11) Kiedorf, G.; Hoang, D.; Müller, a.; Jörke, a.; Markert, J.; Arellano-Garcia, H.; Seidel-Morgenstern, a.; Hamel, C. Kinetics of 1dodecene hydroformylation in a thermomorphic solvent system using a rhodium-biphephos catalyst. Chem. Eng. Sci. 2014, 115, 31−48. (12) Zagajewski, M.; Behr, a.; Sasse, P.; Wittmann, J. Continuously operated miniplant for the rhodium catalyzed hydroformylation of 1dodecene in a thermomorphic multicomponent solvent system (TMS). Chem. Eng. Sci. 2014, 115, 88−94.

Subscripts and Superscripts

A = relative to phase A assoc = contribution due to association forces B = relative to phase B c = critical value cal = calculated value dcl-head = head segment of dodecanal dipole = contribution due to dipole/dipole interactions disp = contribution due to dispersive attraction exp = experimental value Q

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

Article

Industrial & Engineering Chemistry Research

(34) Sauer, E.; Stavrou, M.; Gross, J. Comparison between a Homoand a Heterosegmented Group Contribution Approach Based on the Perturbed-Chain Polar Statistical Associating Fluid Theory Equation of State. Ind. Eng. Chem. Res. 2014, 53, 14854−14864. (35) Dahl, S. r.; Michelsen, M. L. High-pressure vapor-liquid equilibrium with a UNIFAC-based equation of state. AIChE J. 1990, 36, 1829−1836. (36) Jakob, A.; Grensemann, H.; Lohmann, J.; Gmehling, J. Further Development of Modified UNIFAC (Dortmund): Revision and Extension 5. Ind. Eng. Chem. Res. 2006, 45, 7924−7933. (37) Bischof, C. H.; Bücker, H. M.; Lang, B.; Rasch, A.; Vehreschild, A. Combining Source Transformation and Operator Overloading Techniques to Compute Derivatives for MATLAB Programs. Proceedings of the Second IEEE International Workshop on Source Code Analysis and Manipulation (SCAM 2002). Los Alamitos, CA, USA, 2002; pp 65−72. (38) Michelsen, M. L.; Mollerup, J. M. Thermodynamic Models: Fundamentals & Computational Aspects; Tie-Line Publications: Denmark, 2004. (39) López C, D. C.; Barz, T.; Körkel, S.; Wozny, G. Nonlinear illposed problem analysis in model-based parameter estimation and experimental design. Comput. Chem. Eng. 2015, 77, 24−42. (40) Frenkel, M.; Chirico, R. D.; Diky, V.; Yan, X.; Dong, Q.; Muzny, C. ThermoData Engine (TDE): software implementation of the dynamic data evaluation concept. J. Chem. Inf. Model. 2005, 45, 816− 38. (41) Xu, G.; Brennecke, J. F.; Stadtherr, M. A. Reliable Computation of Phase Stability and Equilibrium from the SAFT Equation of State. Ind. Eng. Chem. Res. 2002, 41, 938−952. (42) AspenTech, Aspen Plus V8.2; www.aspentech.com (accessed October 17, 2014).

(13) Vogelpohl, C.; Brandenbusch, C.; Sadowski, G. High-pressure gas solubility in multicomponent solvent systems for hydroformylation. Part II: Syngas solubility. J. Supercrit. Fluids 2014, 88, 74−84. (14) Schäfer, E. Simultaneous Calculation of Fluid Phase Equilibria and Interfacial Properties of Aprotic Solvent Mixtures Incorporating Experimental Studies; Schriftenreihe Thermodynamik; Dr. Hut: 2013. (15) Hentschel, B.; Kiedorf, G.; Gerlach, M.; Hamel, C.; SeidelMorgenstern, A.; Freund, H.; Sundmacher, K. Model-Based Identification and Experimental Validation of the Optimal Reaction Route for the Hydroformylation of 1-Dodecene. Ind. Eng. Chem. Res. 2015, 54, 1755. (16) Hentschel, B.; Peschel, A.; Freund, H.; Sundmacher, K. Simultaneous design of the optimal reaction and process concept for multiphase systems. Chem. Eng. Sci. 2014, 115, 69−87. (17) Kontogeorgis, G.; Coutsikos, P. Thirty Years with EoS/GE Models - What Have We Learned? Ind. Eng. Chem. Res. 2012, 51, 4119−4142. (18) Zagajewski, M.; Dreimann, J.; Behr, A. Verfahrensentwicklung vom Labor zur Miniplant: Hydroformylierung von 1-Dodecen in thermomorphen Lösungsmittelsystemen. Chem. Ing. Tech. 2014, 86, 449−457. (19) Markert, J.; Brunsch, Y.; Munkelt, T.; Kiedorf, G.; Behr, a.; Hamel, C.; Seidel-Morgenstern, a. Analysis of the reaction network for the Rh-catalyzed hydroformylation of 1-dodecene in a thermomorphic multicomponent solvent system. Appl. Catal., A 2013, 462−463, 287− 295. (20) Schäfer, E.; Sadowski, G. Liquid-Liquid Equilibria of Systems with Linear Aldehydes. Experimental Data and Modeling with PCPSAFT. Ind. Eng. Chem. Res. 2012, 51, 14525−14534. (21) Gross, J.; Sadowski, G. Perturbed-Chain SAFT: An Equation of State Based on a Perturbation Theory for Chain Molecules. Ind. Eng. Chem. Res. 2001, 40, 1244−1260. (22) Gross, J.; Vrabec, J. An equation-of-state contribution for polar components: Dipolar molecules. AIChE J. 2006, 52, 1194−1204. (23) Gross, J.; Spuhl, O.; Tumakaka, F.; Sadowski, G. Modeling Copolymer Systems Using the Perturbed-Chain SAFT Equation of State. Ind. Eng. Chem. Res. 2003, 42, 1266−1274. (24) Vogelpohl, C.; Brandenbusch, C.; Sadowski, G. High-pressure gas solubility in multicomponent solvent systems for hydroformylation. Part I: Carbon monoxide solubility. J. Supercrit. Fluids 2013, 81, 23−32. (25) Al-Jimaz, A. S.; Fandary, M. S.; Al-Kandary, J. a.; Fahim, M. a. Phase equilibria for the extraction of sec-butylbenzene from dodecane with N,N-dimethylformamide. Fluid Phase Equilib. 2006, 240, 79−86. (26) Michelsen, M. L.; Kistenmacher, H. On composition-dependent interaction coefficeints. Fluid Phase Equilib. 1990, 58, 229−230. (27) Kuntsche, S.; Barz, T.; Kraus, R.; Arellano-Garcia, H.; Wozny, G. MOSAIC a web-based modeling environment for code generation. Comput. Chem. Eng. 2011, 35, 2257−2273. (28) gPROMS; Process Systems Enterprise: 1997−2015; www. pseenterprize.com/gproms. (29) Merchan, V. A.; Tolksdorf, G.; Kraus, R.; Wozny, G. Extending Documentation-Based Models towards an Efficient Integration into Commercial Process Simulators. Chem. Ing. Tech. 2014, 86, 1117− 1129. (30) Müller, E. A.; Gubbins, K. E. Molecular-Based Equations of State for Associating Fluids: A Review of SAFT and Related Approaches. Ind. Eng. Chem. Res. 2001, 40, 2193−2211. (31) Kontogeorgis, G. M.; Voutsas, E. C.; Yakoumis, I. V.; Tassios, D. P. An Equation of State for Associating Fluids. Ind. Eng. Chem. Res. 1996, 35, 4310−4318. (32) Chapman, W.; Gubbins, K.; Jackson, G.; Radosz, M. SAFT: Equation-of-state solution model for associating fluids. Fluid Phase Equilib. 1989, 52, 31−38. (33) Tumakaka, F. Modellierung von Phasengleichgewichten in Copolymer-Systemen mit polaren und assozierenden Monomeren (Modeling of Phase Equilibria in Polymer Systems and Association with Polar Monomers). Ph.D. Thesis, Universität Dortmund, 2004. R

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