Simulation Study of the Ethoxylation of Octanol in a Microstructured

File failed to load: https://cdn.mathjax.org/mathjax/contrib/a11y/accessibility-menu.js .... Publication Date (Web): November 16, 2016. Copyright © 2...
0 downloads 0 Views 6MB Size
Article pubs.acs.org/IECR

Simulation Study of the Ethoxylation of Octanol in a Microstructured Reactor Pascal D. Hermann,† Toine Cents,‡ Elias Klemm,† and Dirk Ziegenbalg*,† †

Institute of Chemical Technology, University of Stuttgart, 70569 Stuttgart, Germany Sasol Technology Netherlands B.V., 7544 GG Enschede, The Netherlands



S Supporting Information *

ABSTRACT: The influence of physical properties and the resulting limitations of heat and mass transport for the ethoxylation of octanol in a microstructured reactor were investigated by CFD simulations. The reaction was performed under supercritical conditions in a onephase system at pressures between 90 and 100 bar and temperatures between 180 and 240 °C. A 2D CFD model was applied to determine the kinetic parameters by fitting the model to experimental data. Furthermore, the sensitivity of the simulations regarding different physical properties was determined. The influence of diffusion was studied in detail, and it was found that a low diffusion coefficient results in a radial gradient in viscosity which leads to segregated streamlines and a bending of the streamlines toward the center. This “bottleneck” effect causes a large increase of the velocity in the center of the pipe, with the risk of the breakthrough of ethylene oxide. This effect occurs only in the case of higher ethoxylation degrees.

1. INTRODUCTION Surfactants made from ethoxylated fatty alcohols find widespread use in applications such as industrial detergents and processing aids.1 The world market for these ethoxylated fatty alcohols was valued to be $5.125 billion with a production of 3.27 million tons in 2014.2 The reaction of alcohol with ethylene oxide is typically promoted by alkaline catalysts such as KOH or NaOH. The process is usually performed in a semibatchwise operation by slowly adding gaseous ethylene oxide to a cooled stirred tank reactor loaded with the liquid alcohol/alcoholate substrate.3 However, this type of reactor can pose mass-transfer limitations as a result of the dissolution process of the gaseous ethylene oxide in the reaction mixture. This is slow compared to the rapid consumption of the dissolved ethylene oxide, which reduces the productivity of the process.3,4 Furthermore, the heat-transfer capability of this type of reactor is limited. Thus, it is necessary to slow the reaction to prevent a thermal runaway of this highly exothermic reaction. All together, the productivity in a stirred tank reactor is limited by both the mass-transfer limitations and the heat-removal capabilities of the vessel. To intensify the process and improve the safety, Rupp et al. previously transferred the process to a continuous operation in a microcapillary reactor.5−7 A microstructured reactor generally shows improved heat-transfer capabilities and better mechanical resistance to pressure, enabling operation in a one-phase system under elevated pressures and temperatures. Nonetheless, their experiments revealed problems controlling the reaction in the case of high ethoxylation degrees. In addition to the formation of hotspots, breakthrough of ethylene oxide occurs, leading to safety issues.7 © 2016 American Chemical Society

Hence, a better understanding of this breakthrough is crucial for the safe operation of such a process. Because of the submillimeter size of the channels in the reactor, local measurements of the liquid temperature are difficult. Therefore, numerical methods were chosen to gain deeper insight into the submillimeter scale. In previous works, Rupp et al.5−7 used a plug-flow reactor (PFR) model to simulate the product distribution and determine the kinetic parameters. However, this PFR model was not able to fully explain the experimental observations for long-chain products. The computational efforts for a PFR model are low, and such computations often allow for good predictions of reaction conversions and selectivities. Yet, a PFR model excludes the effects of radial and axial transport phenomena and, thus, the occurrence of radial gradients in concentration and physical parameters.8 However, it is possible to describe the radial residence time distribution by simplified methods such as the Taylor−Aris dispersion. Such models can explain a “breakthrough” of reactants in the center of reactors, but this exclude radial mass-transport effects.9 Computational fluid dynamic (CFD) simulations, with detailed descriptions of radial and axial transport phenomena, provide a powerful tool for gaining deeper insight into the roles of heat and mass transport. Because of the need for this deeper insight, a two-dimensional (2D) CFD model was used in this work to investigate the Received: Revised: Accepted: Published: 12675

October 24, 2016 November 11, 2016 November 16, 2016 November 16, 2016 DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research influence of radial heat and mass transport. The physical properties of the pure compounds were described by correlations based on experimental data, when available. Otherwise, predicted data based on equation-of-state models were used. The properties of the reaction mixture were calculated from the properties of the pure compounds by literature correlations. The kinetic parameters of the reactions were numerically fitted to experimental data. To test the robustness of the physical property correlations, sensitivity studies were performed. Additionally, the influence of diffusion was studied in detail, and the experimentally observed breakthrough of ethylene oxide was explained.

2. MODELING AND SIMULATION 2.1. Reactor Modeling. The simulations were performed in ANSYS CFX, release 15.0. The set of equations used were the unsteady Navier−Stokes equations in their conservative form. The following governing equations for mass, momentum, and energy conservation were used:10 Continuity: ρ

∂ρ + ∇(ρu) = 0 ∂t

Figure 1. 9° cylinder segment model.

=

n1 n product

MP1 + ... +

nj n product

MPj

ni MPi n product

(4)

j

(2)

∑ xiMP

M̅ P =

i

(5)

i=0

It is possible to introduce a pseudoproduct species Ps. The parameter s is the chain length of the highest simulated polymerized species, which averages the molar masses of all higher polymerized species. Inserting this species into eq 5 leads to

∂ρhtot ∂ρ −ρ + ∇(ρuhtot) = ∇(λ∇T ) + ∇(uτstress) + SE ∂t ∂t

s−1

(3)

∑ xiMP + xsMP

M̅ P =

In this equation, htot is the total enthalpy per mass unit; λ is the thermal conductivity; and SE is the reaction heat, calculated as SE = ΔHR

MP0 +

where M̅ P is the average molar mass of all products; n is the amount of product species; MPi is the molar mass of product species i; nproduct is the total amount of products, including octanol; and j is the highest polymerized species. MP0 is the molar mass of the initial monomer, octanol in this work. Substituting xi = ni/nproduct into eq 4 gives

(1)

2 The stress tensor is τstress = μ⎡⎣∇u + (∇u)T − 3 δ∇u⎤⎦, p is the pressure, and SM is the source of momentum. Energy:

ρ

∑ i=0

Momentum: ∂ρu + ∇(ρu ⊗ u) = −∇p + ∇τstress + SM ∂t

n product j

In this equation, ρ is the density, u is the velocity, and t is the time.

ρ

n0

M̅ P =

i

s

(6)

i=0

Because the sum of all mole fractions is unity, that is

1 dmEO MEO−1 V dt

j

1=

s−1

∑ xi = ∑ xi + xs i=0

with ΔHR = −95 kJ mol−1.5 The reactor model was a 9° cylinder segment with one cell width and symmetry boundaries on both sides, as shown in Figure 1. This 2D model was used with a varying length and radius. Because of the small diameter, a laminar flow regime was assumed. No-slip boundary conditions were set at the pipe wall. The feed was assumed to be ideally premixed. At the inlet, a fully developed laminar flow was used, and at the outlet, a zerogradient condition was set. This resulted in a 2D model that showed both good performance and precision and was able to resolve radial gradients. It is beneficial to reduce the number of simulated species because the calculation time scales with the number of simulated species. A large number of species can increase the calculation time significantly, especially for numerical optimizations such as the determination of the rate constants with a high number of iterations. The average molar mass of the product for a polymerization can be written as

i=0

this leads to s−1

M̅ P =

k−1

∑ xiMP + (1 − ∑ xi)MP i

k

i=0

i=0

(7)

For ∑s−1 i=0 xi < 1, it is possible to write this equation as M Pk = (M̅ P −

j−1

s−1

∑ xiMPi)(1 −

∑ xi)−1

i=0

i=0

(8)

This equation provides a general equation for the molar mass of a pseudoproduct species that is the sum of all higher polymerized products. It can be used for all polymerizations. For the ethoxylation of octanol, the average molar mass M̅ P of the product can be written as M̅ P = Moctanol + XEOeq EO MEO 0

12676

(9)

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

reaction mixture. Under these conditions, ethylene oxide is in a supercritical state, but no literature data for this supercritical mixture are available. Thus, a large effort was put into correlations for calculating the physical properties of the reaction mixture. The relevant physical properties are the diffusion coefficients of all species in the reaction mixture, the density, the dynamic viscosity, the thermal conductivity, and the heat capacity. All of these properties depend on the temperature, the composition of the reaction mixture, and the change with conversion. Additionally, there are some dependencies on pressure, especially for supercritical EO. To fit the properties of density ρ, thermal conductivity λ, heat capacity Cp, and viscosity η, experimental data of the pure compounds were used, when available. Otherwise, a method based on the Twu−Sim−Tassone equation of state in Aspen Properties was used (HYSGLYCO method).13 The density of the mixture was additionally corrected with experimental data measured at our institute. A more detailed description of the individual correlations and the sources of the data can be found in the Supporting Information. After testing various calculation methods for the diffusion coefficient with a large number of similar results (see Figure 4), the following equation was selected due to its good performance: For VEO/VP ≤ 1 or 2, where VEO is the molar volume of ethylene oxide and VP is the average molar volume of octanol and the product, the following simplified relation is given by Scheibel14

In this equation, eqEO0 is the number of EO equivalents per mole of octanol at the inlet of the reactor, and XEO is the conversion of ethylene oxide. 2.2. Kinetic Modeling. The ethoxylation of octanol is an anionic polymerization. In various publications, it is assumed that only the anionic species is active for the chain-growth reaction. Additionally, there is a proton-transfer reaction that leads to chain growth at a different chain.3,5,11,12 In this work, it is assumed that proton transfer is very fast and, thus, that no limitations on proton transfer exist. A secondorder reaction for each step is assumed with different rate constants for the first two steps (k1, k2) and a constant rate constant for the third and all subsequent steps (k3)

To adapt the concentration change in the CFD model, it is necessary to transfer the change in concentration into equations for mass change: dmPi dt

= M Pi

dci dt

(10)

For a second-order reaction, the change in concentration can be written as

DEO = 17.510 × 10−15

dci = k p c1c 2 i dt

dt

= −k1M P0ρ2

wP0 wEO M P0 MEO

(11)

In this equation, P0 represents octanol, P1 is the product with a chain length of one, EO is ethylene oxide, M is the molar mass, w is the mass fraction, and k1 is the rate constant of the first step. For 1 ≤ i < s, where s is the selected highest simulated species and is thus defined as the pseudoproduct species, this leads for each species i to the following rate of mass change: dmPi dt

= kiM Piρ2

wP w wPi−1 wEO − ki + 1M Piρ2 i EO M Pi MEO M Pi−1 MEO

(12)

For the pseudoproduct species, this results in dmPs dt

= ksM Psρ2

wP w wPs−1 wEO + ksMEOρ2 s EO M Ps MEO M Ps−1 MEO

(13)

with Ms being the molar mass of the pseudoproduct species, calculated according to eq 8. The rate of mass change for ethylene oxide is: s ⎛ dmEO w⎞ = −MEOρ2 ∑ ⎜ki i ⎟ dt Mi ⎠ i=0 ⎝

(15)

where DEO is the diffusion coefficient of ethylene oxide in m2 s−1, T is the temperature in K, η is the dynamic viscosity of the solvent in Pa s, and VPbp is the molar volume of the solute at the boiling point in cm3 mol−1. The molar volume at the boiling point was estimated with a correlation given by Wilke and Chang.15 Unless otherwise stated, the diffusion coefficient of the product was assumed to be identical to that of ethylene oxide. This simplification was necessary to reduce the computational time and improve the numerical stability. 2.4. Simulations. Three different geometries were used in this work: The experimental data used for the kinetic experiments were measured in a 4-m capillary reactor with an inner diameter of 250 μm, and these dimensions were consequently used for the 2D model employed during the estimation of the kinetic parameters.7 For the sensitivity study, a simulated reactor length of 6.2 m with an inner radius of 250 μm was used. Additionally, a simulated reactor length of 1 m with the same radius was used for the diffusion study. Except for the density, the physical properties of the reaction mixture were calculated from the properties of the pure compounds, rather than from experimental data. Thus, it was very important to test the sensitivity of the simulation with regard to the correct description of the physical properties. This was done by multiplying the mixture properties by a factor f. To test the effects of thermal conductivity and heat capacity as well, a larger inner diameter of 500 μm for the simulated reactor and a constant wall temperature were used. The simulated reactor length was 6.2 m. The number of simulated product species s was set to 14; thus, the product with a chain length of 14 was used as the pseudoproduct species. The factors

For octanol, this results with ci = ρwiMi−1 and eq 10 in the following equation for the rate of mass change: dmP0

T N cm ηVPbp1/3 K mol1/3

(14)

2.3. Physical Properties. The chosen reaction conditions, namely, high pressures and temperatures (high-p and -T), raise challenges to the description of the physical properties of the 12677

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

where M̅ w is the weight-average molar mass and M̅ n the number-average molar mass. The conversion of EO, XEO, was calculated with the equation

were changed between 0.01 and 100 for one property by keeping the factors for the others at 1. The temperature used for the sensitivity study was 220 °C, the amount of ethylene oxide used was 9 equivalents per mole of octanol, and the mass flow rate at the inlet was kept constant at 0.37 kg h−1. These conditions were chosen to ensure a conversion of EO between 0.8 and 0.85, which showed the largest radial gradients in conversion and thus were used to investigate the maximum sensitivity. In this work, the focus was on the conversion of ethylene oxide, the polydispersity index (PDI), the maximum radial temperature gradient, and the average temperature in the reactor. The polydispersity index or, more correctly, the molar mass dispersity was used for comparison because of its good usability as an index for the breadth of the product distribution. The PDI was calculated as: s

PDI =

XEO = 1 −

(16)

Table 1. Arrhenius Plot Parameters K1 error K2 error K3 error

EA (kJ mol−1)

ln(k0) (m3 mol−1 s−1)

71.1 1.8 66.0 2.8 69.5 7.3

11.02 0.44 10.86 0.67 12.03 1.75

(17)

where wEO is the mass fraction of EO at the end of the reactor and wEOstart is the mass fraction of EO at the inlet. Both PDI and XEO were averaged by the local mass flow at the outlet of the reactor. The physical properties tested were the density, dynamic viscosity, thermal conductivity, heat capacity, and diffusion coefficient. The simulations were stopped after reaching an imbalance in mass, velocity, momentum, and energy below 1% and a residual below 10−4 between two successive iteration steps. To ensure numerical stability, it was necessary to set a maximum temperature of 600 K. All temperatures higher than 600 K were automatically set to 600 K. Because of the laminar flow regime, diffusion was assumed to be the only radial mass-transport mechanism in a straight microtube reactor. Thus, it was important to additionally test the influence of diffusion on the reaction. The simulated reactor for the diffusion study had a radius of 250 μm and a length of 1 m. A space time of 60 s and a constant wall temperature of 220 °C were employed. To test the influence of the magnitude of diffusion on the reaction and its relevance for performing the reaction in a microreactor with its generally smaller diffusion pathways, a constant diffusion coefficient was assumed.

s

M̅ w −1 = (∑ wM i i)(∑ xiMi) M̅ n i=0 i=0

wEO wEOstart

Figure 2. Parity plots for (a) octanol, (b) product 1, (c) product 2, and (d) pseudoproduct species for all higher ethoxylated products. The black dotted lines gives the ±10% error. The capillary reactor used for the experiments and the 2D model used for the simulations had a length of 4 m with an inner radius of 125 μm. The temperatures were 180, 200, 220, and 240 °C with 3, 6, and 9 equivalents of EO per mole of octanol and a concentration of 0.66 mol % octanolate in octanol. The experiments were performed by Rupp.7 12678

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 3. Sensitivity analysis of physical properties: (a) conversion of EO, (b) molar mass dispersity, (c) maximum temperature, and (d) average temperature. f gives the factor that was used to scale the property. The geometry used was the 2D model with a reactor length of 6.2 m and an inner radius of 250 μm. The mass flow rate used was 0.37 kg h−1, the temperature was set to 220 °C, and 9 equivalents of EO was used per mole of octanol.

The tested diffusion coefficients were varied between 10−8 m2 s−1, which is in the typical magnitude for a supercritical fluid, and 10−12 m2 s−1, which is the typical magnitude for a polymer. This study was performed to test the influence of the diffusion on the reaction in general. Additionally, the study was performed to evaluate whether it is necessary to make the model more complex by employing a conversion- and temperature-dependent function for calculating the diffusion coefficient or whether it is sufficient to assume a constant diffusion coefficient, which significantly reduces the computational effort. Additionally, some published correlations for calculating the diffusion coefficient were tested; these were the ones published by Wilke and Chang,16 Einstein et al.,17 and Reddy and Doraiswamy,18 as well as two different correlations from Scheibel et al.15 These correlations calculate the local diffusion coefficient, which is influenced by temperature and the local average chain length. Except for the simple Scheibel model, the diffusion coefficient of the product is calculated separately. In the simple Scheibel model, it is assumed that the diffussion coefficient of the product is identical to the diffusion coefficient of ethylene oxide.

performed with an inner diameter of 250 μm.7 In this work, the kinetic parameters were numerically determined by utilizing the CFD model. This was done by minimizing the deviation of the simulated data from the experimental data. More precisely, the sum of the mean absolute difference, SMD, was minimized with a Python-based script using the Nelder−Mead algorithm.19 The sum of the mean absolute difference, SMD, was calculated as s

SMD =

⎡1

l

∑ ⎢ ∑ (|xPsim i=0

⎢⎣ l

i

n=1

n

⎤ ⎥ − x Pexp | ) i n ⎥ ⎦

(18)

with

x Pi =

nPi s ∑ j = 0 nPj

where l is the number of experiments, xPisimn is the simulated mole fraction of product i after removal of EO for experiment n and species i, and xPiexpn is the experimental mole fraction. The model for the kinetic determination used s = 3; thus, the 3-foldethoxylated product is a pseudospecies. The activation energies and preexponetial factors obtained from the Arrhenius plot can be found in Table 1. The deviations from the experimental values for the mole fractions of the four species used are shown in Figure 2. Considering the simplifications of the reactor model and the assumptions made concerning the physical properties, the simulations showed a remarkable precision in the description of the conversion and product distribution. The prediction accuracy of the resulting model can be given by the average standard

3. RESULTS AND DISSCUSSION 3.1. Estimation of Kinetic Parameters. The experimental data were measured by Rupp et al. in a laboratory-scale reactor.7 The experiments were performed at four temperatures (180, 200, 220, 240 °C) with 3, 6, and 9 equivalents of EO and a concentration of 0.66 mol % octanolate in octanol. The reactions were conducted in a tempered capillary reactor with a length of 4 m. To ensure isothermal behavior, the experiments were 12679

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 4. Influence of the diffusion coefficient for a constant diffusion coefficient and different published correlations for calculating the diffusion coefficient: (a) conversion of EO, (b) molar mass dispersity PDI, and (c) average velocity. The geometry used was the 2D model with a reactor length of 6.2 m and a radius of 250 μm. The mass flow rate used was 0.37 kg h−1, the temperature was set to 220 °C, and 9 equivalents of EO was used per mole of octanol.

of the laminar flow, with higher velocities in the center and lower velocities near the wall. A larger diffusion coefficient ( f > 1.1) and, thus, faster diffusion leads to the opposite effect, increasing conversion and decreasing the molar mass dispersity. The density ρ has a surprisingly low influence on the conversion of ethylene oxide, considering that the space time is mainly a function of the density. For instance, a factor of 0.1, corresponding to a 10 times lower space time, will decrease the conversion by only 7%. Additionally (not plotted), a lower density leads to a larger radial conversion gradient. A high density has only a small influence on conversion, with a factor of 10 increase slightly reducing the conversion. The molar mass dispersity will increase for small densities and decrease for high densities, similar to the influence of the diffusion coefficient. This is due to the fact that low densities increase the axial velocity (constant mass flow) and, thus, the time for radial diffusion is limited leading to a broader product distribution. In contrast, a higher density leads to a lower axial velocity, with more time for radial diffusion and thus a narrower product distribution and a lower PDI. The density has a negligible influence on the temperature of the simulation. The thermal conductivity λ exhibited the largest influence on conversion. Low values led to full conversion because of the higher temperature of the reaction mixture. The average temperature increased significantly for factors f of 0.5 and smaller. For factors of 0.25 and lower, the maximum temperature reached the simulation limit with a radial temperature gradient of 107 K, as shown in Figure 3c. The radial heat removal by heat conduction was calculated with the thermal conductivity.

deviation from the experimental data. The resulting standard deviation in the mole fractions of the individual product species was 0.017. This corresponds to a relative standard deviation of 9.8%. The experimental error calculated from 10 different reproduction experiments corresponds to a relative standard deviation of 8.6% for the individual species. Thus, the prediction accuracy of the simulation is very close to the experimental error. The good agreement in the product distribution between the simulations and the experimental data indicates that the descriptions of both the physical properties and the kinetics are correct at temperatures between 180 and 240 °C and pressures in the range of 90−100 bar. 3.2. Influence of Physical Properties. The results of the sensitivity studies in terms of the factor f are shown in Figure 3, with f = 1 corresponding to the unaltered correlation as listed in the Supporting Information. The influence of the dynamic viscosity η is negligible for the complete range of factors f under the chosen conditions. However, influences for an increasing dynamic viscosity can be expected for higher ethoxylation degrees because of its influence on the diffusion, as shown in section 3.3. The diffusion coefficient is reciprocally coupled to the dynamic viscosity. Thus, a higher viscosity will decrease the rate of diffusion. Small diffusion coefficients D ( f < 0.9) decrease the conversion and increase the molar mass dispersity. This results from the fact that diffusion is the only means of radial mass transport in a strictly laminar flow regime. A lower radial masstransport rate will increase radial gradients. The broader molar mass dispersity is due to the broad residence time distribution 12680

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 5. Radial streamline position for an unrealistic constant diffusion coefficient: (a) diffusion coefficient of 10−8 m 2 s−1, (b) diffusion coefficient of 10−9 m2 s−1, (c) diffusion coefficient of 10−10 m2 s−1, and (d) diffusion coefficient of 10−11 m2 s −1. The geometry used was the 2D model with a reactor length of 1 m and a radius of 250 μm. The chosen space time was 120 s, the temperature was set to 220 °C, and 9 equivalents of EO was used per mole of octanol.

In a strictly laminar flow regime, this is the sole mechanism of radial heat removal. A higher thermal conductivity will lead to a better removal of the reaction heat and, thus, lower radial temperature gradients. As expected, a higher thermal conductivity decreases the conversion because of the lower temperature caused by the improved heat-removal capability. The heat capacity, Cp, has only a limited influence on conversion, PDI, and average temperature. However, a factor lower than 0.9 will cause a high local temperature maximum, reaching the simulation’s maximum radial temperature gradient of 107 K (maximum absolute temperature of 600 K) for factors of 0.25 and smaller. This occurs because of the lower capability of storing thermal energy in the fluid. In summary, heat removal has a strong influence, and thus thermal conductivity is an essential parameter, even in a microreactor. Radial mass-transport limitations seem to have only a limited influence under the chosen conditions for the sensitivity study. However, simulations with more equivalents of EO showed a growing influence for increasing the conversion and chain length, as discussed in the next subsection. Generally, the effects of all physical properties seem to be insignificant between factors of 0.9 and 1.1. Thus, the description of the physical properties shows a good robustness.

3.3. Influence of Diffusion. 3.3.1. Influence of the Magnitude of the Diffusion Coefficient. As shown in Figure 4, a decreasing diffusion coefficient lowers the conversion of ethylene oxide. According to expectations, the PDI increases with decreasing diffusion coefficient, because very small diffusion coefficients prevent radial mixing and cause a segregation of streamlines. In this case, the laminar flow regime corresponds to short residence times in the center and long residence times near the pipe wall. An unexpected effect is the increased mean velocity for lower constant diffusion coefficients. Thus, it is necessary to include a more detailed description of the diffusion coefficient. All tested correlations for the diffusion showed results similar to those obtained with a constant coefficient D between 1 × 10−9 and 10−8 m2 s−1 under the chosen reaction conditions. The study of the diffusion coefficient was repeated with a radius of 250 μm and a length of 1 m at a space time of 120 s, exceeding a conversion of 0.95 after 0.5 m for all tested diffusion correlations. As shown in Figure 5, lowering the diffusion coefficient led to a bending of the streamlines toward the center for diffusion coefficients of 10−9 m2 s−1 and lower. 3.3.2. Reason for Streamline Bending. As shown in Figure 6, using a diffusion coefficient of 10−8 m2 s−1 led to almost no radial gradients in conversion. A diffusion coefficient of 10−11 m2 s−1 led to a significant radial gradient in the 12681

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 6. Differences in properties on the streamline between an unrealistic constant diffusion coefficient of 10−8 m2 s−1 (left) and an unrealistic constant diffusion coefficient of 10−11 m2 s−1 (right). Red are the values on the streamlines near the pipe center, and blue are the values on the streamlines near the pipe wall. (a,b) Conversion on the streamlines, (c,d) velocity on the streamlines, and (e,f) dynamic viscosity on the streamlines. The geometry used was the 2D model with a reactor length of 1 m and a radius of 250 μm. The chosen space time was 120 s, the temperature was set to 220 °C, and 9 equivalents of EO was used per mole of octanol.

the run with a fixed dynamic viscosity of 0.2 mPa s are shown in Figure 7a. A fixed dynamic viscosity leads to the opposite behavior, namely, a slight streamline bending toward the pipe wall. This can be explained by the radial mass-transport limitation for this diffusion coefficient and thus a longer residence time, higher conversion, and higher density near the pipe wall. A higher density leads to a lower specific volume compared to the pipe center and thus a slight flow toward the wall from the center. A fixed density of 800 kg m−3, as shown in Figure 7b, leads to a streamline bending toward the center. These results provide

conversion of EO, with a high conversion at the pipe wall (blue) and a low conversion in the pipe center (red). Additionally, an increased velocity at the center streamlines was observed. This led to a significant shift of the mass flow toward the center of the pipe, causing the observed bending of the streamlines. This “bottleneck” effect induces an increased average velocity in the reactor and causes a breakthrough of EO. To clarify which physical property is responsible for this bottleneck effect, two isothermal runs were performed, so that temperature-related effects could be excluded. The diffusion coefficient was kept constant at 10−10 m2 s−1. The results from 12682

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 7. Radial streamline position under isothermal conditions for a constant diffusion coefficient of 10−10 m2 s−1 with a (a) fixed viscosity of 0.2 mPa s and (b) fixed density of 800 kg m−3. The geometry used was the 2D model with a reactor length of 1 m and a radius of 250 μm. The chosen space time was 120 s, the temperature was set to 220 °C, and 9 equivalents of EO was used per mole of octanol.

and the diameter of the reactor. Long space times generally reduce the influence of these effects, because the lower average velocity increases the time for radial diffusion. A larger radius will increase the mass-transport limitation effects because of the longer diffusion pathways. Finally, the red region indicates conditions with very slow diffusion and, thus, significant radial mass-transport limitations. As shown, this effect will become relevant for EO-to-octanol ratios above 10 and conversions above approximately 0.7. To test the influence of mass-transport limitations and the bottleneck effect on the results of the simulations, additional runs with 100 equivalents of EO per mole of octanol were performed with a simulated reactor length of 6.2 m and a radius of 500 μm. The pressure was set to 100 bar, and a constant temperature of 220 °C was set at the wall, with a space time of 45 s and 5 mol % octanolate per mole of octanol. The results are depicted in Figure 9. As predicted, the simulation showed streamline bending behavior (Figure 9a). The velocity first decreased with growing chain length but again increased in the center at 4 m and exhibited a maximum at 5 m (Figure 9b). In the first 3 m of the reactor, nearly no radial gradient in the average chain length and thus the dynamic viscosity is visible, as shown in Figure 9c,d. Generally, the radial gradient increases with the reactor length. As shown in Figure 9e, the density shows no drop that could induce the increased center velocity. The increased center temperature (see Figure 9f) counters this effect through the increased reaction rate in the center. The diffusion coefficient decreases with increasing chain length and reaches the range of mass-transfer limitations as determined in the sensitivity study. This is depicted in Figure 10a,b. Thus, the reason for the streamline bending is a low diffusion coefficient that leads to the segregation of the streamlines and thus large radial gradients. This is additionally enforced by the increased center velocity, leading to a bottleneck effect. However, this effect is not unique to ethoxylations. Polymers commonly have an increased viscosity for longer chain lengths, and thus, polymerization reactions will be affected by masstransport limitations. If the polymerization is performed in a tube-shaped reactor, it will be hard to achieve full conversion especially for long-chain products. It is known that the diffusion rate decreases during polymerization. Parida et al.20 stated that the increased viscosity slows the diffusion during polymerization in a microreactor and results in poor reaction control for

Figure 8. Diffusion coefficient DEO as a function of conversion and equivalents of EO depending on the diffusion coefficient at 220 °C and 100 bar. The green area is not limited by radial mass transport. In the orange region, mass-transport limitations depend on space time and radius. In the red area, radial mass transport is limited.

evidence that the bottleneck effect is caused by the gradient in the dynamic viscosity. A limitation in radial mass transport because of a small diffusion coefficient leads to a higher conversion near the pipe wall, resulting in a higher viscosity near the pipe wall. In the laminar flow regime, wall friction is proportional to dynamic viscosity, and thus, a higher viscosity increases the friction at the pipe wall. Additionally, slow radial mass transport results in a lower conversion at the center and, thus, a lower flow resistance at the center because of the higher viscosity at the wall. In sum, these effects cause the bending of the streamlines toward the center and hence to an increased center velocity. As a consequence, the residence time in the center decreases, increasing the risk of breakthrough of ethylene oxide. 3.3.3. Relevance for Ethoxylation. Thus, it was important to test the relevance of these effects on the ethoxylation with different EO-to-octanol ratios with a full conversion- and temperature-dependent description of the diffusion coefficient. The easiest method was to check the calculated diffusion coefficient of EO as shown in Figure 8 for 220 °C at 100 bar. The green and yellow regions in this figure indicate the conditions without relevant radial mass-transport limitations. The orange region indicates the transition region, where the effects of radial transport limitations depend on the space time 12683

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 9. Results for a space time of 45 s, a wall temperature of 220 °C, an EO content of 100 equivalents, and an octanolate loading of 5 mol % per mole of octanol. (a) Radial position of the streamlines; (b) velocity on the streamlines, with red being the center streamlines and blue near the pipe wall in panel a; (c) average chain length eg of the product; (d) dynamic viscosity η; (e) density ρ; and (f) temperature T.

long-chain products. Parida et al.20 also showed experimentally that improved mixing in a microreactor with a coil-flow inverter can reduce the PDI. Both effects are in line with the results presented herein. Cortese et al.21 also performed CFD simulations of an anionic polymerization in a tube-shaped reactor. The authors also stated that the increased viscosity of the reaction mixture slows the diffusion, leading to a segregation of the streamlines. Cortese et al.21 also observed an increased velocity on the center streamlines and stated the radial viscosity gradient as the reason. Although they did not mention the radial shift of the mass flow by streamline bending, the increased

center velocity validates a similar effect as observed in this work. Schulze et al.22 showed that using a model of a laminar flow reactor that is comparable to completely segregated streamlines leads to a much lower PDI than experimentally observed. All together, these publications state that performing a polymerization with high-chain-length products is hard to control in a tube-shaped microreactor. The bottleneck effect observed in this work for low diffusion coefficients is one reason. The other reason is the segregation of the streamlines due to masstransport limitations and thus the residence time distribution in a laminar flow regime. Both the segregation of the streamlines 12684

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research

Figure 10. Results for a space time of 45 s, a wall temperature of 220 °C, an EO content of 100 equivalents, and an octanolate loading of 5 mol % per mole of octanol. (a) Average chain length eg of the product and (b) calculated diffusion coefficient DEO for EO.



and the mass flow shift by streamline bending with increased center velocity explain why radial mixing leads to much lower PDIs, as observed in several publications.20,23 Radial mixing improves the radial mass transport and thus decreases the radial gradients.

Corresponding Author

*E-mail: [email protected]. ORCID

Dirk Ziegenbalg: 0000-0001-6104-4953 Notes

4. CONCLUSIONS

The authors declare no competing financial interest.



The influence of changes in physical properties on the results of simulations of ethoxylation is crucial. The increasing viscosity will induce a limitation on the radial mass transport, increase the molar mass dispersity, and reduce the degree of conversion. This effect is amplified by a possible bending of streamlines with the resulting bottleneck effect and the occurrence of hotspots. Both can lead to a breakthrough of ethylene oxide that was experimentally observed for long-chain products. This effect will also be observed for other homogeneous polymerizations. The supercritical conditions at the inlet of the reactor with low viscosity improve the radial mass transport by diffusion. Thus, the influence of radial mass-transport limitations for this specific reaction can be expected to be low compared to other polymerizations. It can be concluded that achieving high productivities with a homogeneous polymerization in a tube-shaped microreactor will be difficult for long-chain products. The increasing viscosity will slow radial mass transport by diffusion, with the possible coupled effect of streamline bending preventing full conversion and increasing the molar mass dispersion of the product. It is necessary to reduce the radial gradients in the area of masstransport limitations to minimize these effects and prevent a bottleneck effect. There are different possibilities for achieving this. One is to reduce the diameter or change the shape of the reactor with the objective of achieving smaller diffusion pathways at the cost of a higher pressure drop and investment costs. Another is to reduce the mass flow or dilute the reaction mixture with a solvent at the expense of productivity.



AUTHOR INFORMATION



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.iecr.6b04110.

SYMBOLS AND ABBREVIATIONS CFD = computational fluid dynamic eg = ethoxylation grade EO = ethylene oxide eq = equivalents c = concentration (mol m−3) Cp = heat capacity at constant pressure (J kg−1 K−1) D = diffusion coefficient (m2 s−1) i = chain length i j = highest ethoxylated product in the reaction mixture k = reaction rate (mol s−1 m−3) l = number of experiments used ṁ = mass flow rate (kg h−1) M = molar mass (g mol−1) n = mol amount of substance Pi = product with chain length i PDI = molar mass dispersity PFR = plug-flow reactor s = defined highest simulated species, averages s times and higher ethoxylated products SMD = sum of average mean difference u = velocity (m s−1) w = mass fraction x = mole fraction XEO = conversion of ethylene oxide λ = thermal conductivity (W K−1 s−1) η = dynamic viscosity (Pa s) ρ = density (kg m−3) Δ = variable difference ∇ = nabla operator REFERENCES

(1) Evetts, S.; Kovalski, C.; Levin, M.; Stafford, M. High-temperature stability of alcohol ethoxylates. J. Am. Oil Chem. Soc. 1995, 72, 811− 816. (2) Alcohol Ethoxylate Market by Application (Household & Personal Care, Pharmaceutical, Oilfield, Agrochemicals, and Others) & Geog-

Description of the correlations and sources of the data for the physical properties (PDF) 12685

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686

Article

Industrial & Engineering Chemistry Research raphiesGlobal Market Trends, Forecasts to 2019. 2015; available at http://www.researchandmarkets.com/reports/3292449/alcoholethoxylate-market-by-application (accessed 16.02.2016). (3) Di Serio, M.; Tesser, R.; Santacesaria, E. Comparison of Different Reactor Types Used in the Manufacture of Ethoxylated, Propoxylated Products. Ind. Eng. Chem. Res. 2005, 44, 9482−9489. (4) Santacesaria, E.; Diserio, M.; Tesser, R. Role of ethylene oxide solubility in the ethoxylation processes. Catal. Today 1995, 24, 23−28. Catalysis in Multiphase Reactors. (5) Rupp, M.; Ruback, W.; Klemm, E. Octanol ethoxylation in microchannels. Chem. Eng. Process. 2013, 74, 19−26. (6) Rupp, M.; Ruback, W.; Klemm, E. Alcohol ethoxylation kinetics: Proton transfer influence on product distribution in microchannels. Chem. Eng. Process. 2013, 74, 187−192. (7) Rupp, M. Ueber die Ethoxylierung von Octanol im Mikrostrukturreaktor. Ph.D. Thesis, Universität Stuttgart, Stuttgart, Germany, 2015. (8) Zhai, X.; Ding, S.; Cheng, Y.; Jin, Y.; Cheng, Y. CFD simulation with detailed chemistry of steam reforming of methane for hydrogen production in an integrated micro-reactor. Int. J. Hydrogen Energy 2010, 35, 5383−5392. (9) Baldauf, W.; Knapp, H. Measurements of diffusivities in liquids by the dispersion method. Chem. Eng. Sci. 1983, 38, 1031−1037. (10) ANSYS CFX-Solver Theory Guide; ANSYS, Inc.: Canonsburg, PA, 2013. (11) Gee, G.; Higginson, W. C. E.; Levesley, P.; Taylor, K. J. 266. Polymerisation of epoxides. Part I. Some kinetic aspects of the addition of alcohols to epoxides catalysed by sodium alkoxides. J. Chem. Soc. 1959, 1338. (12) Chiu, Y.; Naser, J.; Easton, A.; Ngian, K.; Pratt, K. Kinetics of a catalyzed semi-batch ethoxylation of nonylphenol. Chem. Eng. Sci. 2010, 65, 1167−1172. (13) Aspen Properties Reference Manual, version 7.3; AspenTech: Cambridge, MA, 2011. (14) Hayduk, W.; Laudie, H. Prediction of diffusion coefficients for nonelectrolytes in dilute aqueous solutions. AIChE J. 1974, 20, 611− 615. (15) Logan, B. E. Environmental Transport Processes, 1st ed.; Wiley: New York, 1999; p 67. (16) Wilke, C. R.; Chang, P. Correlation of diffusion coefficients in dilute solutions. AIChE J. 1955, 1, 264−270. (17) Cussler, E. L. Diffusion: Mass Transfer in Fluid Systems, 3rd ed.; Cambridge University Press: Cambridge, U.K., 2010. (18) Reddy, K. A.; Doraiswamy, L. K. Estimating Liquid Diffusivity. Ind. Eng. Chem. Fundam. 1967, 6, 77−79. (19) Nelder, J. A.; Mead, R. A Simplex Method for Function Minimization. Computer Journal 1965, 7, 308−313. (20) Parida, D.; Serra, C. A.; Garg, D. K.; Hoarau, Y.; Bally, F.; Muller, R.; Bouquey, M. Coil Flow Inversion as a Route To Control Polymerization in Microreactors. Macromolecules 2014, 47, 3282− 3287. (21) Cortese, B.; Noel, T.; de Croon, M. H. J. M.; Schulze, S.; Klemm, E.; Hessel, V. Modeling of Anionic Polymerization in Flow With Coupled Variations of Concentration, Viscosity, and Diffusivity. Macromol. React. Eng. 2012, 6, 507−515. (22) Schulze, S.; Cortese, B.; Rupp, M.; de Croon, M. H.; Hessel, V.; Couet, J.; Lang, J.; Klemm, E. Investigations on the anionic polymerization of butadiene in capillaries by kinetic measurements and reactor simulation. Green Process. Synth. 2013, 2, 381−395. (23) Ziegenbalg, D.; Kompter, C.; Schoenfeld, F.; Kralisch, D. Evaluation of different micromixers by CFD simulations for the anionic polymerisation of styrene. Green Process. Synth. 2012, 1, 211− 214.

12686

DOI: 10.1021/acs.iecr.6b04110 Ind. Eng. Chem. Res. 2016, 55, 12675−12686