- Email: [email protected]
Numerical simulation of mass transfer in a liquid–liquid membrane contactor for laminar flow conditions
Available online at www.sciencedirect.com
Computers and Chemical Engineering 32 (2008) 1325–1333
Numerical simulation of mass transfer in a liquid–liquid membrane contactor for laminar flow conditions Roland Kieffer, Catherine Charcosset ∗ , Franc¸ois Puel, Denis Mangin Laboratoire d’Automatique et de G´enie des Proc´ed´es, UMR CNRS 5007, UCBLyon 1, ESCPE-Lyon, 43 Bd du 11 Novembre 1918, 69 622 Villeurbanne Cedex, France Received 6 November 2006; received in revised form 31 May 2007; accepted 6 June 2007 Available online 20 June 2007
Abstract Liquid–liquid phase membrane contactors are increasingly being used for mixing and reaction. The principle is the following: component A flows through the membrane device inlet to mix/react with component B which comes from the membrane pores. This study presents a numerical simulation using computational fluid dynamics (CFD) of momentum and mass transfer in a tubular membrane contactor for laminar flow conditions. The velocity and concentration profiles of components A–C are obtained by resolution of the Navier-Stokes and convection-diffusion equations. The numerical simulations show that mixing between A and B is obtained by diffusion along the streamlines separating both components. The mixing/reaction zone width is within the region of a few hundred of microns, and depends on the diffusion coefficients of A and B. Hollow fiber membrane devices are found to be of particular interest because their inner diameter is close to the mixing zone width. © 2007 Elsevier Ltd. All rights reserved. Keywords: Computational fluid dynamics; Membrane contactor; Membrane reactor; Micromixing; Microreactor; Mixing
1. Introduction Membrane contactors represent an emerging technology in which the membrane is used as a tool for inter-phase mass transfer operations. These membrane systems provide a high interfacial area between two phases to achieve high overall rates of mass transfer (Drioli, Criscuoli, & Curcio, 2003; Sirkar, Shanbhag, & Kovvali, 1999). Membrane contactors have attracted considerable attention recently for liquid–gas phase processes such as simultaneous extraction and injection of gases in natural water production, and membrane distillation for desalination and wastewater treatment (Drioli et al., 2003). Liquid–liquid membrane contactors are reported for the preparation of nanoprecipitates (Chen, Luo, Xu, & Wang, 2004; Fei, Jia, Yuelian, & Liu, 2003; Jia & Liu, 2002; Jia, Liu, & He, 2003), polyaluminium chloride particles (Fei et al., 2003), l-asparagine crystals (Zarkadas & Sirkar, 2006) and polymeric nanoparticles (Charcosset & Fessi, 2005).
∗
Corresponding author. Tel.: +33 4 72 43 18 67; fax: +33 4 72 43 16 99. E-mail address: [email protected] (C. Charcosset).
0098-1354/$ – see front matter © 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.compchemeng.2007.06.013
The principle of liquid–liquid membrane contactors is the following. Component A flows through the inlet of the lumen side of the membrane contactor and component B comes from the membrane pores (shell side). Components A and B mix and react inside the lumen side of the membrane device. This membrane process (for a hollow fiber membrane device) is an alternative to micromixing and microchannel reactors. Standard mixing devices include the stirred tank reactor, the T-tube or Y-tube, and more sophisticated ones, like the two-impinging-jets mixing device, the sliding surface mixing device and the vortex reactor (B´enet, Muhr, Plasari, & Rousseaux, 2002). The common feature of micromixing and microreactor designs is small channels with dimensions below 1 mm. Application fields of micromixers include both modern, specialised issues such as sample preparation for chemical analysis and traditional, widespread usable mixing tasks such as reaction, gas absorption, emulsification, foaming and blending (Hessel, L¨owe, & Sch¨onfeld, 2005). CFD has been largely used as a tool to model membrane separation processes (Belfort, 1989; Ghidossi, Veyret, & Moulin, 2006). The Navier-Stokes equations were solved to estimate the pressure drop across the membrane module and the shear rate on the membrane surface. Several studies used the concentration polarization model, and other studies were carried
1326
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
Nomenclature CA CB CA,in CB,in D Di k L Ni Ni,in Ni,out p Qin,L Qin,S Qout r R RC Reout v v vin,L vin,S vout Xi z
concentration of component A (mol/m3 ) concentration of component B (mol/m3 ) concentration of component A at the inlet of the lumen side (mol/m3 ) concentration of component B at the inlet of the shell side (mol/m3 ) inner diameter (m) diffusion coefficient of component i (i = A–C) (m2 /s) reaction rate constant (m3 /mol s) membrane length (m) molecular flux of component i (i = A and B) (mol/m2 s) total molecular flux of component i calculated at the inlets (mol/s) total molecular flux of component i calculated at the outlet of the lumen side (mol/s) pressure (Pa) flow rate at the inlet of the lumen side (m3 /s) flow rate at the inlet of the shell side (m3 /s) flow rate at the outlet of the lumen side (m3 /s) radial coordinate (m) inner radius (m) second-order reaction term for the preparation of C (mol/m3 s) Reynolds number calculated at the outlet of the lumen side velocity vector (m/s) norm of the velocity vector (m/s) mean velocity at the inlet of the lumen side (m/s) velocity at the inlet of the shell side (m/s) mean velocity at the outlet of the lumen side (m/s) conversion rate of component i (i = A and B) axial coordinate (m)
Greek letters η dynamic viscosity (Pa s) ρ fluid density (kg/m3 )
out inside the membrane porous matrix. Numerical simulations were performed for various membrane filtration processes: microfiltration and ultrafiltration (Ghidossi et al., 2006; Pellerin, Michelitsch, Darcovich, Lin, & Tam, 1995; Wiley & Fletcher, 2003), nanofiltration and reverse osmosis (i.e. Geraldes, Semi˜ao, & Pinho, 2000). Various membrane module configurations were also investigated: flat sheet (Geraldes et al., 2000; Pellerin et al., 1995; Wiley & Fletcher, 2003), rotating circular membrane geometry (Serra & Wiesner, 2000), tubular membranes with inserts (Bellhouse, Costigan, Abhinava, & Merry, 2001) and centrifugal membrane separation (Pharoah, Djilali, & Vickers, 2000). CFD was also used to understand the decrease in membrane fouling by using different instabilities or turbulence promoters, i.e. for the simulation of flow patterns through spac-
ers in membrane channels (i.e. Cao, Wiley, & Fane, 2001; Dendukuri, Karode, & Kumar, 2005). The present study presents a numerical simulation using CFD of momentum and mass transfer in the lumen side of a membrane contactor for laminar flow conditions. The aim of the simulation is to calculate the concentration profiles of components A–C. The influence of various process parameters (diffusion coefficients, concentrations, velocities and inner radius) on the mass transfer is investigated. The numerical results are discussed in terms of mixing and reaction between A and B. 2. Theory 2.1. Model of velocity and concentrations Fig. 1 shows a hollow fiber of length L and radius R (2D geometry) in a cylindrical coordinate system. The real 3-D geometry can be obtained from the 2-D geometry under the hypothesis of axial symmetry. The membrane may be a tubular or a hollow fiber membrane device, according to the inner radius data. Component A is introduced at the inlet of the lumen side, with concentration CA,in . Component B flows from the inlet of the shell side, with concentration CB,in . The solvents of both components A and B are assumed identical or at least to have similar dynamic viscosity and density. The inlet velocity at the lumen side, vin,L , represents the mean velocity. The inlet velocity at the shell side, vin,S , is assumed to be uniform along the membrane length and equal to the flow rate at the shell side divided by the membrane surface (2πRL). This assumption implies that the porosity of the membrane surface is rather high, which is verified for membranes used for contactor applications. An other assumption is that there is no pore blockage nor gradual fiber blockage by any chemical species. In laminar flow, the Navier-Stokes equations apply: ∂v − ηv + ρ(v · ∇)v + ∇p = 0 ∂t ∇ ·v=0 ρ
(1)
where v, p, ρ and η denote, respectively, the velocity vector, the pressure, the density of the fluid and the dynamic viscosity. In the following, the norm of the velocity vector, v, will be noted v. The transfer equations are the convection-diffusion equations: ∂Ci + ∇ · (−Di ∇Ci + Ci v) = 0, ∂t
i = A–C
(2)
Fig. 1. Schematic representation of the lumen side of the membrane contactor (axial symmetry).
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
with Ci and Di the concentration and diffusion coefficients of component i. The boundary conditions are expressed as follows:
2.2. Numerical values of the model parameters Laminar flow is obtained on the all membrane length, if the Reynolds number calculated at the outlet of the lumen side is such that Reout < 2100, transition flow if the Reynolds number is between 2100 and 4000, and turbulent flow for higher data. The Reynolds number is given by: Reout =
(3) At the inlet of the lumen side, the flow is assumed to be fully developed and a parabolic profile for the axial velocity, vz , is specified, with vin,L being the mean velocity. At the inlet of the shell side, the radial velocity, vr , is set to −vin,S . At the outlet of the lumen side, the relative pressure is set to zero (the solution flows at the outlet of the lumen side at atmospheric pressure). The molecular fluxes of components A and B, (Ni = Cvi , i = A and B), are specified at the inlet of the lumen and shell sides. At r = 0, axial symmetry conditions are specified. For a reaction between A and B which results in the production of C, a reaction term RC must be added on the right side of Eq. (2), negative for A and B and positive for C. In this work, the simulations are realized for a second-order reaction. Although, this reaction is very simple, the simulations obtained give the general trends of the concentration profiles and the influence of parameters. For a second-order reaction, the reaction term RC is expressed as: RC = kCA CB
(4)
where k (m3 /mol s) is the rate constant. The performance of the membrane contactor is evaluated in terms of conversion rate Xi of component i, defined as (Schweich, 2001): Xi = 1 −
Ni,out , Ni,in
i = A and B
(5)
where Ni,in and Ni,out are the total molecular flows at the inlets and outlet, respectively. The conversion rate Xi represents the amount of component i which has been consumed during the reaction. Ni,out is calculated by integrating the local molecular flow at the outlet of the lumen side (z = L): Ni,out = ci (r)v(r) dS (6) z=L
where Ci and v are, respectively, the concentration of component i and the norm of the velocity vector at the outlet of the lumen side. The inlet molecular flows Ni,in of components i are calculated from the inlet concentrations and velocities data (NA,in = CA,in vin,L πR2 , NB,in = CB,in vin,S 2πRL). The software COMSOL MultiphysicsTM (Grenoble, France) is used to solve the differential equations set, by the finite element method and the direct linear system solvers UMFPACK and SPOOLES.
1327
ρvout 2R η
(7)
where vout is the mean fluid velocity at the outlet of the lumen side. The density and the viscosity data are those of water (ρ = 1 × 103 kg/m3 and η = 1 × 10−3 Pa s). The mean velocity at the outlet is obtained from a mass balance: vout =
2L vin,S + vin,L R
(8)
From Eqs. (7) and (8), we therefore obtain: Reout =
2ρ (2Lvin,S + Rvin,L ) η
(9)
Various inner radius are used for the numerical simulations: (1) 3 × 10−3 m, the inner radius of tubular membranes sold by Orelis (Miribel, France) to perform nanofiltration, ultrafiltration and microfiltration operations, (2) 3 × 10−4 m, the inner radius of hollow fiber modules, commercialized by Spectrumlab (Rancho Dominguez, USA) for ultrafiltration and microfiltration and (3) 1 × 10−4 m, the inner radius of hollow fiber modules sold by Liqui-Cel (Charlotte, USA) for membrane contactors mainly used to remove dissolved gasses. The membrane length used for the numerical simulations is L = 0.4 m. The numerical simulations are realized for a large range of inlet velocities at the lumen and shell sides, which can be obtained with these membrane devices. The permeating flowrates Qin,S were measured experimentally for the various devices, and the velocities at the inlet of the shell side were deduced from these data (vin,S = Qin,S /2πRL). A range of inlet velocities at the shell side from 1.3 to 39 × 10−4 m/s was chosen, and the value of 3.75 × 10−4 m/s was used for many simulations. For an inlet velocity at the shell side equal to 3.75 × 10−4 m/s, the flow is laminar inside the membrane device for an inlet velocity at the lumen side vin,L below 0.25 m/s for the 3 × 10−3 m inner diameter, below 2.5 m/s for the 3 × 10−4 m inner diameter, and below 7.5 m/s for the 1 × 10−4 m inner diameter. Two typical values of diffusion coefficients are used in the simulation: that of the sulphate ion, 2.13 × 10−9 m2 /s, and hydronium ion H3 O+ , 9.3 × 10−9 m2 /s, at 25 ◦ C (Lide, 2003). 3. Results 3.1. Velocity profiles The velocity profiles in a membrane contactor are simulated. A tubular membrane geometry (inner radius 3 × 10−3 m) is used, rather than a hollow fiber configuration (3 × 10−4 or 1 × 10−4 m), because the influence of process parameters
1328
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
This equation predicts that the axial mean velocity increases linearly with membrane length. The maximum velocity at the outlet of the lumen side (0.38 m/s) is slightly different from 4 × vin,L (0.4 m/s) which would be obtained if the axial velocity profile was perfectly parabolic. However, the mass balance is well verified at the outlet of the lumen side of the contactor (Qout = Qin,L + Qin,S = 5.66 × 10−6 m3 /s). The radial velocity, vr , is shown in Fig. 2c. It does not change with membrane length, and is equal to zero at the membrane centerline. Moreover, the slight increase of the radial velocity observed near the membrane surface can be found by a mass balance in the lumen side.
Fig. 2. Velocity profiles in the membrane contactor at different membrane sections (1) z = 0, (2) z = L/4, (3) z = L/2, (4) z = 3L/4, (5) z = L with L = 0.4 m; (a) two-dimensional, (b) axial velocity and (c) radial velocity (vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s).
appears then more clearly. It must be underlined that in all figures the scales of the r- and z-axis are different, to see more clearly the various profiles. Fig. 2 shows the velocity profiles in a membrane contactor for the inlet velocities vin,L = 0.1 m/s and vin,S = 3.75 × 10−4 m/s. The flow rates are identical Qin,L = Qin,S = 2.83 × 10−6 m3 /s. The norm of the velocity vector, v, is shown in Fig. 2a, the axial velocity, vz , in Fig. 2b, and the radial velocity, vr , in Fig. 2c. The axial velocity profile is almost parabolic with a mean velocity increasing with membrane length because of continuous fluid permeation. The axial velocity profile is close to the analytical solution obtained from a mass balance in the z-direction, that is to say: r 2 2vin,S vz = 2 vin,L + (10) z 1− R R
Fig. 3. Concentration profiles of A and B in the membrane contactor: (a) twodimensional (DA = DB = 9.3 × 10−9 m2 /s), (b) DA = DB = 2.13 × 10−9 m2 /s (at the outlet section z = 0.4 m) and (c) DA = DB = 9.3 × 10−9 m2 /s (at the outlet section z = 0.4 m). Parameters values vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s, CA,in = CB,in = 100 mol/m3 .
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
3.2. Concentration profiles The steady state for the concentrations is obtained at time t = L/vin,L . All the following simulations are presented in steady state. The concentrations profiles of components A and B are shown in Fig. 3, two-dimensional (Fig. 3a), and at the outlet section of the lumen side (z = 0.4 m), for diffusion coefficients equal to DA = DB = 2.13 × 10−9 m2 /s and DA = DB = 9.3 × 10−9 m2 /s (Fig. 3b and c, respectively). The velocities are the following: vin,L = 0.1 m/s and vin,S = 3.75 × 10−4 m/s (the flow rates are identical Qin,L = Qin,S = 2.83 × 10−6 m3 /s). The initial concentrations of the two solutions are assumed to be constant and equal to CA,in = CB,in = 100 mol/m3 . Fig. 3a shows that components A and B separate along the streamlines coming from the inlet of the lumen side (z = 0). Mixing is obtained by diffusion of both components along the streamlines. Fig. 3b and c show that concentration A is equal to zero near the shell inlet, and increases to reach CA,in at the membrane centreline (axis r = 0). Concentration B is equal to CB,in near the shell inlet, and then decreases to zero. When concentrations of A and B are simultaneously different from zero, the two components mix and may react together. In the following, we will call the zone where concentrations of A and B are simultaneously different from zero, the “mixing zone”. The width of the
1329
mixing/reaction zone plays a major role on the performance of the process. The influence of diffusion coefficients on the concentration profiles of A and B are shown in Fig. 3b and c, for diffusion coefficients equal to DA = DB =2.13 × 10−9 m2 /s and DA =DB =9.3 × 10−9 m2 /s, respectively. The mixing zone width is close to 400 m for DA =DB =2.13 × 10−9 m2 /s, and to 850 m for DA = DB = 9.3 × 10−9 m2 /s. In a membrane contactor process, mixing relies on molecular diffusion; therefore, the performance is strongly influenced by the diffusion coefficients of the molecular species. The mixing zone width increases with an increase in diffusion coefficients of both components. If the diffusion coefficients of A, DA , and B, DB , are not equal, which could be the case in reality, the mixing/reaction zone width will be determined by the diffusion coefficient of both components. From the general point of view, the mixing zone width will be higher than the width obtained for D = min(DA , DB ) and smaller than the width obtained for D = max(DA , DB ). 3.3. Influence of process parameters The following numerical simulations are realized to determine the effect of process parameters on the mixing zone width and its position inside the membrane. The effect of the
Fig. 4. Influence of process parameters on the concentration profiles of B at the outlet of the lumen side and for DB = 9.3 × 10−9 m/s. (a) Influence of concentration CB,in (vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s); (b) influence of the membrane length (CA,in = CB,in = 100 mol/m3 , vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s); (c) influence of the inlet velocity at the lumen side vin,L (vin,S = 3.75 × 10−4 m/s, CB,in = 100 mol/m3 ); (d) influence of the inlet velocity at the shell side vin,S (vin,L = 0.1 m/s, CB,in = 100 mol/m3 ).
1330
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
inlet concentrations, the length of the membrane and the inlet velocities are investigated. The influence of the density and the dynamic viscosity of the solvents could be evaluated in a further study. The geometry used is again a tubular geometry with an inner radius equal to 3 × 10−3 m. The influence of process parameters is shown in Fig. 4 at the outlet section of the lumen side (z = 0.4 m) and for diffusion coefficients DA = DB = 9.3 × 10−9 m2 /s. For lower diffusion coefficients (i.e. DA = DB = 2.13 × 10−9 m2 /s), the effect of the parameters would be the same with a smaller mixing zone width. Fig. 4a shows the effect of the inlet concentration, CB,in , on the concentration profile of B in the lumen side of the membrane contactor. The inlet concentration varies from 5 to 100 mol/m3 . The numerical simulations show that the width and position of the mixing zone do not changed with the inlet concentration of B. Only the concentration at the membrane surface is modified. When the dimensionless concentration CB /CB,in is plotted, the curves are found to be exactly the same. The fact that the inlet concentration of B does not modify the dimensionless profile is related to the Fick’s law diffusion equation, which resolution gives a similar result: the dimensionless concentration profile does not depend on the inlet concentration. The influence of the inlet concentration, CA,in , on the concentration profile of A is checked to be similar to the influence of CB,in on the concentration profile of B. Fig. 4b shows the variation of the concentration of component B in the lumen side of the membrane contactor at different sections along the membrane length (L/4, L/2, 3L/4, L). It appears that the mixing/reaction zone moves to the membrane centreline along the lumen side of the membrane contactor. Moreover, the width of the mixing zone increases slightly with length due to diffusion during longer times. Fig. 4c shows the effect of the inlet velocity at the lumen side, vin,L , on the concentration profile of B in the lumen side the membrane contactor. The velocity vin,L varies from 0.01 to 0.15 m/s (×15). The numerical simulations predict that an increase in vin,L moves the mixing zone near the shell side, and decreases slightly its width. Fig. 4d shows the effect of the permeating velocity, vin,S , on the concentration profile of B in the lumen side. The velocity vin,L varies by a factor of 30, from 1.3 × 10−4 to 3.9 × 10−3 m/s. An increase in vin,S leads to move the mixing zone to the membrane centreline and to a smaller mixing width. The effect of the membrane pore diameter is expected to be similar to that of the permeating velocity, as both parameters are related (a decrease in permeating velocity is observed for a decrease in membrane pore size for similar transmembrane pressure). 3.4. Reaction between A and B Fig. 5 shows the results of the simulations which involve a second-order reaction between A and B for the production of C. The inlet concentrations of A and B, CA,in and CB,in , are the same (CA,in = CB,in = 100 mol/m3 ), as well as the diffusion coefficients (DA = DB = DC = 9.3 × 10−9 m2 /s). The velocities used for the numerical simulations are: vin,L = 0.1 m/s and vin,s = 3.75 × 10−4 m/s (the flow rates are identical Qin,L = Qin,S = 2.83 × 10−6 m3 /s).
Fig. 5. Concentration profiles of A–C for a second-order reaction at z = 0.4 m. (a) Concentrations of A–C with and without reaction (k = 1 m3 /mol s) and (b) influence of rate constant. Parameters values vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s, CA,in = CB,in = 100 mol/m3 , DA = DB = DC = 9.3 × 10−9 m2 /s.
The concentrations profiles of A–C are shown in Fig. 5a, in the absence of reaction, and for a second-order reaction with a reaction kinetic constant k = 1 m3 /mol s. The reaction between A and B occurs in the mixing zone, which is a key characteristic of this membrane process. When the reaction occurs, the concentrations of A and B are lower than those obtained in the absence of reaction, due to their consumption. Fig. 5b shows the influence of the reaction rate constant, k, on the concentration profiles of C calculated at the outlet of the lumen side (z = 0.4 m). The rate constant varies from 0.001 to 1 m3 /mol. The numerical simulations show that an increasing rate constant increases the C concentration in the mixing zone, but does not change the mixing zone width. It is also interesting to compare the respective kinetics of the two consecutive mechanisms: (i) the mixing which brings the reactants in contact and (ii) the reaction itself. Thus, Fig. 5a indicates that all components A and B are consumed in the mixing zone, for k = 1 m3 /mol s. There is no zone where A and B are found simultaneously, which means that the reaction is very fast and that mixing is the limiting phenomenon. For k = 0.1 m3 /mol s (Fig. 5b), this is still the case because the concentration profile of C is close to the one obtained with k = 1 m3 /mol s. On the contrary, for k lower than 0.1 m3 /mol s, both reaction and mixing influence the conversion. In an application of the membrane contactor like precipitation between two solutions flowing from the inlets of the lumen and
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
1331
Fig. 6. Influence of inner radius: (1) 3 × 10−3 m, (2) 3 × 10−4 m, (3) 1 × 10−4 m: (a) two-dimensional concentration profiles of C and (b) concentration profiles of A–C at the outlet section of the lumen side (z = 0.4 m). Parameters values: vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 , 3.75 × 10−5 and 1.25 × 10−5 m/s, respectively, for 3 × 10−3 , 3 × 10−4 and 1 × 10−4 m inner radius, CA,in = CB,in = 100 mol/m3 , DA = DB = DC = 2.13 × 10−9 m2 /s, k = 1 m3 /mol s.
shell sides, the mixing/reaction strongly depends on the reaction constants, on the diffusion coefficients of the species and on the process parameters. For the optimization of the process, the values of the process parameters, like the inlet velocities, must be chosen according to the reaction characteristics, like rate constants. 3.5. Influence of the inner radius of the membrane The influence of the inner radius on the concentration profiles of A–C is investigated for the following conditions The initial concentrations, CA,in and CB,in , are the same (CA,in = CB,in = 100 mol/m3 ), as well as the diffusion coefficients (DA = DB = DC = 2.13 × 10−9 m2 /s). The mean velocity at the inlet of the lumen side is the same for the various inner radius (vin,L = 0.1 m/s). The flow rates at the inlet of the lumen and
shell sides are assumed equal (Qin,L = Qin,S = 2.83 × 10−6 m3 /s), therefore the velocity at the shell side, vin,L , is calculated from the mean inlet velocity at the lumen side vin,L as: vin,S = vin,L × R/2L. The membrane length is equal to L = 0.4 m. The simulations are presented in Fig. 6, for inner radius: 3 × 10−3 , 3 × 10−4 and 1 × 10−4 m. The two-dimensional concentration profiles of C are shown in Fig. 6a1–a3, and the concentration profiles of A–C at the outlet of the lumen side (z = 0.4 m) in Fig. 6b1–b3, respectively. Again, the scales of the r- and z-axis are different to show clearly the profiles inside the lumen side. Mixing is caused by diffusion along the streamlines separating components A and B. For a hollow fiber membrane device, the mixing zone width is close to the inner membrane diameter data; therefore mixing and reaction occur inside all the lumen side of the membrane contactor (Fig. 6a2 and a3). In case of a
1332
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
tubular membrane, mixing and reaction occur in a very small part of the lumen side (Fig. 6a1). Fig. 6b1–b3 show the concentration profiles of A–C at the outlet section of the lumen side. The conversion rates of A and B are calculated for the various radius. For the 3 × 10−3 m tubular membrane, the conversion rates are XA = 0.055 and XB = 0.058 (this small difference is due to the precision of the simulations). For the 3 × 10−4 m hollow fiber, XA = 0.50 and XB = 0.50, and for the 1 × 10−4 m hollow fiber XA = 0.87 and XB = 0.87. The conversions of both components A and B are identical because the inlet flow rates at the lumen and shell sides are the same. The conversion rates increase when the inner radius decreases. The very high conversion rates obtained with the hollow fiber membrane devices confirm the potential of liquid–liquid hollow fiber membrane contactors to perform mixing and reaction. These results are in accordance with experimental data recently obtained with membrane contactors. Zarkadas and Sirkar (2006) realized experiments with porous hollow fiber devices to evaluate mixing between isopropanaol (introduced at the membrane shell) and water (lumen inlet). It was shown that radial mixing in porous hollow fiber devices was significant, because the time required for the permeate to be convectively transported to the hollow fiber centreline was one to two times smaller than the mean residence time inside the device. These authors underlined that this membrane process will lead to a new method of mixing characterized by low energy input and the absence of mechanical components inside the mixing device (i.e. stirrers and motors). Jia, Zhao, Liu, He, and Liu (2006) also recently studied the experimental micromixing efficiency of a hollow fiber membrane contactor by employing competing parallel reactions. They conclude that the membrane contactor exhibits excellent micromixing effects compared with other reactors. 4. Conclusion This paper presents a numerical simulation of momentum and mass transfer in a membrane contactor for liquid–liquid phase processes and laminar flow conditions. The numerical simulations show that mixing is obtained by diffusion of A and B along the streamlines separating both components. The mixing zone width is within the region of a few hundred of microns, depends on the diffusion coefficients of A and B, and on the inlet velocities. A hollow fiber membrane device is particularly interesting for mixing/reaction because the inner diameter is close to the mixing zone width. In this work, the concentration profiles in the lumen side of the membrane contactor and the influence of the various process parameters are simulated for a simple reaction. However, the simulations give the general trends which would be obtained for more sophisticated reactions. In particular, the presented results are in accordance with our current work on the simulation of precipitation in the lumen side of a membrane contactor. The membrane contactor provides an interesting alternative to other micromixer geometries and is an example of membrane technologies for process intensification (Charpentier, 2002; Stankiewicz, 2003). Like in micromixing devices where
no turbulence occurs, the mixing in the membrane contactor relies on molecular diffusion. Accordingly, diffusive mixing can be optimized by maximisation of the constituting factors. For micromixers, it is pointed out that “the art of micromixing” translates to an efficient maximisation of interfacial surface area and concentration gradient (Hessel et al., 2005). This is also the case for hollow fiber membrane contactors. Another advantage of membrane contactors is the scale-up ability of membrane devices (Charcosset & Fessi, 2005; Zarkadas & Sirkar, 2006). The molecular flow rates may easily be increased by adding more fibers in the device, increasing the fiber length, and using devices in parallel, to achieve the desired conversion rate and productivity.
References Belfort, G. (1989). Fluid mechanics in membrane filtration: Recent developments. Journal of Membrane Science, 40, 123–147. Bellhouse, B. J., Costigan, G., Abhinava, K., & Merry, A. (2001). The performance of helical screw-thread inserts in tubular membranes. Separation and Purification Technology, 22–23, 89–113. B´enet, N., Muhr, H., Plasari, E., & Rousseaux, J. M. (2002). New technologies for the precipitation of solid particles with controlled properties. Powder Technology, 128, 93–98. Cao, Z., Wiley, D. E., & Fane, A. G. (2001). CFD simulations of net-type turbulence promoters in a narrow channel. Journal of Membrane Science, 185, 157–176. Charcosset, C., & Fessi, H. (2005). Preparation of nanoparticles with a membrane contactor. Journal of Membrane Science, 266, 115–120. Charpentier, J. C. (2002). The triplet “molecular processes-product-process” engineering: The future of chemical engineering? Chemical Engineering Science, 57, 4667–4690. Chen, G. G., Luo, G. S., Xu, J. H., & Wang, J. D. (2004). Membrane dispersion precipitation method to prepare nanoparticles. Powder Technology, 139, 180–185. Dendukuri, D., Karode, S. K., & Kumar, A. (2005). Flow visualization through spacer filled channels by computational fluid dynamics-II: Improved feed spacer designs. Journal of Membrane Science, 249, 41–49. Drioli, E., Criscuoli, A., & Curcio, E. (2003). Membrane contactors and catalytic membrane reactors in process intensification. Chemical Engineering and Technology, 26, 975–981. Fei, H., Jia, Z., Yuelian, P., & Liu, Z. (2003). Synthesis of polyaluminium chloride with a membrane reactor: Effects of operating modes. Journal of Membrane Science, 227, 15–21. Geraldes, V., Semi˜ao, V., & Pinho, M. N. (2000). Numerical modelling of mass transfer in slits with semi-permeable membrane walls. Engineering Computations, 17(3), 192–217. Ghidossi, R., Veyret, D., & Moulin, P. (2006). Computational fluid dynamics applied to membranes: State of the art and opportunities. Chemical Engineering and Processing, 45, 437–454. Hessel, V., L¨owe, H., & Sch¨onfeld, F. (2005). Micromixers—A review on passive and active mixing principles. Chemical Engineering Science, 60, 2479–2501. Jia, Z., & Liu, Z. (2002). Synthesis of nanosized BaSO4 particles with a membrane reactor: Effects of operating parameters. Journal of Membrane Science, 209, 153–161. Jia, Z., Liu, Z., & He, F. (2003). Synthesis of nanosized BaSO4 and CaCO3 particles with a membrane reactor: Effects of additives on particles. Journal of Colloid and Interface Science, 66, 322–327. Jia, Z., Zhao, Y., Liu, L., He, F., & Liu, Z. (2006). A membrane reactor intensifying micromixing: Effects of parameters on segregation index. Journal of Membrane Science, 276, 295–300. Lide, D. R. (Ed.). (2003). Handbook of chemistry and physics. Cleveland, Ohio, USA: CRC Press.
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333 Pellerin, E., Michelitsch, E., Darcovich, K., Lin, S., & Tam, C. M. (1995). Turbulent transport in membrane modules by CFD simulation in two dimensions. Journal of Membrane Science, 100, 139–153. Pharoah, J. G., Djilali, N., & Vickers, G. W. (2000). Fluid mechanics and mass transport in centrifugal membrane separation. Journal Membrane Science, 165, 277–289. Schweich, D. (Ed.). (2001). G´enie de la r´eaction chimique. New York: Technique et Documentation. Serra, C. A., & Wiesner, M. R. (2000). A comparison of rotating and stationary membrane disk filters using computational fluid dynamics. Journal of Membrane Science, 165, 19–29.
1333
Sirkar, K. K., Shanbhag, P. V., & Kovvali, A. S. (1999). Membrane in a reactor: A functional perspective. Industrial and Engineering Chemistry Research, 38, 3715–3737. Stankiewicz, A. (2003). Reactive separations for process intensification: An industrial perspective. Chemical Engineering and Processing, 42, 137–144. Wiley, D. E., & Fletcher, D. F. (2003). Techniques for computational fluid dynamics modelling of flow in membrane channels. Journal of Membrane Science, 211, 127–137. Zarkadas, D. M., & Sirkar, K. K. (2006). Antisolvent crystallization in porous hollow fiber devices. Chemical Engineering Science, 61, 5030–5048.
Computers and Chemical Engineering 32 (2008) 1325–1333
Numerical simulation of mass transfer in a liquid–liquid membrane contactor for laminar flow conditions Roland Kieffer, Catherine Charcosset ∗ , Franc¸ois Puel, Denis Mangin Laboratoire d’Automatique et de G´enie des Proc´ed´es, UMR CNRS 5007, UCBLyon 1, ESCPE-Lyon, 43 Bd du 11 Novembre 1918, 69 622 Villeurbanne Cedex, France Received 6 November 2006; received in revised form 31 May 2007; accepted 6 June 2007 Available online 20 June 2007
Abstract Liquid–liquid phase membrane contactors are increasingly being used for mixing and reaction. The principle is the following: component A flows through the membrane device inlet to mix/react with component B which comes from the membrane pores. This study presents a numerical simulation using computational fluid dynamics (CFD) of momentum and mass transfer in a tubular membrane contactor for laminar flow conditions. The velocity and concentration profiles of components A–C are obtained by resolution of the Navier-Stokes and convection-diffusion equations. The numerical simulations show that mixing between A and B is obtained by diffusion along the streamlines separating both components. The mixing/reaction zone width is within the region of a few hundred of microns, and depends on the diffusion coefficients of A and B. Hollow fiber membrane devices are found to be of particular interest because their inner diameter is close to the mixing zone width. © 2007 Elsevier Ltd. All rights reserved. Keywords: Computational fluid dynamics; Membrane contactor; Membrane reactor; Micromixing; Microreactor; Mixing
1. Introduction Membrane contactors represent an emerging technology in which the membrane is used as a tool for inter-phase mass transfer operations. These membrane systems provide a high interfacial area between two phases to achieve high overall rates of mass transfer (Drioli, Criscuoli, & Curcio, 2003; Sirkar, Shanbhag, & Kovvali, 1999). Membrane contactors have attracted considerable attention recently for liquid–gas phase processes such as simultaneous extraction and injection of gases in natural water production, and membrane distillation for desalination and wastewater treatment (Drioli et al., 2003). Liquid–liquid membrane contactors are reported for the preparation of nanoprecipitates (Chen, Luo, Xu, & Wang, 2004; Fei, Jia, Yuelian, & Liu, 2003; Jia & Liu, 2002; Jia, Liu, & He, 2003), polyaluminium chloride particles (Fei et al., 2003), l-asparagine crystals (Zarkadas & Sirkar, 2006) and polymeric nanoparticles (Charcosset & Fessi, 2005).
∗
Corresponding author. Tel.: +33 4 72 43 18 67; fax: +33 4 72 43 16 99. E-mail address: [email protected] (C. Charcosset).
0098-1354/$ – see front matter © 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.compchemeng.2007.06.013
The principle of liquid–liquid membrane contactors is the following. Component A flows through the inlet of the lumen side of the membrane contactor and component B comes from the membrane pores (shell side). Components A and B mix and react inside the lumen side of the membrane device. This membrane process (for a hollow fiber membrane device) is an alternative to micromixing and microchannel reactors. Standard mixing devices include the stirred tank reactor, the T-tube or Y-tube, and more sophisticated ones, like the two-impinging-jets mixing device, the sliding surface mixing device and the vortex reactor (B´enet, Muhr, Plasari, & Rousseaux, 2002). The common feature of micromixing and microreactor designs is small channels with dimensions below 1 mm. Application fields of micromixers include both modern, specialised issues such as sample preparation for chemical analysis and traditional, widespread usable mixing tasks such as reaction, gas absorption, emulsification, foaming and blending (Hessel, L¨owe, & Sch¨onfeld, 2005). CFD has been largely used as a tool to model membrane separation processes (Belfort, 1989; Ghidossi, Veyret, & Moulin, 2006). The Navier-Stokes equations were solved to estimate the pressure drop across the membrane module and the shear rate on the membrane surface. Several studies used the concentration polarization model, and other studies were carried
1326
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
Nomenclature CA CB CA,in CB,in D Di k L Ni Ni,in Ni,out p Qin,L Qin,S Qout r R RC Reout v v vin,L vin,S vout Xi z
concentration of component A (mol/m3 ) concentration of component B (mol/m3 ) concentration of component A at the inlet of the lumen side (mol/m3 ) concentration of component B at the inlet of the shell side (mol/m3 ) inner diameter (m) diffusion coefficient of component i (i = A–C) (m2 /s) reaction rate constant (m3 /mol s) membrane length (m) molecular flux of component i (i = A and B) (mol/m2 s) total molecular flux of component i calculated at the inlets (mol/s) total molecular flux of component i calculated at the outlet of the lumen side (mol/s) pressure (Pa) flow rate at the inlet of the lumen side (m3 /s) flow rate at the inlet of the shell side (m3 /s) flow rate at the outlet of the lumen side (m3 /s) radial coordinate (m) inner radius (m) second-order reaction term for the preparation of C (mol/m3 s) Reynolds number calculated at the outlet of the lumen side velocity vector (m/s) norm of the velocity vector (m/s) mean velocity at the inlet of the lumen side (m/s) velocity at the inlet of the shell side (m/s) mean velocity at the outlet of the lumen side (m/s) conversion rate of component i (i = A and B) axial coordinate (m)
Greek letters η dynamic viscosity (Pa s) ρ fluid density (kg/m3 )
out inside the membrane porous matrix. Numerical simulations were performed for various membrane filtration processes: microfiltration and ultrafiltration (Ghidossi et al., 2006; Pellerin, Michelitsch, Darcovich, Lin, & Tam, 1995; Wiley & Fletcher, 2003), nanofiltration and reverse osmosis (i.e. Geraldes, Semi˜ao, & Pinho, 2000). Various membrane module configurations were also investigated: flat sheet (Geraldes et al., 2000; Pellerin et al., 1995; Wiley & Fletcher, 2003), rotating circular membrane geometry (Serra & Wiesner, 2000), tubular membranes with inserts (Bellhouse, Costigan, Abhinava, & Merry, 2001) and centrifugal membrane separation (Pharoah, Djilali, & Vickers, 2000). CFD was also used to understand the decrease in membrane fouling by using different instabilities or turbulence promoters, i.e. for the simulation of flow patterns through spac-
ers in membrane channels (i.e. Cao, Wiley, & Fane, 2001; Dendukuri, Karode, & Kumar, 2005). The present study presents a numerical simulation using CFD of momentum and mass transfer in the lumen side of a membrane contactor for laminar flow conditions. The aim of the simulation is to calculate the concentration profiles of components A–C. The influence of various process parameters (diffusion coefficients, concentrations, velocities and inner radius) on the mass transfer is investigated. The numerical results are discussed in terms of mixing and reaction between A and B. 2. Theory 2.1. Model of velocity and concentrations Fig. 1 shows a hollow fiber of length L and radius R (2D geometry) in a cylindrical coordinate system. The real 3-D geometry can be obtained from the 2-D geometry under the hypothesis of axial symmetry. The membrane may be a tubular or a hollow fiber membrane device, according to the inner radius data. Component A is introduced at the inlet of the lumen side, with concentration CA,in . Component B flows from the inlet of the shell side, with concentration CB,in . The solvents of both components A and B are assumed identical or at least to have similar dynamic viscosity and density. The inlet velocity at the lumen side, vin,L , represents the mean velocity. The inlet velocity at the shell side, vin,S , is assumed to be uniform along the membrane length and equal to the flow rate at the shell side divided by the membrane surface (2πRL). This assumption implies that the porosity of the membrane surface is rather high, which is verified for membranes used for contactor applications. An other assumption is that there is no pore blockage nor gradual fiber blockage by any chemical species. In laminar flow, the Navier-Stokes equations apply: ∂v − ηv + ρ(v · ∇)v + ∇p = 0 ∂t ∇ ·v=0 ρ
(1)
where v, p, ρ and η denote, respectively, the velocity vector, the pressure, the density of the fluid and the dynamic viscosity. In the following, the norm of the velocity vector, v, will be noted v. The transfer equations are the convection-diffusion equations: ∂Ci + ∇ · (−Di ∇Ci + Ci v) = 0, ∂t
i = A–C
(2)
Fig. 1. Schematic representation of the lumen side of the membrane contactor (axial symmetry).
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
with Ci and Di the concentration and diffusion coefficients of component i. The boundary conditions are expressed as follows:
2.2. Numerical values of the model parameters Laminar flow is obtained on the all membrane length, if the Reynolds number calculated at the outlet of the lumen side is such that Reout < 2100, transition flow if the Reynolds number is between 2100 and 4000, and turbulent flow for higher data. The Reynolds number is given by: Reout =
(3) At the inlet of the lumen side, the flow is assumed to be fully developed and a parabolic profile for the axial velocity, vz , is specified, with vin,L being the mean velocity. At the inlet of the shell side, the radial velocity, vr , is set to −vin,S . At the outlet of the lumen side, the relative pressure is set to zero (the solution flows at the outlet of the lumen side at atmospheric pressure). The molecular fluxes of components A and B, (Ni = Cvi , i = A and B), are specified at the inlet of the lumen and shell sides. At r = 0, axial symmetry conditions are specified. For a reaction between A and B which results in the production of C, a reaction term RC must be added on the right side of Eq. (2), negative for A and B and positive for C. In this work, the simulations are realized for a second-order reaction. Although, this reaction is very simple, the simulations obtained give the general trends of the concentration profiles and the influence of parameters. For a second-order reaction, the reaction term RC is expressed as: RC = kCA CB
(4)
where k (m3 /mol s) is the rate constant. The performance of the membrane contactor is evaluated in terms of conversion rate Xi of component i, defined as (Schweich, 2001): Xi = 1 −
Ni,out , Ni,in
i = A and B
(5)
where Ni,in and Ni,out are the total molecular flows at the inlets and outlet, respectively. The conversion rate Xi represents the amount of component i which has been consumed during the reaction. Ni,out is calculated by integrating the local molecular flow at the outlet of the lumen side (z = L): Ni,out = ci (r)v(r) dS (6) z=L
where Ci and v are, respectively, the concentration of component i and the norm of the velocity vector at the outlet of the lumen side. The inlet molecular flows Ni,in of components i are calculated from the inlet concentrations and velocities data (NA,in = CA,in vin,L πR2 , NB,in = CB,in vin,S 2πRL). The software COMSOL MultiphysicsTM (Grenoble, France) is used to solve the differential equations set, by the finite element method and the direct linear system solvers UMFPACK and SPOOLES.
1327
ρvout 2R η
(7)
where vout is the mean fluid velocity at the outlet of the lumen side. The density and the viscosity data are those of water (ρ = 1 × 103 kg/m3 and η = 1 × 10−3 Pa s). The mean velocity at the outlet is obtained from a mass balance: vout =
2L vin,S + vin,L R
(8)
From Eqs. (7) and (8), we therefore obtain: Reout =
2ρ (2Lvin,S + Rvin,L ) η
(9)
Various inner radius are used for the numerical simulations: (1) 3 × 10−3 m, the inner radius of tubular membranes sold by Orelis (Miribel, France) to perform nanofiltration, ultrafiltration and microfiltration operations, (2) 3 × 10−4 m, the inner radius of hollow fiber modules, commercialized by Spectrumlab (Rancho Dominguez, USA) for ultrafiltration and microfiltration and (3) 1 × 10−4 m, the inner radius of hollow fiber modules sold by Liqui-Cel (Charlotte, USA) for membrane contactors mainly used to remove dissolved gasses. The membrane length used for the numerical simulations is L = 0.4 m. The numerical simulations are realized for a large range of inlet velocities at the lumen and shell sides, which can be obtained with these membrane devices. The permeating flowrates Qin,S were measured experimentally for the various devices, and the velocities at the inlet of the shell side were deduced from these data (vin,S = Qin,S /2πRL). A range of inlet velocities at the shell side from 1.3 to 39 × 10−4 m/s was chosen, and the value of 3.75 × 10−4 m/s was used for many simulations. For an inlet velocity at the shell side equal to 3.75 × 10−4 m/s, the flow is laminar inside the membrane device for an inlet velocity at the lumen side vin,L below 0.25 m/s for the 3 × 10−3 m inner diameter, below 2.5 m/s for the 3 × 10−4 m inner diameter, and below 7.5 m/s for the 1 × 10−4 m inner diameter. Two typical values of diffusion coefficients are used in the simulation: that of the sulphate ion, 2.13 × 10−9 m2 /s, and hydronium ion H3 O+ , 9.3 × 10−9 m2 /s, at 25 ◦ C (Lide, 2003). 3. Results 3.1. Velocity profiles The velocity profiles in a membrane contactor are simulated. A tubular membrane geometry (inner radius 3 × 10−3 m) is used, rather than a hollow fiber configuration (3 × 10−4 or 1 × 10−4 m), because the influence of process parameters
1328
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
This equation predicts that the axial mean velocity increases linearly with membrane length. The maximum velocity at the outlet of the lumen side (0.38 m/s) is slightly different from 4 × vin,L (0.4 m/s) which would be obtained if the axial velocity profile was perfectly parabolic. However, the mass balance is well verified at the outlet of the lumen side of the contactor (Qout = Qin,L + Qin,S = 5.66 × 10−6 m3 /s). The radial velocity, vr , is shown in Fig. 2c. It does not change with membrane length, and is equal to zero at the membrane centerline. Moreover, the slight increase of the radial velocity observed near the membrane surface can be found by a mass balance in the lumen side.
Fig. 2. Velocity profiles in the membrane contactor at different membrane sections (1) z = 0, (2) z = L/4, (3) z = L/2, (4) z = 3L/4, (5) z = L with L = 0.4 m; (a) two-dimensional, (b) axial velocity and (c) radial velocity (vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s).
appears then more clearly. It must be underlined that in all figures the scales of the r- and z-axis are different, to see more clearly the various profiles. Fig. 2 shows the velocity profiles in a membrane contactor for the inlet velocities vin,L = 0.1 m/s and vin,S = 3.75 × 10−4 m/s. The flow rates are identical Qin,L = Qin,S = 2.83 × 10−6 m3 /s. The norm of the velocity vector, v, is shown in Fig. 2a, the axial velocity, vz , in Fig. 2b, and the radial velocity, vr , in Fig. 2c. The axial velocity profile is almost parabolic with a mean velocity increasing with membrane length because of continuous fluid permeation. The axial velocity profile is close to the analytical solution obtained from a mass balance in the z-direction, that is to say: r 2 2vin,S vz = 2 vin,L + (10) z 1− R R
Fig. 3. Concentration profiles of A and B in the membrane contactor: (a) twodimensional (DA = DB = 9.3 × 10−9 m2 /s), (b) DA = DB = 2.13 × 10−9 m2 /s (at the outlet section z = 0.4 m) and (c) DA = DB = 9.3 × 10−9 m2 /s (at the outlet section z = 0.4 m). Parameters values vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s, CA,in = CB,in = 100 mol/m3 .
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
3.2. Concentration profiles The steady state for the concentrations is obtained at time t = L/vin,L . All the following simulations are presented in steady state. The concentrations profiles of components A and B are shown in Fig. 3, two-dimensional (Fig. 3a), and at the outlet section of the lumen side (z = 0.4 m), for diffusion coefficients equal to DA = DB = 2.13 × 10−9 m2 /s and DA = DB = 9.3 × 10−9 m2 /s (Fig. 3b and c, respectively). The velocities are the following: vin,L = 0.1 m/s and vin,S = 3.75 × 10−4 m/s (the flow rates are identical Qin,L = Qin,S = 2.83 × 10−6 m3 /s). The initial concentrations of the two solutions are assumed to be constant and equal to CA,in = CB,in = 100 mol/m3 . Fig. 3a shows that components A and B separate along the streamlines coming from the inlet of the lumen side (z = 0). Mixing is obtained by diffusion of both components along the streamlines. Fig. 3b and c show that concentration A is equal to zero near the shell inlet, and increases to reach CA,in at the membrane centreline (axis r = 0). Concentration B is equal to CB,in near the shell inlet, and then decreases to zero. When concentrations of A and B are simultaneously different from zero, the two components mix and may react together. In the following, we will call the zone where concentrations of A and B are simultaneously different from zero, the “mixing zone”. The width of the
1329
mixing/reaction zone plays a major role on the performance of the process. The influence of diffusion coefficients on the concentration profiles of A and B are shown in Fig. 3b and c, for diffusion coefficients equal to DA = DB =2.13 × 10−9 m2 /s and DA =DB =9.3 × 10−9 m2 /s, respectively. The mixing zone width is close to 400 m for DA =DB =2.13 × 10−9 m2 /s, and to 850 m for DA = DB = 9.3 × 10−9 m2 /s. In a membrane contactor process, mixing relies on molecular diffusion; therefore, the performance is strongly influenced by the diffusion coefficients of the molecular species. The mixing zone width increases with an increase in diffusion coefficients of both components. If the diffusion coefficients of A, DA , and B, DB , are not equal, which could be the case in reality, the mixing/reaction zone width will be determined by the diffusion coefficient of both components. From the general point of view, the mixing zone width will be higher than the width obtained for D = min(DA , DB ) and smaller than the width obtained for D = max(DA , DB ). 3.3. Influence of process parameters The following numerical simulations are realized to determine the effect of process parameters on the mixing zone width and its position inside the membrane. The effect of the
Fig. 4. Influence of process parameters on the concentration profiles of B at the outlet of the lumen side and for DB = 9.3 × 10−9 m/s. (a) Influence of concentration CB,in (vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s); (b) influence of the membrane length (CA,in = CB,in = 100 mol/m3 , vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s); (c) influence of the inlet velocity at the lumen side vin,L (vin,S = 3.75 × 10−4 m/s, CB,in = 100 mol/m3 ); (d) influence of the inlet velocity at the shell side vin,S (vin,L = 0.1 m/s, CB,in = 100 mol/m3 ).
1330
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
inlet concentrations, the length of the membrane and the inlet velocities are investigated. The influence of the density and the dynamic viscosity of the solvents could be evaluated in a further study. The geometry used is again a tubular geometry with an inner radius equal to 3 × 10−3 m. The influence of process parameters is shown in Fig. 4 at the outlet section of the lumen side (z = 0.4 m) and for diffusion coefficients DA = DB = 9.3 × 10−9 m2 /s. For lower diffusion coefficients (i.e. DA = DB = 2.13 × 10−9 m2 /s), the effect of the parameters would be the same with a smaller mixing zone width. Fig. 4a shows the effect of the inlet concentration, CB,in , on the concentration profile of B in the lumen side of the membrane contactor. The inlet concentration varies from 5 to 100 mol/m3 . The numerical simulations show that the width and position of the mixing zone do not changed with the inlet concentration of B. Only the concentration at the membrane surface is modified. When the dimensionless concentration CB /CB,in is plotted, the curves are found to be exactly the same. The fact that the inlet concentration of B does not modify the dimensionless profile is related to the Fick’s law diffusion equation, which resolution gives a similar result: the dimensionless concentration profile does not depend on the inlet concentration. The influence of the inlet concentration, CA,in , on the concentration profile of A is checked to be similar to the influence of CB,in on the concentration profile of B. Fig. 4b shows the variation of the concentration of component B in the lumen side of the membrane contactor at different sections along the membrane length (L/4, L/2, 3L/4, L). It appears that the mixing/reaction zone moves to the membrane centreline along the lumen side of the membrane contactor. Moreover, the width of the mixing zone increases slightly with length due to diffusion during longer times. Fig. 4c shows the effect of the inlet velocity at the lumen side, vin,L , on the concentration profile of B in the lumen side the membrane contactor. The velocity vin,L varies from 0.01 to 0.15 m/s (×15). The numerical simulations predict that an increase in vin,L moves the mixing zone near the shell side, and decreases slightly its width. Fig. 4d shows the effect of the permeating velocity, vin,S , on the concentration profile of B in the lumen side. The velocity vin,L varies by a factor of 30, from 1.3 × 10−4 to 3.9 × 10−3 m/s. An increase in vin,S leads to move the mixing zone to the membrane centreline and to a smaller mixing width. The effect of the membrane pore diameter is expected to be similar to that of the permeating velocity, as both parameters are related (a decrease in permeating velocity is observed for a decrease in membrane pore size for similar transmembrane pressure). 3.4. Reaction between A and B Fig. 5 shows the results of the simulations which involve a second-order reaction between A and B for the production of C. The inlet concentrations of A and B, CA,in and CB,in , are the same (CA,in = CB,in = 100 mol/m3 ), as well as the diffusion coefficients (DA = DB = DC = 9.3 × 10−9 m2 /s). The velocities used for the numerical simulations are: vin,L = 0.1 m/s and vin,s = 3.75 × 10−4 m/s (the flow rates are identical Qin,L = Qin,S = 2.83 × 10−6 m3 /s).
Fig. 5. Concentration profiles of A–C for a second-order reaction at z = 0.4 m. (a) Concentrations of A–C with and without reaction (k = 1 m3 /mol s) and (b) influence of rate constant. Parameters values vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 m/s, CA,in = CB,in = 100 mol/m3 , DA = DB = DC = 9.3 × 10−9 m2 /s.
The concentrations profiles of A–C are shown in Fig. 5a, in the absence of reaction, and for a second-order reaction with a reaction kinetic constant k = 1 m3 /mol s. The reaction between A and B occurs in the mixing zone, which is a key characteristic of this membrane process. When the reaction occurs, the concentrations of A and B are lower than those obtained in the absence of reaction, due to their consumption. Fig. 5b shows the influence of the reaction rate constant, k, on the concentration profiles of C calculated at the outlet of the lumen side (z = 0.4 m). The rate constant varies from 0.001 to 1 m3 /mol. The numerical simulations show that an increasing rate constant increases the C concentration in the mixing zone, but does not change the mixing zone width. It is also interesting to compare the respective kinetics of the two consecutive mechanisms: (i) the mixing which brings the reactants in contact and (ii) the reaction itself. Thus, Fig. 5a indicates that all components A and B are consumed in the mixing zone, for k = 1 m3 /mol s. There is no zone where A and B are found simultaneously, which means that the reaction is very fast and that mixing is the limiting phenomenon. For k = 0.1 m3 /mol s (Fig. 5b), this is still the case because the concentration profile of C is close to the one obtained with k = 1 m3 /mol s. On the contrary, for k lower than 0.1 m3 /mol s, both reaction and mixing influence the conversion. In an application of the membrane contactor like precipitation between two solutions flowing from the inlets of the lumen and
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
1331
Fig. 6. Influence of inner radius: (1) 3 × 10−3 m, (2) 3 × 10−4 m, (3) 1 × 10−4 m: (a) two-dimensional concentration profiles of C and (b) concentration profiles of A–C at the outlet section of the lumen side (z = 0.4 m). Parameters values: vin,L = 0.1 m/s, vin,S = 3.75 × 10−4 , 3.75 × 10−5 and 1.25 × 10−5 m/s, respectively, for 3 × 10−3 , 3 × 10−4 and 1 × 10−4 m inner radius, CA,in = CB,in = 100 mol/m3 , DA = DB = DC = 2.13 × 10−9 m2 /s, k = 1 m3 /mol s.
shell sides, the mixing/reaction strongly depends on the reaction constants, on the diffusion coefficients of the species and on the process parameters. For the optimization of the process, the values of the process parameters, like the inlet velocities, must be chosen according to the reaction characteristics, like rate constants. 3.5. Influence of the inner radius of the membrane The influence of the inner radius on the concentration profiles of A–C is investigated for the following conditions The initial concentrations, CA,in and CB,in , are the same (CA,in = CB,in = 100 mol/m3 ), as well as the diffusion coefficients (DA = DB = DC = 2.13 × 10−9 m2 /s). The mean velocity at the inlet of the lumen side is the same for the various inner radius (vin,L = 0.1 m/s). The flow rates at the inlet of the lumen and
shell sides are assumed equal (Qin,L = Qin,S = 2.83 × 10−6 m3 /s), therefore the velocity at the shell side, vin,L , is calculated from the mean inlet velocity at the lumen side vin,L as: vin,S = vin,L × R/2L. The membrane length is equal to L = 0.4 m. The simulations are presented in Fig. 6, for inner radius: 3 × 10−3 , 3 × 10−4 and 1 × 10−4 m. The two-dimensional concentration profiles of C are shown in Fig. 6a1–a3, and the concentration profiles of A–C at the outlet of the lumen side (z = 0.4 m) in Fig. 6b1–b3, respectively. Again, the scales of the r- and z-axis are different to show clearly the profiles inside the lumen side. Mixing is caused by diffusion along the streamlines separating components A and B. For a hollow fiber membrane device, the mixing zone width is close to the inner membrane diameter data; therefore mixing and reaction occur inside all the lumen side of the membrane contactor (Fig. 6a2 and a3). In case of a
1332
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333
tubular membrane, mixing and reaction occur in a very small part of the lumen side (Fig. 6a1). Fig. 6b1–b3 show the concentration profiles of A–C at the outlet section of the lumen side. The conversion rates of A and B are calculated for the various radius. For the 3 × 10−3 m tubular membrane, the conversion rates are XA = 0.055 and XB = 0.058 (this small difference is due to the precision of the simulations). For the 3 × 10−4 m hollow fiber, XA = 0.50 and XB = 0.50, and for the 1 × 10−4 m hollow fiber XA = 0.87 and XB = 0.87. The conversions of both components A and B are identical because the inlet flow rates at the lumen and shell sides are the same. The conversion rates increase when the inner radius decreases. The very high conversion rates obtained with the hollow fiber membrane devices confirm the potential of liquid–liquid hollow fiber membrane contactors to perform mixing and reaction. These results are in accordance with experimental data recently obtained with membrane contactors. Zarkadas and Sirkar (2006) realized experiments with porous hollow fiber devices to evaluate mixing between isopropanaol (introduced at the membrane shell) and water (lumen inlet). It was shown that radial mixing in porous hollow fiber devices was significant, because the time required for the permeate to be convectively transported to the hollow fiber centreline was one to two times smaller than the mean residence time inside the device. These authors underlined that this membrane process will lead to a new method of mixing characterized by low energy input and the absence of mechanical components inside the mixing device (i.e. stirrers and motors). Jia, Zhao, Liu, He, and Liu (2006) also recently studied the experimental micromixing efficiency of a hollow fiber membrane contactor by employing competing parallel reactions. They conclude that the membrane contactor exhibits excellent micromixing effects compared with other reactors. 4. Conclusion This paper presents a numerical simulation of momentum and mass transfer in a membrane contactor for liquid–liquid phase processes and laminar flow conditions. The numerical simulations show that mixing is obtained by diffusion of A and B along the streamlines separating both components. The mixing zone width is within the region of a few hundred of microns, depends on the diffusion coefficients of A and B, and on the inlet velocities. A hollow fiber membrane device is particularly interesting for mixing/reaction because the inner diameter is close to the mixing zone width. In this work, the concentration profiles in the lumen side of the membrane contactor and the influence of the various process parameters are simulated for a simple reaction. However, the simulations give the general trends which would be obtained for more sophisticated reactions. In particular, the presented results are in accordance with our current work on the simulation of precipitation in the lumen side of a membrane contactor. The membrane contactor provides an interesting alternative to other micromixer geometries and is an example of membrane technologies for process intensification (Charpentier, 2002; Stankiewicz, 2003). Like in micromixing devices where
no turbulence occurs, the mixing in the membrane contactor relies on molecular diffusion. Accordingly, diffusive mixing can be optimized by maximisation of the constituting factors. For micromixers, it is pointed out that “the art of micromixing” translates to an efficient maximisation of interfacial surface area and concentration gradient (Hessel et al., 2005). This is also the case for hollow fiber membrane contactors. Another advantage of membrane contactors is the scale-up ability of membrane devices (Charcosset & Fessi, 2005; Zarkadas & Sirkar, 2006). The molecular flow rates may easily be increased by adding more fibers in the device, increasing the fiber length, and using devices in parallel, to achieve the desired conversion rate and productivity.
References Belfort, G. (1989). Fluid mechanics in membrane filtration: Recent developments. Journal of Membrane Science, 40, 123–147. Bellhouse, B. J., Costigan, G., Abhinava, K., & Merry, A. (2001). The performance of helical screw-thread inserts in tubular membranes. Separation and Purification Technology, 22–23, 89–113. B´enet, N., Muhr, H., Plasari, E., & Rousseaux, J. M. (2002). New technologies for the precipitation of solid particles with controlled properties. Powder Technology, 128, 93–98. Cao, Z., Wiley, D. E., & Fane, A. G. (2001). CFD simulations of net-type turbulence promoters in a narrow channel. Journal of Membrane Science, 185, 157–176. Charcosset, C., & Fessi, H. (2005). Preparation of nanoparticles with a membrane contactor. Journal of Membrane Science, 266, 115–120. Charpentier, J. C. (2002). The triplet “molecular processes-product-process” engineering: The future of chemical engineering? Chemical Engineering Science, 57, 4667–4690. Chen, G. G., Luo, G. S., Xu, J. H., & Wang, J. D. (2004). Membrane dispersion precipitation method to prepare nanoparticles. Powder Technology, 139, 180–185. Dendukuri, D., Karode, S. K., & Kumar, A. (2005). Flow visualization through spacer filled channels by computational fluid dynamics-II: Improved feed spacer designs. Journal of Membrane Science, 249, 41–49. Drioli, E., Criscuoli, A., & Curcio, E. (2003). Membrane contactors and catalytic membrane reactors in process intensification. Chemical Engineering and Technology, 26, 975–981. Fei, H., Jia, Z., Yuelian, P., & Liu, Z. (2003). Synthesis of polyaluminium chloride with a membrane reactor: Effects of operating modes. Journal of Membrane Science, 227, 15–21. Geraldes, V., Semi˜ao, V., & Pinho, M. N. (2000). Numerical modelling of mass transfer in slits with semi-permeable membrane walls. Engineering Computations, 17(3), 192–217. Ghidossi, R., Veyret, D., & Moulin, P. (2006). Computational fluid dynamics applied to membranes: State of the art and opportunities. Chemical Engineering and Processing, 45, 437–454. Hessel, V., L¨owe, H., & Sch¨onfeld, F. (2005). Micromixers—A review on passive and active mixing principles. Chemical Engineering Science, 60, 2479–2501. Jia, Z., & Liu, Z. (2002). Synthesis of nanosized BaSO4 particles with a membrane reactor: Effects of operating parameters. Journal of Membrane Science, 209, 153–161. Jia, Z., Liu, Z., & He, F. (2003). Synthesis of nanosized BaSO4 and CaCO3 particles with a membrane reactor: Effects of additives on particles. Journal of Colloid and Interface Science, 66, 322–327. Jia, Z., Zhao, Y., Liu, L., He, F., & Liu, Z. (2006). A membrane reactor intensifying micromixing: Effects of parameters on segregation index. Journal of Membrane Science, 276, 295–300. Lide, D. R. (Ed.). (2003). Handbook of chemistry and physics. Cleveland, Ohio, USA: CRC Press.
R. Kieffer et al. / Computers and Chemical Engineering 32 (2008) 1325–1333 Pellerin, E., Michelitsch, E., Darcovich, K., Lin, S., & Tam, C. M. (1995). Turbulent transport in membrane modules by CFD simulation in two dimensions. Journal of Membrane Science, 100, 139–153. Pharoah, J. G., Djilali, N., & Vickers, G. W. (2000). Fluid mechanics and mass transport in centrifugal membrane separation. Journal Membrane Science, 165, 277–289. Schweich, D. (Ed.). (2001). G´enie de la r´eaction chimique. New York: Technique et Documentation. Serra, C. A., & Wiesner, M. R. (2000). A comparison of rotating and stationary membrane disk filters using computational fluid dynamics. Journal of Membrane Science, 165, 19–29.
1333
Sirkar, K. K., Shanbhag, P. V., & Kovvali, A. S. (1999). Membrane in a reactor: A functional perspective. Industrial and Engineering Chemistry Research, 38, 3715–3737. Stankiewicz, A. (2003). Reactive separations for process intensification: An industrial perspective. Chemical Engineering and Processing, 42, 137–144. Wiley, D. E., & Fletcher, D. F. (2003). Techniques for computational fluid dynamics modelling of flow in membrane channels. Journal of Membrane Science, 211, 127–137. Zarkadas, D. M., & Sirkar, K. K. (2006). Antisolvent crystallization in porous hollow fiber devices. Chemical Engineering Science, 61, 5030–5048.