3D multiphase flow simulation of Marangoni convection on reactive absorption of CO2 by monoethanolamine in microchannel

Shuai Chen,Jiahong Lan,Yu Zhang,Jia Guo,Zhikai Cao,Yong Sha

Department of Chemical and Biochemical Engineering,College of Chemistry and Chemical Engineering,Xiamen University,Xiamen 361005,China

Keywords:Numerical simulation Microchannel Falling film Marangoni convection MEA-CO2

ABSTRACT A multiphase flow 3D numerical simulation method employing the coupled volume of fluid (VOF) and level set model is established to study the reactive absorption of CO2 by the monoethanolamine (MEA)aqueous solution in a falling film microchannel.Based on the flow-reaction-mass transfer model of the MEA-CO2 system in the falling film microchannel,the enhancement effect of the Marangoni convection in this reactive absorption process is analyzed.The enhancement factor of the Marangoni convection obtained in this work is in good agreement with experimental results in the literature.With consideration of the absorption ratio as well as the enhancement effect of the Marangoni convection,the influence of different MEA concentrations on absorption of CO2 is investigated.Furthermore,the appropriate MEA concentration for absorption enhanced by the Marangoni convection is acquired.

1.Introduction

Absorption of CO2by the monoethanolamine (MEA) aqueous solution is the most utilized CO2removal method,and it has obvious advantages in absorption efficiency[1].Usually,it is conducted by means of the column equipment [2–4],such as the packed column which utilizes the falling film along the packing surface to carry out absorption [5].

Nowadays,the microchannel falling film can be used in various gas–liquid reaction processes,such as chemical absorption,sulfonation,chlorination,etc.[6–8].The recently developed microchannel technology can make the liquid film thickness to 10 μm and the interface area beyond 1.0×104m2.m-3.In the last decade,lots of progress have been achieved in the microchannel technology,and some industrial microchannel devices has been also utilized for industrial productions[9].The microchannel technology has a pronounced enhancement effect on the mass and heat transfer process[10–12].CO2,as an exhaust gas of many chemical production processes,many works focus on its absorption with high efficiency.The microchannel technique,such as falling film absorption,can significantly improve the CO2absorption effect[13].

In addition,many studies showed that the Marangoni convection could have a significant enhancement effect on mass transfer in traditional falling film devices [14,15].For the microchannel,experiments proved that the Marangoni convection can also affect the mass transfer process[16–18].Sobieszuk et al.[13]studied the MEA–CO2reactive absorption system in the microchannel through experiments,and results showed that the Marangoni convection did exist in the microchannel.The enhanced effect of the Marangoni convection on the CO2reactive absorption can reach as more as 5–6 times under certain operating conditions.

However,the experimental study can only characterize the reaction-mass transfer process in the microchannel from the macroscopic view of the mass transfer coefficient.It is difficult to obtain the detailed flow details,concentration,temperature distribution and other information in the microchannel[19].As a direct research method,the numerical simulation method can obtain the local mass transfer information,which is difficult to be obtained in experiments[20–22].Therefore,it is necessary to carry out further numerical simulation investigation in order to utilize the Marangoni convection to enhance the absorption of CO2by the MEA aqueous solution in the falling film microchannel.

In this work,for a falling film microchannel,a multiphase flow numerical simulation method with coupled volume of fluid (VOF)method and level set method is established to investigate the effect of the Marangoni convection on the MEA–CO2reactive absorption process.The reactive absorption process of CO2by the MEA aqueous solution is simulated and analyzed.The appropriate MEA concentration for enhanced absorption is acquired.

2.Numerical Model

2.1.Physical model

The simulation object in this work is a single microchannel,and its structure size is shown in Fig.1.The simulation utilizes the same size as the single microchannel in the experiment of Sobieszuk et al.[12].For the MEA–CO2reactive absorption process,their work showed that the Marangoni convection,which can significantly enhance absorption of CO2in the MEA aqueous solution,can occur along with the falling film in the microchannel.

The length,width and thickness of the liquid film in the microchannel are 78,0.6 and 0.3 mm respectively.The thickness of the gas film is 3.0 mm.The gas and liquid in the channel are co-current flowing.The inlet is the velocity inlet,and the outlet is the pressure outlet.The non-slip wall is used in four around walls.The width of the single microchannel is generally between 0.3 and 1.0 mm in the literature,so 0.6 mm width utilized in this work is a relatively common width [18].

2.2.Governing equation

The multiphase flow model adopts the coupling method of the VOF model [23] and level set model [24].The VOF method tracks the volume fraction of different phases by solving the VOF function transport equation,and tracks the interface by the piecewise linear interface calculation method[25].However,since the volume fraction in the VOF model is a discontinuous function,specific errors could occur in calculation of the surface tension.It could generate spurious flows at the interface and have adverse effects on accuracy of the numerical simulation.Therefore,the VOF model and level set model are coupled to establish a multiphase flow model in this work.For the numerical implementation,firstly,the unit normal vector of the phase interface is calculated by using the level set function,and a new phase interface is obtained through the VOF function coupled with the piecewise linear reconstruction method.Then,the level set function is initialized again based on the new phase interface,and the phase conservation can be guaranteed.Detailed calculation procedures can be found in the work of Chakraborty et al.[26].

The transport equation of VOF function and level set function are as follows:

Fig.1.Single microchannel.

where α is the VOF function,defined as the volume fraction of the target fluid in the cell.0<α<1 means in the interface cell,and α=0 or 1 means in the single phase cell.φ is the level set function,and it represents the positive and negative distance to the interface.φ=0 means in the phase interface,φ>0 means in the main phase,φ<0 means in the secondary phase.The volume fraction of the target fluid in the cell is obtained by Eq.(1),and the interface normal vector is obtained by Eq.(2).

The Navier-Stokes equation and the continuity equation for the entire computational domain are as follows:

where αqis the volume fraction of phase q;Sqis the mass source phase of phase q;F represents the surface tension per unit volume of the fluid;g is the acceleration of gravity.

The continuum surface force model[27]is utilized to represent the surface tensionF.It divides the surface tension into forces parallel and perpendicular to the phase interface.The force perpendicular to the phase interfaceis calculated by the following equation:

where a is the grid size.

In different regions,density ρ and viscosity μ are calculated by the Heaviside function:

It is known that the Marangoni convection is generated by the forceparallel to the interface,and it is calculated by Eq.(8):

where ∇Sis the gradient along the interface.

Since the absorption process of CO2by the MEA aqueous solution involves the reaction in the liquid phase and the interphase mass transfer process between the gas phase and liquid phase,it is necessary to obtain the reaction and interphase mass transfer source in the component transport equation.

The phase-averaged species equation,written for an arbitrary phase q,is as following.

where Ji，qis the diffusion flux of component i in the phase q;Si，qis the mass source of component i in the phase q.

For the reaction mechanism between MEA and CO2,Sada et al.[28]considered it as a second-order irreversible reaction as shown as Eq.(10).

where R is the reaction rate which is zero in the gas.k2is the second-order rate constant.

For the interphase mass transfer of CO2,assuming that the stable reactive absorption state is reached in the liquid film,the following equation can be obtained by the material balance of CO2in the liquid film.

Its boundary conditions are as follows.

According to Eqs.(11)–(13),the total mass transfer flux of CO2in the liquid phase can be obtained by solving Eq.(14) [29] as follows.

At the gas–liquid interface,fluxes of CO2through the gas phase and liquid phase should be consistent under the steady state,and it is shown as Eq.(15) as below.where Δx is the cell size perpendicular to the interface direction;is the concentration of CO2in the gas side at the phase interface.

Under the steady state condition,interfacial concentrations of CO2in the gas and liquid side are subject to Henry’s law:

where H is the Henry coefficient.is the concentration of carbon dioxide in the gas bulk,andis the concentration of carbon dioxide in the liquid bulk.Hence,through solving Eqs.(15)and(16)together,the interfacial CO2concentration on both sides of the gas and liquid phase along the interface and the mass transfer flux through the interface can be obtained.The corresponding mass source of CO2in the liquid can be expressed as below.

2.3.Physical property

Due to change of gas–liquid properties with the composition in the reactive absorption process of CO2by the MEA aqueous solution,the density,viscosity,diffusion coefficient,Henry coefficient of CO2in the MEA aqueous solution and reaction constant k2are modified according to the experimental correlation [30].

In addition,the liquid film in the microchannel is very thin,and generally a cooling system is utilized to surround the microchannel device,such as the experiment of Sobieszuk et al.[13].In this situation,the thermal effect of the reactive absorption has little influence on temperature.It is reasonable to assume that the temperature in the liquid film can keep constant.Therefore,the constant temperature condition in the wall can be set.Furthermore,the influence of the Rayleigh effect can be ignored in the falling film [9].Under the given 0.6 mm scale of the channel width in the model,the influence of the Laplace pressure on the phase equilibrium can be ignored[31],and the gas–liquid equilibrium can be represented by the common method.

The surface tension of the MEA–CO2solution can be calculated according to following experimental regression [32]:

where σaqis the surface tension of unloaded MEA aqueous solutions and σionstands for the contribution from ions.β is the CO2loading in the solution,which represents the total amount of CO2molecule absorbed by 1 mol MEA in the solution.ωMEAis the mass fraction of the molecular MEA.

The first item on the right of Eq.(18)represents the surface tension of the initial aqueous solution.The second item represents the contribution of ions to the surface tension due to the absorption of CO2.The third item represents influence of the hydrolysis of ions in the absorption solution as well as the action of ions with residual MEA and water on the surface tension.

For the aqueous MEA solution–CO2system,according to the established reaction mechanism and the gas–liquid phase equilibrium,the average error between the calculated CO2loadings in the phase equilibrium state and experimental value is no more than 5%[33,34].It shows that the relevant part of the model in this work is reliable.

For the absorption process of CO2by the MEA aqueous solution in the falling film microchannel,due to the narrow channel,the liquid film can form a meniscus in the channel under the solid–liquid interaction,so the gas–liquid mass transfer area can increase due to the presence of the meniscus.The shape of this meniscus depends on the contact angle of the liquid phase on the wall of the microchannel.According to the experiment of Sobieszuk et al.[13],the contact angle used in the numerical simulation is set to 15°,which also is satisfied with the general situation [18].

2.4.Grid division and numerical methods

In this work,structured hexahedral meshes are selected,and the total number of meshes is 561,600.Along the X-direction,the grid size in the 0.3 mm liquid film region close to the wall is 0.02 mm,and the grid size in the gas phase region gradually changes from 0.02 to 0.2 mm.The uniform sizes of the grid in the Y-direction and Z-direction are 0.25 and 0.02 mm respectively.

FLUENT 14.0 is used as the solver.The mass transfer equation is written by user define function (UDF) in FLUENT.The momentum source of the Marangoni force,the mass sources of the interface mass transfer and reaction are written by UDF.The physical properties,such as the density,viscosity,etc.,are also rewritten and calculated by UDF.The pressure velocity coupling equation uses the PISO method,and the spatial difference of the pressure term uses PRESTO! algorithm.The gas–liquid interface is obtained by the geo-reconstruct algorithm with high accuracy.The calculation precisions of the equations are 10–5while the unsteady time step is set to10–5s.

The Marangoni convection is unsteady in nature,so average values of the liquid outlet flowrate and concentration are taken as the reference.While they remain unchanged in the calculation process,the stable reactive absorption state can be acquired.Several mesh density schemes are adopted.With the increase in the mesh density,the difference of steady-state solutions under different mesh densities can disappear gradually,and the aforementioned grid scheme is confirmed by this strategy.

3.Model Validation

The enhancement effect of the Marangoni convection on the mass transfer is usually expressed by the enhancement factor of the mass transfer coefficient.It is the ratio of the mass transfer coefficient with and without the Marangoni convection.Sobieszuk et al.[13]only gave the enhancement factor of the Marangoni convection as the experimental result.Therefore,this work compares the computational enhancement factor with their experimental value in order to check whether the computational model agrees with the experiment under the same parameter conditions in the single channel.The reactive absorption is carried out at 298.15 K and 101.325 kPa.The volume fraction of CO2in the inlet is 23%.After conversion of original experimental flow rates,in the single channel,the gas inlet flow rate is set to 1.14×10–7m3.s-1while the liquid inlet flow rate is set to 1.31×10–8m3.s-1.

The total mass transfer flux of CO2can be calculated by Eq.(19)according to the simulation results.

where QV，inand QV，outare volume flowrates at the inlet and outlet of the gas phase respectively;b is the channel width;L is the channel length.

Then,the total mass transfer coefficient K is obtained as below.

The gasphase mass transfer coefficient kvcan becalculated by Eq.(21) [35] as below.

With and without consideration of the Marangoni convection,the liquid phase mass transfer coefficient is represented as kL，Maand kL，0respectively,so the Marangoni convection enhancement factor Ekon the liquid mass transfer coefficient is shown as Eq.(23).The enhancement factor given in this paper is the total value after absorption of the whole liquid film through the channel.

Fig.2.Enhancement factor at CO2 volume fraction of 23% in inlet.

As shown in Fig.2,when the volume fraction of CO2in the inlet is 23%,the Marangoni convection enhancement factors obtained by numerical simulation are basically in agreement with experimental results.The maximum deviation is 6.67% when the MEA solution is 1.5 mol.L–1.This proves the reliability of the numerical model and calculation established in this work.The higher the concentration of MEA at the liquid inlet is,the bigger the Marangoni enhancement factor is.The Marangoni enhancement factor can be about 2.5 at a high MEA concentration.With the increase of the MEA concentration,the reactive absorption becomes stronger,and the concentration of reaction products in the solution becomes higher.The change of the surface tension in the solution is sensitive to the change of the solution composition,so the Marangoni convection is more likely to occur.However,when the MEA concentration is bigger than 2 mol.L–1,the Marangoni enhancement factor fails to increase.It should be noted that the CO2concentration in the gas phase is rapidly consumed when the initial MEA concentration is bigger than 2 mol.L–1.Along the flowing direction of the liquid film,CO2is absorbed mainly in the very upper falling film.Along the whole falling film surface,the variation of the liquid composition,the corresponding concentration gradient as well as the surface tension gradient is not significant in the high MEA concentration,so it could result in a relatively weak intensity of the Marangoni convection.In addition,the liquid viscosity is also bigger at a higher MEA concentration.It also could result in a weak Marangoni convection as well as a smaller Marangoni enhancement factor.

4.Results and Discussion

4.1.Flow pattern and its effects

The numerical simulation of the CO2reactive absorption process in the falling film MEA aqueous solution is carried out under temperature,pressure,gas–liquid flow and gas concentration consistent with those in Section 3,and the initial concentration of the MEA aqueous solution is 1 mol.L–1.The streamline distribution in different cross-sections of the falling liquid film is shown in Fig.3,and it can represent the existence of the Marangoni convection in different heights of the falling film.It should be noted that the macroscopic flow direction of the liquid film is the negative direction of the Y-axis,and the Marangoni convection represented in Fig.3 actually exists on the X–Z plane.The velocity of twodimensional convective cells can be obtained by acquiring the Xand Z velocities in the X–Z plane.

At Y=75 mm,i.e.,near the liquid inlet,the cross-section of the falling film is near regular meniscus shape,and there is no trace of the Marangoni convection because there is no significant change of the liquid composition at the very beginning of absorption of CO2.With the increase of absorption of CO2in the falling film,at Y=50 mm,several Marangoni convective cells emerge near the gas–liquid interface,and the maximum velocity of convective cells can be about 0.007 m.s-1.Along the channel,with the gradual increase of reactive absorption and development of the Marangoni convection,the Marangoni convection cells gradually migrate towards two concave sides of the meniscus,and the convection gradually disappears at the center of the concave meniscus.At Y=5 mm,there are several convective cells in the liquid layer near the wall,and the maximum velocity of convective cells also decreases to 0.003 m.s-1.

Fig.3.Streamline and velocity distributions in different cross-sections.

As shown in Fig.3,the Marangoni convection occurs near the gas–liquid interface except near the liquid phase entrance.The presence of convective cells allows the inner fresh MEA aqueous solution to be convective to vicinity of the interface rather than diffusive,and it can significantly enhance the absorption amount of CO2.In addition,due to existence of the Marangoni convection,the significant deformation of the gas–liquid interface happens,and it can also enhance the absorption amount of CO2due to expansion of the gas–liquid contact area.

The corresponding quantitative results can be observed from the interfacial mass transfer flux shown in Fig.4 and the CO2loading shown in Fig.5.For the convenience of further analysis,Figs.4 and 5 also represent results without consideration of the Marangoni convection in the same conditions.

Fig.4.Mass transfer flux between gas and liquid interface.

At the entrance of the falling film,since the gas–liquid phase just starts to contact,the mass transfer driving force and the mass transfer flux is the largest as shown in Fig.4.Fig.4 is the projection of the gas–liquid interface on the Y–Z plane.However,near the entrance of the falling film,the total amount of CO2entering the liquid film and reacting with MEA is little at the very beginning of absorption,and the change of the liquid composition is limit.It has little influence on the surface tension.Therefore,near the entrance,there is little difference for the mass transfer flux between with and without the Marangoni convection.However,along the downward flow path,with the gradual adsorption of CO2,the amount of CO2entering the liquid film gradually increases.The composition in the liquid phase changes significantly,and the surface tension of the liquid phase also begins to change significantly.In the middle and lower region of the falling liquid film,the mass transfer flux with the Marangoni convection is obviously bigger than that without the Marangoni convection.The enhancement effect of the Marangoni convection on the mass transfer is significant in this region.

As shown in Fig.5,since the velocity along the flow direction of the liquid film near the wall is small,CO2and reaction products are continuously accumulated on the liquid film surface near the wall.The CO2loading near the wall region can approach to the saturation value of 0.67 mol.mol-1,and it is consistent with the data in the literature [32].If there is no existence of the Marangoni convection,the surface of the falling film is smooth and has a regular meniscus shape.In this situation,the CO2and reaction products are mainly concentrated on the surface of the liquid film,and little CO2can reach the inside of the falling film.However,with consideration of the Marangoni convection,the surface of the liquid film deforms to a certain degree.Furthermore,due to occurrence of the Marangoni convection,CO2and reaction products on the surface of the liquid film can be brought to the bottom layer of the liquid film while the fresh MEA inside the liquid film is brought to the vicinity of the interface.Therefore,an enhanced mass transfer flux can emerge in the middle and lower region of the falling film.

Fig.5.CO2 loading in liquid phase.

4.2.Suitable MEA concentration for absorption of flue gas

Since the Marangoni convection can obviously improve the absorption effect of CO2in the MEA aqueous solution,it is necessary to consider the effect of Marangoni convection on industrial flue gas conditions.In addition,as mentioned in Section 3,there is a complicated relationship between the Marangoni enhancement effect and the concentration of the MEA aqueous solution,so further investigation should be carried out.

Generally,the volume fraction of CO2in the coal-fired flue gas is 10%–15% [36].This means,in the industrial condition,CO2is absorbed in the low concentration.From the perspective of increasing CO2absorption rate,the high MEA concentration should be better.However,a high MEA concentration can cause corrosion of the equipment.At present,when the corrosion inhibitor is added,the MEA aqueous solution concentration is generally no more than 5 mol.L–1(mass fraction 30%) [37].

Fig.6.Marangoni enhancement effect and absorption ratio at different MEA concentrations.

The volume fraction of the gas phase CO2in the inlet is set to 15%,and other parameter conditions are the same as those in Section 3.In order to directly investigate the absorption effect of CO2,the absorption ratio of CO2is introduced.It is calculated by Eq.(24).Meanwhile,the enhancement effect of the Marangoni convection on the absorption ratio is represented by Eaas shown in Eq.(25).It is expressed by the ratio of absorption ratios with and without the Marangoni convection.Results are shown in Fig.6.

where XCO2is the absorption ratio of CO2and xCO2is the volume fraction of CO2.

As shown in Fig.6,in comparison with the case without the Marangoni convection,a higher CO2absorption ratio can be achieved when the Marangoni convection is considered.When the MEA concentration is lower than 2 mol.L–1,the higher the MEA concentration is,the stronger the enhancement effect of the Marangoni convection is.When the MEA concentration is more than 2 mol.L–1,the absorption ratio with the Marangoni convection is about 1.7 times bigger than that without the Marangoni convection.Although the Marangoni enhancement effect begins to decrease when the MEA concentration is bigger than 2 mol.L–1,even in this high concentration region,the absorption ratio of CO2is still obviously bigger than that without the Marangoni convection.As mentioned before,in the high MEA concentration,along the whole falling film surface,the variation of the liquid composition and the surface tension gradient is not significant due to fast adsorption.In addition,the liquid viscosity is bigger at a higher MEA concentration.These could result in a weak enhancement effect of the Marangoni convection.

Concerning the specific aim of the absorption ratio of CO2,a suitable MEA concentration can be chosen according to Fig.6.Obviously,a less MEA concentration can be utilized when the Marangoni convection is considered.It may have a significant economic advantage for application of the microchannel reactor in absorption of CO2.

It should be mentioned that results from this part represent enhanced absorption of CO2just for the microchannel configure and operation condition in this work.But for other configures and operation conditions,the similar methodology and procedure can be carried out in order to utilize the Marangoni convection to enhance absorption of CO2in the MEA aqueous solution.

5.Conclusions

For the microchannel falling film reactor,the 3D numerical simulation method established in this work can represent absorption of CO2in the MEA aqueous solution as well as the Marangoni convection.Results show that the Marangoni convection can occur in this reactive absorption process,and it can significantly promote surface renewal of the falling film and the interfacial mass transfer flux.Due to the influence of Marangoni convection,absorption of CO2can be significantly improved.With regard to the certain aim of the absorption ratio of CO2,a less MEA concentration can be utilized when the Marangoni convection is considered,and it may have a significant economic advantage.The quantitative results from this work is helpful for utilizing the Marangoni convection to enhance absorption of CO2in the MEA aqueous solution in the microchannel.

Nomenclature

A interface area,m2

a grid size,m

b channel width,m

c concentration,mol.m-3

D diffusion coefficient,m2.s-1

EaMarangoni convection enhancement factor on absorption ratio

EkMarangoni convection enhancement factor on liquid mass transfer coefficient

F surface tension per unit volume,N.m-3

g gravity acceleration,9.81 m.s-2

H Henry coefficient

HφHeaviside function

J diffusion flux,kg.m-2.s-1

K total mass transfer coefficient,m.s-1

kLliquidphase mass transfer coefficient,m.s-1

kVgasphase mass transfer coefficient,m.s-1

k2second-order rate constant,m3.mol-1.s-1

L channellength,m

N mass transfer flux,mol.m-2.s-1

p pressure,Pa

Q volume flowrate,m3.s-1

R reaction rate,mol.m-3.s-1

SCO2，Lmass source of CO2in the liquid,kg.m-3.s-1

T temperature,K

u velocity,m.s-1

XCO2absorption ratio of CO2

xCO2volume fraction of CO2

α volume fraction

β CO2loading,mol CO2.(mol MEA)–1

κ curvature,m-1

μ viscosity,Pa.s

σ surface tension,N.m-1

ρ density,kg.m-3

φ level set function,m

ω mass fraction

Subscripts

i component

I interface

L liquid phase

q phase q

V vapor phase

in inlet

out outlet

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

The authors acknowledge financial support provided by National Natural Science Foundation of China (21978243).