- Email: [email protected]
Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/he
Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel Mohamed El Amine Ben Amara*, Sassi Ben Nasrallah Energy and Thermal Systems Laboratory, National Engineering School of Monastir, University of Monastir, Avenue Ibn El Jazzar 5019, Monastir, Tunisia
article info
abstract
Article history:
Water management is emerging as one of the problems related to Proton Exchange
Received 15 June 2014
Membrane Fuel Cell. In fact, the humidity in a PEMFC plays a key role on its performance.
Received in revised form
The membrane must be sufficiently moistened to ensure the transport of protons. How-
12 August 2014
ever, liquid water may form and block the transport of gas to the electrodes. This can
Accepted 11 September 2014
generate a sharp decrease in current produced by the cell.
Available online xxx
In this paper, the droplet behavior in a proton exchange membrane (PEM) fuel cell micro-channel was simulated by using the lattice Boltzmann method (LBM) based on the
Keywords:
Shan-Chen Pseudo-potential model. A three-dimensional case was considered and a D3Q19
Lattice Boltzmann
scheme was utilized to keep track of the deformation of the liquidegas interface. Visual-
Two phase
ization of droplet shape is obtained for different capillary numbers and the hysteresis
Fuel cell
between the advancing and receding contact angle is clearly observed. Also flow structures
Micro-channel
in the micro-channel were illustrated. The effect of wettability on droplet displacement
Droplet
behavior is also explored. It was found that hydrophobic micro-channel is better than the hydrophilic micro-channel for droplets evacuation. This work presents a basic understanding for the droplet behavior in a fuel cell micro-channel and the effect of important parameters on its dynamics. Copyright © 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
Introduction A fuel cell is a device that generates electricity by a chemical reaction. It is the object of intensive research by the scientific community and industry [2]. There are different types of fuel cells. A fuel cell proton exchange membrane (PEM) is composed of an anode, a cathode, two bipolar plates, two catalysts and a membrane. The electricity generation is done
through oxidation on an electrode of a fuel in this case hydrogen (H2 / 2Hþ þ 2e) coupled with the reduction on the other electrode of an oxidant, such as the oxygen which produces water (1/2O2 þ 2Hþ þ 2e / H2O) The hydrogen and oxygen are fed to the anode and cathode by the flow channels and the gas-diffusion-layers (GDL) which is an essential element for the proper functioning of the fuel cell; it should allow the gas transfer on the electrode's entire surface, evacuate water from the active layers, evacuate the
* Corresponding author. Tel.: þ216 97575662. E-mail addresses: [email protected], [email protected] (M.E.A. Ben Amara). http://dx.doi.org/10.1016/j.ijhydene.2014.09.077 0360-3199/Copyright © 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
2
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Nomenclature cs ei fi eq fi G Gads i lu mu P Pin Pout R r T ts u U wi
speed of sound, lu ts1 particle speed, lu ts1 distribution function equilibrium distribution function cohesion parameter adhesion parameters lattice streaming vector direction lattice unit mass unit pressure, mu lu1 ts2 pressure inside the droplet, mu lu1 ts2 pressure outside of the drop, mu lu1 ts2 droplet radius, lu specific gas constant temperature, K time step macroscopic velocity characteristic velocity, lu ts1 weighting factor
Greek symbols U collision operator n kinematic viscosity, lu2 ts1 Dx lattice spacing, lu Dt time step, ts r density, mu lu3 j potential function t relaxation time, ts q contact angle m viscosity, mu lu1 ts1 s liquid surface tension, lu mu ts2 Abbreviations Bo bond number Ca capillary number Re Reynolds number
heat generated by the electrochemical reactions, improve the strength of the membrane electrode assembly and spreading electrons injected from the electrode in the plane of the active layer. The diffusion layer is a porous medium made of carbon fiber and usually treated with a hydrophobic agent to improve its water removal properties. It must have sufficient pores to provide reactant gas access from flow-field channels to catalyst layers and provide passage for removal of product water from catalyst-layer area to flow-field in channels which must have a well-defined shape. The water produced in the porous membrane is shunted out to the exterior of the fuel cell [17]. The polymer membrane must be hydrated to ensure high conductivity of protons. The bipolar plate containing the channels must have an adequate geometry for improving the performance of a fuel cell. PEM operating at reduced temperature and pressure levels, which triggers the development of two-phase flows inside the bipolar plates channels. There are many types of bipolar plates: the original GlobeTech geometry, the serpentine geometry,
the spiral geometry and the discontinuous channels geometry as cited by Yuan et al. [20]. The water management [10] in fuel cells is currently one of the critical issues. Thus, it is important to understand how water behaves inside the fuel cell. Removing the water produced by the electrochemical reactions at the cathode from anode side GDL is also essential to achieve continuous operation. The removal of the water produced by the fuel cell is carried out by mechanical entrainment at the same time as the evacuation of unconsumed reactant gas. The water present in the fuel cell can come from two sources: a portion of this water is fed through the wet gas stream entering the cell across the catalyst layer to the electrode-membrane interface, this water is then absorbed by the electrolyte and contributing to its hydration, when the second portion of water is produced by the chemical reactions of the fuels, in fact, after migrating through the membrane, the protons react with the electrons and oxygen gas, introduced on the cathode side, producing water. Water vapor condensation is responsible for the growth of water droplets in the feed channels of the fuel cell that are discharged by the gas flow. This phenomenon is important to analyze, several studies address the mechanisms of gas dynamics and its effects on the mobility of water droplets in the micro-channels. However, there are few direct observations of two-phase flow during fuel cell operation; some authors employ transparent fuel cell, magnetic resonance imaging, or X-rays to visualize the water. Nevertheless, these diagnostic tools require complex and expensive equipment. It is imperative to prevent fuel cell flooding, for this it is possible to optimize the stack design and the operating conditions. For example varying the shape of the bipolar plate channel can be a good solution for flooding problem [14]. Proton exchange membrane fuel cells micro-channels are generally rectangular in shape but there are other shapes such as trapezoidal, triangular, and circular shapes. The change of the micro-channel’s geometry may affect the water removal ability; for circular micro-channels, the condensed water forms a liquid film at the bottom of the channel and for tapered micro-channels the water forms small droplets. Liu et al. [15] studied the reactant gas transport and the cell performance of a proton exchange membrane fuel cell (PEMFC) with a tapered flow channel design, a two dimensional numerical modeling was done with the reduction in the channel depth along the streamwise direction. The authors concluded that the application of tapered micro-channels improve water management and fuel cell performance. The effects may be increased with decreasing the taper ratio of the fuel channel. Metz [18] proposes a flooding passive management using capillary forces to move the excess water. This requires a special structure of gas distribution micro-channels. He shows that the use of tapered channels may allow permanent gas intake even when large quantities of water are produced. The channels are designed so that to force water to rise along the hydrophilic walls then extended in a second microchannel located at the top of the conical section. For such a system to be effective it is essential to add an absorbent material to the upper end of each channel.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Buie [4] chooses an active water management by incorporating an electroosmotic pump to remove product water from the cathode area. Minor [19] has studied and modeled the droplet dynamic behavior in air flow field channels using experimental studies. Micro digital particle image velocimetry techniques were used to provide quantitative visualizations of the flow inside the liquid phase for the case of air flow around a droplet adhered to the wall of a rectangular gas channel model. The results of the study yield very interesting results such as the observation of a variety of rotational secondary flow patterns within the droplet. Using a transparent fuel cell, Chen [6] shows that several water flow regimes inside the microchannels can be observed. When the air flow increases, we can go from a regime with droplets stalling to film flow before reaching a chaotic flow regime. Hussaini et al. [13] also show in their work the existence of these different flow regimes. Anderson [1] also reviewed numerical simulations on the two-phase flow in gas flow channels of proton exchange membrane fuel cells. In his review, he exhibits in situ and ex situ experimental setups utilized to visualize and quantify two-phase flow phenomena in terms of flow regime maps, flow maldistribution, and pressure drop measurements. He summarizes the current literature on CFD simulations of gaseliquid two-phase flow in PEM fuel cells including the multi-fluid model, mixture model, volume of fraction method (VOF), and Lattice Boltzmann method (LBM). In his work, Cho et al. [8] did both an experimental and numerical study of droplet deformation and detachment in micro gas flow channels. The volume of fluid (VOF) scheme was used to simulate the water droplet dynamics and a high resolution CCD camera is employed to capture the droplet shape-change and detachment. They found that the local gas pressure is increased at the stagnation point in the front of the droplet and small droplets are more sensitive to viscous force. The volume-of-fluid model (VOF) was also applied by Zhu et al. [21] to investigate the effect of micro-channel geometry on droplet dynamics. Many geometric configurations were tested and it was shown that for the wettable configuration, the movement of water changes drastically when a droplet surface is in contact with the hydrophilic walls. In this work, a three-dimensional, two-phase model has been developed to simulate droplet emerging from pores of a porous fuel cell electrode into the gas flow channel as sketched in Fig. 1. Numerical algorithm based on lattice Boltzmann method is developed to solve the Boltzmann equation with the twophase ShaneChen model and appropriate boundary conditions. The model is capable of predicting the droplet dynamics in a micro-channel of a proton exchange membrane fuel cell.
Lattice Boltzmann method The Numerical Approach used in this work is based on the lattice Boltzmann method which is an alternative to the traditional approaches. Recently the lattice Boltzmann method has met with significant success for the numerical simulation of many
3
Fig. 1 e Sessile droplet inside a micro-channel.
industrial and scientific problems [3,5,7,9,11,12,16]. Traditional numerical methods solve the macroscopic transport equations of fluid flow, mass and heat transfer by directly discrediting them. Contrary to classical methods (difference methods (FDM), finite volume methods (FVM)), the LB method utilizes micro and meso-scale theoretical based on kinetic equations to recover the macroscopic NaviereStokes equation for fluid motion. LBM inherits many of the advantages of molecular dynamics and kinetic theories, due to its microscopic origin. The advantages of LBM include simple calculation procedure, simple and efficient implementation for parallel computation, easy and robust handling of complex geometries. The starting point of the lattice Boltzmann method is to solve, on a discrete lattice, the following Boltzmann equation for the discrete velocity distribution: / vfi þ ei V fi ¼ U fi vt
(1)
where ei is the speed of a particle at a position xi and time t and U(fi) is the collision operator controlling the rate of change in / the distribution function f during the collision, the term ! e i $V fi models the change in the distribution function due to the spread of the particles during their movement. The collision function represents the collision of fluid molecules at each node and has the following form (Bhatnagar et al.) [3]: eq Ui ¼ fi x; t fi x; t t
(2)
eq fi ðx; tÞ
is the equilibrium distribution function and t is where the relaxation time which is related to the viscosity of the fluid (n ¼ (2t1)/6, where n is the kinematic viscosity). The equilibrium distribution functions for different models were derived by He and Luo [12]. The function is given in the following form for the three-dimensional LB model with nine microscopic velocity vectors (D3Q19) (Fig. 2): " ! e i $! e i $! u ðx; tÞ ð! u ðx; tÞÞ2 eq x ; tÞ ¼ wi rðx; tÞ 1 þ þ fi ð! 2 cs 2c4s # ! u ðx; tÞ$! u ðx; tÞ 2c2s
(3)
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
4
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
uðx; tÞ ¼
19 X
fi ðx; tÞ ! e i =rðx; tÞ
(7)
i¼1
One of the main advantages of the lattice Boltzmann method is the easy introduction of boundary conditions; the most popular boundary condition for (LBE) method is the bounce back scheme, in this rule, particles which are incident upon a solid boundary are reversed and leave in the direction res et al. [9]). This method from which they came (D'Humie makes the resolution of complex solid boundaries straightforward. Fig. 3 shows the bounce back behavior of a solid wall. As the lattice Boltzmann method is a kinetic method, macroscopic boundary conditions do not have direct equivalents. They have to be replaced by appropriate microscopic rules which induce the desired macroscopic behavior. For the bounce back rule on wall nodes: x ; tÞ ¼ fi ð! x ; tÞ x2wall fi ð! Fig. 2 e The lattice velocities of D3Q19 model.
where r and ! u are the density and the macroscopic velocity of the node. wi is the weight factor for ith direction; we also assume that these factors are the same for directions having the same velocity, these variables satisfy the following relationP ship: Wi ¼ 1. For D3Q19, the weighting factors are defined by: i
8 < 1=3 i ¼ 0 Wi ¼ 1=18 i ¼ 1 to 6 : 1=36 i ¼ 7 to 18 cs is the isothermal speed of sound, it is obtained by the ChapmaneEnskog expansion. pffiffiffi Cs ¼ 1= 3 for D2Q9 lattice (Chen et al., [5]; Frisch et al., [10]; He and Luo, [12]). In three-dimensional problems, the most popular scheme is the D3Q19 model (See Fig. 2 below). In this model, the lattice velocities are defined as: ! ei ¼
(8)
with ! e i ¼ ! e i.
8 <0 i ¼ 0 ð±1; 0; 0Þ; ð0; ±1; 0Þ; ð0; 0; ±1Þ i ¼ 1 to 6 : ð±1; ±1; 0Þ; ð±1; 0; ±1Þ; ð0; ±1; ±1Þ i ¼ 7 to 18
Shan and Chen-type lattice Boltzmann Shan and Chen proposed a multiple phases LBM model by introducing an interparticle potential between fluid components and based on the BGK collision model. In this model, one distribution function is introduced for each of the fluid components. In the Shan-Chen model, a force, between the two fluids is introduced that effectively perturbs the equilibrium velocity for each fluid [16]. In D3Q19 model, this force is given by: Fðx; tÞ ¼ Gjðx; tÞ
X
wi jðx þ ei Dt; tÞei
(9)
i
where G is the interaction strength, wi is weight coefficient, and j is the interaction potential. By choosing a potential function of the form: j r ¼ j0 1 eðr0 =rÞ
(10)
The discretized Lattice Boltzmann equation is as follow: eq fi ð! x þ! e i $Dt; t þ DtÞ fi ð! x ; tÞ ¼ 1=t fi ð! x ; tÞ fi ð! x ; tÞ
(5)
There are two basic recurrent steps during simulation of viscous flow in a typical LB algorithm: -Collision step: the arriving particles at the points interact with another and change their velocity directions, thus, at time t the particles at node x come into collision with each x ; tÞ to other which changes the distribution function from fi ð! eq x ; tÞ ¼ fi ð! x ; tÞ þ 1t ðfi ð! x ; tÞ fi ð! x ; tÞÞ. fi* ð! -Streaming step: particles move during the time step Dt, along lattice bonds to the neighboring lattice nodes and the x ; tÞ spreads along the vector ! e i, distribution function fi* ð! x þ! e i $Dt; t þ DtÞ ¼ fi* ð! x ; tÞ. more formally: fi ð! The two macroscopic properties, density (r) and velocity (u) of the nodes, are calculated using the following relations: rðx; tÞ ¼
9 X i¼1
fi ðx; tÞ
(6) Fig. 3 e Rebound conditions at a wall.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
5
Fig. 7 e Impact of Gads on the shape of the droplets: (a) q ≈ 177 .5, (b) q ≈ 90, (c) q ≈ 17 .8.
Fig. 4 e Plot of pressure vs density at different values of G based on ShaneChen model.
Fig. 8 e Plot of the contact angle versus Gads.
Fig. 5 e Instantaneous phase separation.
It is possible to obtain a liquid vapor equilibrium whose coexisting densities can be derived from Maxwell construction (equal area construction). Here, j0 and r0 are arbitrary constants. The equation of state for the ShaneChen model is given by applying the ChapmanneEnskog expansion [15] the pressure is related to the density as: P ¼ rRT þ GRT=2j2 ðrÞ. Fig. 4 shows pressure plots for different values of G. The shapes of these curves are similar to the typical curve of van der Waals equation of state used to predict the interaction potential in the ShaneChan model. According to this
model, the phase transition occurs when the value of G is greater or equal than a critical value Gc ¼ 92.4. Adhesive forces between the fluid and solid phases are introduced into the model by Martys and Chen, [6]: Fads ðx; tÞ ¼ Gads jðx; tÞ
X
wi sðx þ ei Dt; tÞei
(11)
i
Here s ¼ 0, 1 for nodes in the liquid and on solid walls, respectively. Gads represents the particle interaction strength between fluid and solid walls, and varying the parameter allows simulation of the complete range of contact angles. With these definitions, in simulation, the cohesive force and the attractive force are added to the velocities that compute the equilibrium distribution function with the following formula:
Fig. 6 e Pressure differences across the interface for different drop radii from LBM simulations. Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
6
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Results and discussion
Fig. 9 e Geometry of the LBM experiment consisting in a droplet in a closed rectangular micro-channel.
t ueq ¼ u þ ðF þ Fads Þ r
(12)
We can remark that these forces are not explicitly taken into account in the discrete lattice Boltzmann equation.
In this section, some simulation results are repeated to examine the validity of this model and code. First, a simulation of a spinodal phase separation has been conducted in a 51 51 51 lu3 domain initialized by a random distribution of two phases; the system will automatically separates to the liquid phase and the vapor phase due to the instability of van der Waals equation of state. Fig. 5 shows the Phase field morphology at different times. The variation of the density is shown in color scale with the maximum in red (in web version) and the minimum in blue (in web version). Second, the ShaneChan model is applied to verify the Laplace law for static drops. This law states that the pressure difference required to maintain the shape of an interface between two fluids is proportional to the interface curvature. Considering a spherical drop within a three-dimensional domain, the Laplace law can be written as: DP ¼ Pin Pout ¼ s/R.
Fig. 10 e Sketch of three-dimensional of a droplet sitting on a surface in a small gas channel exposed to a gas flow moving right to left at Ca ¼ 0.083. The initial contact angle q at equilibrium is 165.38 . Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Fig. 11 e Streamlines behind the droplet in XZ plane at Ca ¼ 0.083.
Where DP is the pressure difference across the fluid interface s is the surface tension, R is the radii of curvature, Pin and Pout are the pressures inside and outside of the drop, respectively. The simulation was conducted for several initial droplet radii. Then the changes of pressure differences with respect of 1/R are plotted. Drops with different radii are generated inside a computational domain using the ShaneChan model. A plot of P versus 1/R is a straight line and its slope will be s as it is seen in Fig. 6. For Dr ¼ 438.686 mu lu2, G ¼ 120, J0 ¼ 4 and r0 ¼ 200 mu lu2, DP ¼ 14.473 (1/R) þ 0.00181, the surface tension is equal to 14.473 lu mu ts2. We can conclude that Laplace's Law is well satisfied for this model. The equilibrium contact angle, q defined as the angle between the vaporeliquid interface and the solid wall, can be obtained through numerical experiment of droplet placement on solid surfaces with different parameter Gads. Fig. 7 shows droplet shapes for various fluid-solid interaction strengths. Three typed of solid surfaces were defined with Gads ¼ 50, 187.16 and 300 corresponding to q ¼ 177 , 90 and 17 0.8. In Fig. 8 we report contact angle versus fluid-solid interaction parameter Gads.
7
It is found that contact angle varies al-most linearly with the fluid-solid interaction strength Gads. A linear relation between Gads and q can be given by: q ¼ 0.64Gads þ 209.8. We implemented the lattice Boltzmann model for nonideal fluids to simulate the droplet dynamics in a microchannel of a proton exchange membrane fuel. Initially, the three-dimensional computation domain consists of an incompressible fluid droplet in a rectangular micro-channel (See Fig. 9). The two steps “stream and collide” algorithm for a (D3Q19) lattice is used to simulate lattice Boltzmann equation on 150 40 30 site lattices. Bounce back boundary conditions are imposed on the walls and periodic boundary conditions are imposed on the two domain ends. To mimic the flow inside the micro-channel, a constant body force along the z direction was applied to displace both the fluid and the droplet. Driven by this force the droplet move but its speed and shape change with the inlet velocity and the wettability. There are several dimensionless numbers encountered in micro-fluidic applications. For low Reynolds number regime (fluid inertia forces are negligible), the displacement of immiscible fluids in micro-fluidic devices can be described by two main dimensionless numbers: the Capillary number (Ca) and Bond number (Bo). The capillary number reflects the competition between viscous forces and capillary forces. It is given by: Ca ¼ mU/s. Where s is the surface tension which is calculated using Laplace's law test, U is a characteristic velocity and m is the dynamic viscosity of the liquid phase. The Bond number reflects the balance between gravitational and capillary forces and is defined as: Bo ¼ R2Drg/s. where g is the gravitational acceleration. In this work, the focus will be on the influence of capillary number. Calculation of Ca requires the viscosity which is given by:y ¼ c2s ðt 1=2Þ where cs is the sound velocity and t is the BGK relaxation time. We focus on the dynamic characteristics of liquid drop's spreading. We will present essentially the velocity and density fields and streamlines. As a first simulation step, we have simulated droplet motion through a 150 40 30 voxels micro-channel. Normalized density plots at various time steps were generated. In order to focus on the effect of shear the shear induced by the gas on the droplet, the physical parameters are chosen so that the gravity is negligible compared with the surface tension This means that the Bond number is less than 1. The droplet is initially spherical caps centered on the pore due to surface tension, after that the droplet moves down the micro-channel in the direction of the gas flow. The droplet will deform as a
Fig. 12 e Streamlines generated from 2D mean velocity vectors in YZ (a) and XY (b) planes (Ca ¼ 0.083). Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
8
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Fig. 13 e Evolution of the interfaces of droplets for different capillary numbers.
result of these forces and exhibit a contact angle hysteresis defined as the difference between the advancing and receding contact angles. Fig. 10 shows the evolution of the phase distribution patterns under external force (The fluid flow is always laminar Re z 1.9). After a certain time step, the droplet continues to move at constant speed and reaches equilibrium with an unchanged shape. The streamlines are plotted to show the flow field in the micro-channel. Flow recirculations behind the droplet and near the contact point are clearly visible in Fig. 11. Velocity streamlines in XY and YZ planes are presented in Fig. 12. Flow recirculation around the bubble is observed. Fig. 13 illustrates the effect of capillary number on the deformation of droplets with q ¼ 165.38 , this comparison shows that as the capillary number increases, the droplet displacement and deformation also increase. Fig. 14 shows the high deformation of droplet under high capillary number. The droplet is strongly deformed at Ca ¼ 0.2075 but the contact surface area decreases. We note that for low capillary number the droplet cannot be evacuated from the micro-channel. The streamlines plot show a symmetrical vortex around the deformed droplet (See Fig. 15).
Fig. 15 e Vortex around the deformed droplet at Ca ¼ 0.2075.
We consider now the wetting case, Fig. 16 shows the droplet behavior in micro-channel with hydrophilic surface (q < 90 ). Because of shearing forces, the droplet breaks down into two droplets with different sizes: a principal droplet and a satellite droplet which touches the top wall, this droplet is drained and disappears completely. The wettability effect is very strong and can transform droplets in thin liquid films spreading along the solid walls. We can remark that the droplet velocity is relatively high for a hydrophobic channel walls than for a hydrophilic walls.
Conclusion We have developed an LBM algorithm with an external force for two-phase droplet dynamics in three-dimensional microchannel of a proton exchange membrane fuel, in this model the two phase fluid is modeled by ShaneChen scheme. The density is simulated inside the micro-channel. We remark that the droplet is highly deformed into a shape with a tip. It is observed that the droplet deformation increases with increase in capillary number.
Fig. 14 e Droplet behavior under high capillary number Ca ¼ 0.2075. Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
9
Fig. 16 e Time evolution of dynamic behavior of droplet in micro-channel with hydrophilic surface. (q < 90 ).
It is clear that the average contact angle hysteresis, which is the difference between the advancing and receding contact angles, increases during the droplet motion. In the wetting case, droplet adheres more to the channel walls and its evacuation became more difficult than the case of a hydrophobic micro-channel. Recirculation in the wake behind the droplet is observed which depend on the capillary number Ca and the wettability. Other geometries must be tested and porous effect must be included in future work.
references
[1] Anderson R, Zhang L, Ding Y, Blanco M, Bi X, Wilkinson DP. A critical review of two-phase flow in gas flow channels of
[2] [3]
[4]
[5]
[6]
[7]
proton exchange membrane fuel cells. J Power Sources 2010;195:4531e53. Barbir F. PEM fuel cells: theory and practice. Elsevier Academic Press; 2005. Bhatnagar PL, Gross EP, Krook M. A model for collision processes in gases, I. Small amplitude processes in charged and neutral one-component system. Phys Rev 1954;94:511e25. Buie CR, Posner JD, Fabian T, Cha SW, Kim D, Prinz FB, et al. Water management in proton exchange membrane fuel cells using integrated electroosmotic pumping. J Power Sources 2006;161(1):191e202. Chen H, Chen S, Matthaeus WH. Recovery of the NavierStokes equations using a lattice-gas Boltzmann method. Phys Rev A 1992;45:5339e42. Chen J. Experimental study on the two phase flow behavior in PEM fuel cell parallel channels with porous media inserts. J Power Sources 2010;195:1122e9. Chen S, Doolen GD. Lattice Boltzmann method for fluid flows. Annu Rev Fluid Mech 1998;30:329e64.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
10
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
[8] Cho SC, Wang Y, Chen KS. Droplet dynamics in a polymer electrolyte fuel cell gas flow channel: forces, deformation, and detachment. I: theoretical and numerical analyses. J Power Sources 2012;206:119e28. [9] D'Humieres D, Lallemand P. Numerical simulations of hydrodynamics with lattice gas automata in two dimensions. Complex Syst 1987;1:598. [10] Eikerling M. Water management in cathode catalyst layers of PEM fuel cells. J Electrochem Soc 2006;153(3):E58. [11] Frisch U, Hasslacher B, Pomeau Y. Lattice-gas automata for Navier-Stokes equations. Phys Rev Lett 1986;56:1505e8. [12] He X, Luo LS. Theory of the lattice Boltzmann method: from the Boltzmann equation to the lattice Boltzmann equation. Phys Rev E. 1997;56:6811e8. [13] Hussaini IS, Wang CY. Visualization and quantification of cathode channel flooding in PEM fuel cells. J Power Sources 2009;187:444e51. [14] Li H, Tang Y, Wang Z, Shi Z, Wu S, Song D, et al. A review of water flooding issues in the proton exchange membrane fuel cell. J Power Sources 2008;178:103e17. [15] Liu HC, Yan WM, Soong CY, Chenc F, Chud HS. Reactant gas transport and cell performance of proton exchange
[16]
[17] [18]
[19]
[20]
[21]
membrane fuel cells with tapered flow field design. J Power Sources 2006;158:78e87. Martys NS, Chen H. Simulation of multicomponent fluids in complex three-dimensional geometries by the lattice Boltzmann method. Phys Rev E. 1996;53:743e50. Mench M. Fuel cell engines. John Wiley & Sons, Inc; 2008. Metz T, Paust N, Mu¨ller C, Zengerle R, Koltay P. Passive water removal in fuel cells by capillary droplet actuation. Sens Actuators A e Phys July 2007;5 [in press], Corrected Proof. Minor G, Zhu X, Oshkai P, Sui P, Djilali N. Water transport dynamics in fuel cell micro-channels. In: Kakac S, Pramuanjaroenkij A, Vasiliev L, editors. Fundementals and applications. Dordrecht: Springer; 2008. Mini-Micro Fuel 75 Cells. Yuan XZ, Wang H, Zhang J, Wilkinson DP. Bipolar plates for PEM fuel cells-from materials to processing. J New Mater Electrochem Syst 2005;8(4):257e67. Zhu X, Liao Q, Sui PC, Djilali N. Numerical investigation of water droplet dynamics in a low-temperature fuel cell microchannel: effect of channel geometry. J Power Sources 2010;195:801e12.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/he
Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel Mohamed El Amine Ben Amara*, Sassi Ben Nasrallah Energy and Thermal Systems Laboratory, National Engineering School of Monastir, University of Monastir, Avenue Ibn El Jazzar 5019, Monastir, Tunisia
article info
abstract
Article history:
Water management is emerging as one of the problems related to Proton Exchange
Received 15 June 2014
Membrane Fuel Cell. In fact, the humidity in a PEMFC plays a key role on its performance.
Received in revised form
The membrane must be sufficiently moistened to ensure the transport of protons. How-
12 August 2014
ever, liquid water may form and block the transport of gas to the electrodes. This can
Accepted 11 September 2014
generate a sharp decrease in current produced by the cell.
Available online xxx
In this paper, the droplet behavior in a proton exchange membrane (PEM) fuel cell micro-channel was simulated by using the lattice Boltzmann method (LBM) based on the
Keywords:
Shan-Chen Pseudo-potential model. A three-dimensional case was considered and a D3Q19
Lattice Boltzmann
scheme was utilized to keep track of the deformation of the liquidegas interface. Visual-
Two phase
ization of droplet shape is obtained for different capillary numbers and the hysteresis
Fuel cell
between the advancing and receding contact angle is clearly observed. Also flow structures
Micro-channel
in the micro-channel were illustrated. The effect of wettability on droplet displacement
Droplet
behavior is also explored. It was found that hydrophobic micro-channel is better than the hydrophilic micro-channel for droplets evacuation. This work presents a basic understanding for the droplet behavior in a fuel cell micro-channel and the effect of important parameters on its dynamics. Copyright © 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
Introduction A fuel cell is a device that generates electricity by a chemical reaction. It is the object of intensive research by the scientific community and industry [2]. There are different types of fuel cells. A fuel cell proton exchange membrane (PEM) is composed of an anode, a cathode, two bipolar plates, two catalysts and a membrane. The electricity generation is done
through oxidation on an electrode of a fuel in this case hydrogen (H2 / 2Hþ þ 2e) coupled with the reduction on the other electrode of an oxidant, such as the oxygen which produces water (1/2O2 þ 2Hþ þ 2e / H2O) The hydrogen and oxygen are fed to the anode and cathode by the flow channels and the gas-diffusion-layers (GDL) which is an essential element for the proper functioning of the fuel cell; it should allow the gas transfer on the electrode's entire surface, evacuate water from the active layers, evacuate the
* Corresponding author. Tel.: þ216 97575662. E-mail addresses: [email protected], [email protected] (M.E.A. Ben Amara). http://dx.doi.org/10.1016/j.ijhydene.2014.09.077 0360-3199/Copyright © 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
2
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Nomenclature cs ei fi eq fi G Gads i lu mu P Pin Pout R r T ts u U wi
speed of sound, lu ts1 particle speed, lu ts1 distribution function equilibrium distribution function cohesion parameter adhesion parameters lattice streaming vector direction lattice unit mass unit pressure, mu lu1 ts2 pressure inside the droplet, mu lu1 ts2 pressure outside of the drop, mu lu1 ts2 droplet radius, lu specific gas constant temperature, K time step macroscopic velocity characteristic velocity, lu ts1 weighting factor
Greek symbols U collision operator n kinematic viscosity, lu2 ts1 Dx lattice spacing, lu Dt time step, ts r density, mu lu3 j potential function t relaxation time, ts q contact angle m viscosity, mu lu1 ts1 s liquid surface tension, lu mu ts2 Abbreviations Bo bond number Ca capillary number Re Reynolds number
heat generated by the electrochemical reactions, improve the strength of the membrane electrode assembly and spreading electrons injected from the electrode in the plane of the active layer. The diffusion layer is a porous medium made of carbon fiber and usually treated with a hydrophobic agent to improve its water removal properties. It must have sufficient pores to provide reactant gas access from flow-field channels to catalyst layers and provide passage for removal of product water from catalyst-layer area to flow-field in channels which must have a well-defined shape. The water produced in the porous membrane is shunted out to the exterior of the fuel cell [17]. The polymer membrane must be hydrated to ensure high conductivity of protons. The bipolar plate containing the channels must have an adequate geometry for improving the performance of a fuel cell. PEM operating at reduced temperature and pressure levels, which triggers the development of two-phase flows inside the bipolar plates channels. There are many types of bipolar plates: the original GlobeTech geometry, the serpentine geometry,
the spiral geometry and the discontinuous channels geometry as cited by Yuan et al. [20]. The water management [10] in fuel cells is currently one of the critical issues. Thus, it is important to understand how water behaves inside the fuel cell. Removing the water produced by the electrochemical reactions at the cathode from anode side GDL is also essential to achieve continuous operation. The removal of the water produced by the fuel cell is carried out by mechanical entrainment at the same time as the evacuation of unconsumed reactant gas. The water present in the fuel cell can come from two sources: a portion of this water is fed through the wet gas stream entering the cell across the catalyst layer to the electrode-membrane interface, this water is then absorbed by the electrolyte and contributing to its hydration, when the second portion of water is produced by the chemical reactions of the fuels, in fact, after migrating through the membrane, the protons react with the electrons and oxygen gas, introduced on the cathode side, producing water. Water vapor condensation is responsible for the growth of water droplets in the feed channels of the fuel cell that are discharged by the gas flow. This phenomenon is important to analyze, several studies address the mechanisms of gas dynamics and its effects on the mobility of water droplets in the micro-channels. However, there are few direct observations of two-phase flow during fuel cell operation; some authors employ transparent fuel cell, magnetic resonance imaging, or X-rays to visualize the water. Nevertheless, these diagnostic tools require complex and expensive equipment. It is imperative to prevent fuel cell flooding, for this it is possible to optimize the stack design and the operating conditions. For example varying the shape of the bipolar plate channel can be a good solution for flooding problem [14]. Proton exchange membrane fuel cells micro-channels are generally rectangular in shape but there are other shapes such as trapezoidal, triangular, and circular shapes. The change of the micro-channel’s geometry may affect the water removal ability; for circular micro-channels, the condensed water forms a liquid film at the bottom of the channel and for tapered micro-channels the water forms small droplets. Liu et al. [15] studied the reactant gas transport and the cell performance of a proton exchange membrane fuel cell (PEMFC) with a tapered flow channel design, a two dimensional numerical modeling was done with the reduction in the channel depth along the streamwise direction. The authors concluded that the application of tapered micro-channels improve water management and fuel cell performance. The effects may be increased with decreasing the taper ratio of the fuel channel. Metz [18] proposes a flooding passive management using capillary forces to move the excess water. This requires a special structure of gas distribution micro-channels. He shows that the use of tapered channels may allow permanent gas intake even when large quantities of water are produced. The channels are designed so that to force water to rise along the hydrophilic walls then extended in a second microchannel located at the top of the conical section. For such a system to be effective it is essential to add an absorbent material to the upper end of each channel.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Buie [4] chooses an active water management by incorporating an electroosmotic pump to remove product water from the cathode area. Minor [19] has studied and modeled the droplet dynamic behavior in air flow field channels using experimental studies. Micro digital particle image velocimetry techniques were used to provide quantitative visualizations of the flow inside the liquid phase for the case of air flow around a droplet adhered to the wall of a rectangular gas channel model. The results of the study yield very interesting results such as the observation of a variety of rotational secondary flow patterns within the droplet. Using a transparent fuel cell, Chen [6] shows that several water flow regimes inside the microchannels can be observed. When the air flow increases, we can go from a regime with droplets stalling to film flow before reaching a chaotic flow regime. Hussaini et al. [13] also show in their work the existence of these different flow regimes. Anderson [1] also reviewed numerical simulations on the two-phase flow in gas flow channels of proton exchange membrane fuel cells. In his review, he exhibits in situ and ex situ experimental setups utilized to visualize and quantify two-phase flow phenomena in terms of flow regime maps, flow maldistribution, and pressure drop measurements. He summarizes the current literature on CFD simulations of gaseliquid two-phase flow in PEM fuel cells including the multi-fluid model, mixture model, volume of fraction method (VOF), and Lattice Boltzmann method (LBM). In his work, Cho et al. [8] did both an experimental and numerical study of droplet deformation and detachment in micro gas flow channels. The volume of fluid (VOF) scheme was used to simulate the water droplet dynamics and a high resolution CCD camera is employed to capture the droplet shape-change and detachment. They found that the local gas pressure is increased at the stagnation point in the front of the droplet and small droplets are more sensitive to viscous force. The volume-of-fluid model (VOF) was also applied by Zhu et al. [21] to investigate the effect of micro-channel geometry on droplet dynamics. Many geometric configurations were tested and it was shown that for the wettable configuration, the movement of water changes drastically when a droplet surface is in contact with the hydrophilic walls. In this work, a three-dimensional, two-phase model has been developed to simulate droplet emerging from pores of a porous fuel cell electrode into the gas flow channel as sketched in Fig. 1. Numerical algorithm based on lattice Boltzmann method is developed to solve the Boltzmann equation with the twophase ShaneChen model and appropriate boundary conditions. The model is capable of predicting the droplet dynamics in a micro-channel of a proton exchange membrane fuel cell.
Lattice Boltzmann method The Numerical Approach used in this work is based on the lattice Boltzmann method which is an alternative to the traditional approaches. Recently the lattice Boltzmann method has met with significant success for the numerical simulation of many
3
Fig. 1 e Sessile droplet inside a micro-channel.
industrial and scientific problems [3,5,7,9,11,12,16]. Traditional numerical methods solve the macroscopic transport equations of fluid flow, mass and heat transfer by directly discrediting them. Contrary to classical methods (difference methods (FDM), finite volume methods (FVM)), the LB method utilizes micro and meso-scale theoretical based on kinetic equations to recover the macroscopic NaviereStokes equation for fluid motion. LBM inherits many of the advantages of molecular dynamics and kinetic theories, due to its microscopic origin. The advantages of LBM include simple calculation procedure, simple and efficient implementation for parallel computation, easy and robust handling of complex geometries. The starting point of the lattice Boltzmann method is to solve, on a discrete lattice, the following Boltzmann equation for the discrete velocity distribution: / vfi þ ei V fi ¼ U fi vt
(1)
where ei is the speed of a particle at a position xi and time t and U(fi) is the collision operator controlling the rate of change in / the distribution function f during the collision, the term ! e i $V fi models the change in the distribution function due to the spread of the particles during their movement. The collision function represents the collision of fluid molecules at each node and has the following form (Bhatnagar et al.) [3]: eq Ui ¼ fi x; t fi x; t t
(2)
eq fi ðx; tÞ
is the equilibrium distribution function and t is where the relaxation time which is related to the viscosity of the fluid (n ¼ (2t1)/6, where n is the kinematic viscosity). The equilibrium distribution functions for different models were derived by He and Luo [12]. The function is given in the following form for the three-dimensional LB model with nine microscopic velocity vectors (D3Q19) (Fig. 2): " ! e i $! e i $! u ðx; tÞ ð! u ðx; tÞÞ2 eq x ; tÞ ¼ wi rðx; tÞ 1 þ þ fi ð! 2 cs 2c4s # ! u ðx; tÞ$! u ðx; tÞ 2c2s
(3)
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
4
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
uðx; tÞ ¼
19 X
fi ðx; tÞ ! e i =rðx; tÞ
(7)
i¼1
One of the main advantages of the lattice Boltzmann method is the easy introduction of boundary conditions; the most popular boundary condition for (LBE) method is the bounce back scheme, in this rule, particles which are incident upon a solid boundary are reversed and leave in the direction res et al. [9]). This method from which they came (D'Humie makes the resolution of complex solid boundaries straightforward. Fig. 3 shows the bounce back behavior of a solid wall. As the lattice Boltzmann method is a kinetic method, macroscopic boundary conditions do not have direct equivalents. They have to be replaced by appropriate microscopic rules which induce the desired macroscopic behavior. For the bounce back rule on wall nodes: x ; tÞ ¼ fi ð! x ; tÞ x2wall fi ð! Fig. 2 e The lattice velocities of D3Q19 model.
where r and ! u are the density and the macroscopic velocity of the node. wi is the weight factor for ith direction; we also assume that these factors are the same for directions having the same velocity, these variables satisfy the following relationP ship: Wi ¼ 1. For D3Q19, the weighting factors are defined by: i
8 < 1=3 i ¼ 0 Wi ¼ 1=18 i ¼ 1 to 6 : 1=36 i ¼ 7 to 18 cs is the isothermal speed of sound, it is obtained by the ChapmaneEnskog expansion. pffiffiffi Cs ¼ 1= 3 for D2Q9 lattice (Chen et al., [5]; Frisch et al., [10]; He and Luo, [12]). In three-dimensional problems, the most popular scheme is the D3Q19 model (See Fig. 2 below). In this model, the lattice velocities are defined as: ! ei ¼
(8)
with ! e i ¼ ! e i.
8 <0 i ¼ 0 ð±1; 0; 0Þ; ð0; ±1; 0Þ; ð0; 0; ±1Þ i ¼ 1 to 6 : ð±1; ±1; 0Þ; ð±1; 0; ±1Þ; ð0; ±1; ±1Þ i ¼ 7 to 18
Shan and Chen-type lattice Boltzmann Shan and Chen proposed a multiple phases LBM model by introducing an interparticle potential between fluid components and based on the BGK collision model. In this model, one distribution function is introduced for each of the fluid components. In the Shan-Chen model, a force, between the two fluids is introduced that effectively perturbs the equilibrium velocity for each fluid [16]. In D3Q19 model, this force is given by: Fðx; tÞ ¼ Gjðx; tÞ
X
wi jðx þ ei Dt; tÞei
(9)
i
where G is the interaction strength, wi is weight coefficient, and j is the interaction potential. By choosing a potential function of the form: j r ¼ j0 1 eðr0 =rÞ
(10)
The discretized Lattice Boltzmann equation is as follow: eq fi ð! x þ! e i $Dt; t þ DtÞ fi ð! x ; tÞ ¼ 1=t fi ð! x ; tÞ fi ð! x ; tÞ
(5)
There are two basic recurrent steps during simulation of viscous flow in a typical LB algorithm: -Collision step: the arriving particles at the points interact with another and change their velocity directions, thus, at time t the particles at node x come into collision with each x ; tÞ to other which changes the distribution function from fi ð! eq x ; tÞ ¼ fi ð! x ; tÞ þ 1t ðfi ð! x ; tÞ fi ð! x ; tÞÞ. fi* ð! -Streaming step: particles move during the time step Dt, along lattice bonds to the neighboring lattice nodes and the x ; tÞ spreads along the vector ! e i, distribution function fi* ð! x þ! e i $Dt; t þ DtÞ ¼ fi* ð! x ; tÞ. more formally: fi ð! The two macroscopic properties, density (r) and velocity (u) of the nodes, are calculated using the following relations: rðx; tÞ ¼
9 X i¼1
fi ðx; tÞ
(6) Fig. 3 e Rebound conditions at a wall.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
5
Fig. 7 e Impact of Gads on the shape of the droplets: (a) q ≈ 177 .5, (b) q ≈ 90, (c) q ≈ 17 .8.
Fig. 4 e Plot of pressure vs density at different values of G based on ShaneChen model.
Fig. 8 e Plot of the contact angle versus Gads.
Fig. 5 e Instantaneous phase separation.
It is possible to obtain a liquid vapor equilibrium whose coexisting densities can be derived from Maxwell construction (equal area construction). Here, j0 and r0 are arbitrary constants. The equation of state for the ShaneChen model is given by applying the ChapmanneEnskog expansion [15] the pressure is related to the density as: P ¼ rRT þ GRT=2j2 ðrÞ. Fig. 4 shows pressure plots for different values of G. The shapes of these curves are similar to the typical curve of van der Waals equation of state used to predict the interaction potential in the ShaneChan model. According to this
model, the phase transition occurs when the value of G is greater or equal than a critical value Gc ¼ 92.4. Adhesive forces between the fluid and solid phases are introduced into the model by Martys and Chen, [6]: Fads ðx; tÞ ¼ Gads jðx; tÞ
X
wi sðx þ ei Dt; tÞei
(11)
i
Here s ¼ 0, 1 for nodes in the liquid and on solid walls, respectively. Gads represents the particle interaction strength between fluid and solid walls, and varying the parameter allows simulation of the complete range of contact angles. With these definitions, in simulation, the cohesive force and the attractive force are added to the velocities that compute the equilibrium distribution function with the following formula:
Fig. 6 e Pressure differences across the interface for different drop radii from LBM simulations. Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
6
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Results and discussion
Fig. 9 e Geometry of the LBM experiment consisting in a droplet in a closed rectangular micro-channel.
t ueq ¼ u þ ðF þ Fads Þ r
(12)
We can remark that these forces are not explicitly taken into account in the discrete lattice Boltzmann equation.
In this section, some simulation results are repeated to examine the validity of this model and code. First, a simulation of a spinodal phase separation has been conducted in a 51 51 51 lu3 domain initialized by a random distribution of two phases; the system will automatically separates to the liquid phase and the vapor phase due to the instability of van der Waals equation of state. Fig. 5 shows the Phase field morphology at different times. The variation of the density is shown in color scale with the maximum in red (in web version) and the minimum in blue (in web version). Second, the ShaneChan model is applied to verify the Laplace law for static drops. This law states that the pressure difference required to maintain the shape of an interface between two fluids is proportional to the interface curvature. Considering a spherical drop within a three-dimensional domain, the Laplace law can be written as: DP ¼ Pin Pout ¼ s/R.
Fig. 10 e Sketch of three-dimensional of a droplet sitting on a surface in a small gas channel exposed to a gas flow moving right to left at Ca ¼ 0.083. The initial contact angle q at equilibrium is 165.38 . Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Fig. 11 e Streamlines behind the droplet in XZ plane at Ca ¼ 0.083.
Where DP is the pressure difference across the fluid interface s is the surface tension, R is the radii of curvature, Pin and Pout are the pressures inside and outside of the drop, respectively. The simulation was conducted for several initial droplet radii. Then the changes of pressure differences with respect of 1/R are plotted. Drops with different radii are generated inside a computational domain using the ShaneChan model. A plot of P versus 1/R is a straight line and its slope will be s as it is seen in Fig. 6. For Dr ¼ 438.686 mu lu2, G ¼ 120, J0 ¼ 4 and r0 ¼ 200 mu lu2, DP ¼ 14.473 (1/R) þ 0.00181, the surface tension is equal to 14.473 lu mu ts2. We can conclude that Laplace's Law is well satisfied for this model. The equilibrium contact angle, q defined as the angle between the vaporeliquid interface and the solid wall, can be obtained through numerical experiment of droplet placement on solid surfaces with different parameter Gads. Fig. 7 shows droplet shapes for various fluid-solid interaction strengths. Three typed of solid surfaces were defined with Gads ¼ 50, 187.16 and 300 corresponding to q ¼ 177 , 90 and 17 0.8. In Fig. 8 we report contact angle versus fluid-solid interaction parameter Gads.
7
It is found that contact angle varies al-most linearly with the fluid-solid interaction strength Gads. A linear relation between Gads and q can be given by: q ¼ 0.64Gads þ 209.8. We implemented the lattice Boltzmann model for nonideal fluids to simulate the droplet dynamics in a microchannel of a proton exchange membrane fuel. Initially, the three-dimensional computation domain consists of an incompressible fluid droplet in a rectangular micro-channel (See Fig. 9). The two steps “stream and collide” algorithm for a (D3Q19) lattice is used to simulate lattice Boltzmann equation on 150 40 30 site lattices. Bounce back boundary conditions are imposed on the walls and periodic boundary conditions are imposed on the two domain ends. To mimic the flow inside the micro-channel, a constant body force along the z direction was applied to displace both the fluid and the droplet. Driven by this force the droplet move but its speed and shape change with the inlet velocity and the wettability. There are several dimensionless numbers encountered in micro-fluidic applications. For low Reynolds number regime (fluid inertia forces are negligible), the displacement of immiscible fluids in micro-fluidic devices can be described by two main dimensionless numbers: the Capillary number (Ca) and Bond number (Bo). The capillary number reflects the competition between viscous forces and capillary forces. It is given by: Ca ¼ mU/s. Where s is the surface tension which is calculated using Laplace's law test, U is a characteristic velocity and m is the dynamic viscosity of the liquid phase. The Bond number reflects the balance between gravitational and capillary forces and is defined as: Bo ¼ R2Drg/s. where g is the gravitational acceleration. In this work, the focus will be on the influence of capillary number. Calculation of Ca requires the viscosity which is given by:y ¼ c2s ðt 1=2Þ where cs is the sound velocity and t is the BGK relaxation time. We focus on the dynamic characteristics of liquid drop's spreading. We will present essentially the velocity and density fields and streamlines. As a first simulation step, we have simulated droplet motion through a 150 40 30 voxels micro-channel. Normalized density plots at various time steps were generated. In order to focus on the effect of shear the shear induced by the gas on the droplet, the physical parameters are chosen so that the gravity is negligible compared with the surface tension This means that the Bond number is less than 1. The droplet is initially spherical caps centered on the pore due to surface tension, after that the droplet moves down the micro-channel in the direction of the gas flow. The droplet will deform as a
Fig. 12 e Streamlines generated from 2D mean velocity vectors in YZ (a) and XY (b) planes (Ca ¼ 0.083). Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
8
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
Fig. 13 e Evolution of the interfaces of droplets for different capillary numbers.
result of these forces and exhibit a contact angle hysteresis defined as the difference between the advancing and receding contact angles. Fig. 10 shows the evolution of the phase distribution patterns under external force (The fluid flow is always laminar Re z 1.9). After a certain time step, the droplet continues to move at constant speed and reaches equilibrium with an unchanged shape. The streamlines are plotted to show the flow field in the micro-channel. Flow recirculations behind the droplet and near the contact point are clearly visible in Fig. 11. Velocity streamlines in XY and YZ planes are presented in Fig. 12. Flow recirculation around the bubble is observed. Fig. 13 illustrates the effect of capillary number on the deformation of droplets with q ¼ 165.38 , this comparison shows that as the capillary number increases, the droplet displacement and deformation also increase. Fig. 14 shows the high deformation of droplet under high capillary number. The droplet is strongly deformed at Ca ¼ 0.2075 but the contact surface area decreases. We note that for low capillary number the droplet cannot be evacuated from the micro-channel. The streamlines plot show a symmetrical vortex around the deformed droplet (See Fig. 15).
Fig. 15 e Vortex around the deformed droplet at Ca ¼ 0.2075.
We consider now the wetting case, Fig. 16 shows the droplet behavior in micro-channel with hydrophilic surface (q < 90 ). Because of shearing forces, the droplet breaks down into two droplets with different sizes: a principal droplet and a satellite droplet which touches the top wall, this droplet is drained and disappears completely. The wettability effect is very strong and can transform droplets in thin liquid films spreading along the solid walls. We can remark that the droplet velocity is relatively high for a hydrophobic channel walls than for a hydrophilic walls.
Conclusion We have developed an LBM algorithm with an external force for two-phase droplet dynamics in three-dimensional microchannel of a proton exchange membrane fuel, in this model the two phase fluid is modeled by ShaneChen scheme. The density is simulated inside the micro-channel. We remark that the droplet is highly deformed into a shape with a tip. It is observed that the droplet deformation increases with increase in capillary number.
Fig. 14 e Droplet behavior under high capillary number Ca ¼ 0.2075. Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
9
Fig. 16 e Time evolution of dynamic behavior of droplet in micro-channel with hydrophilic surface. (q < 90 ).
It is clear that the average contact angle hysteresis, which is the difference between the advancing and receding contact angles, increases during the droplet motion. In the wetting case, droplet adheres more to the channel walls and its evacuation became more difficult than the case of a hydrophobic micro-channel. Recirculation in the wake behind the droplet is observed which depend on the capillary number Ca and the wettability. Other geometries must be tested and porous effect must be included in future work.
references
[1] Anderson R, Zhang L, Ding Y, Blanco M, Bi X, Wilkinson DP. A critical review of two-phase flow in gas flow channels of
[2] [3]
[4]
[5]
[6]
[7]
proton exchange membrane fuel cells. J Power Sources 2010;195:4531e53. Barbir F. PEM fuel cells: theory and practice. Elsevier Academic Press; 2005. Bhatnagar PL, Gross EP, Krook M. A model for collision processes in gases, I. Small amplitude processes in charged and neutral one-component system. Phys Rev 1954;94:511e25. Buie CR, Posner JD, Fabian T, Cha SW, Kim D, Prinz FB, et al. Water management in proton exchange membrane fuel cells using integrated electroosmotic pumping. J Power Sources 2006;161(1):191e202. Chen H, Chen S, Matthaeus WH. Recovery of the NavierStokes equations using a lattice-gas Boltzmann method. Phys Rev A 1992;45:5339e42. Chen J. Experimental study on the two phase flow behavior in PEM fuel cell parallel channels with porous media inserts. J Power Sources 2010;195:1122e9. Chen S, Doolen GD. Lattice Boltzmann method for fluid flows. Annu Rev Fluid Mech 1998;30:329e64.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077
10
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y x x x ( 2 0 1 4 ) 1 e1 0
[8] Cho SC, Wang Y, Chen KS. Droplet dynamics in a polymer electrolyte fuel cell gas flow channel: forces, deformation, and detachment. I: theoretical and numerical analyses. J Power Sources 2012;206:119e28. [9] D'Humieres D, Lallemand P. Numerical simulations of hydrodynamics with lattice gas automata in two dimensions. Complex Syst 1987;1:598. [10] Eikerling M. Water management in cathode catalyst layers of PEM fuel cells. J Electrochem Soc 2006;153(3):E58. [11] Frisch U, Hasslacher B, Pomeau Y. Lattice-gas automata for Navier-Stokes equations. Phys Rev Lett 1986;56:1505e8. [12] He X, Luo LS. Theory of the lattice Boltzmann method: from the Boltzmann equation to the lattice Boltzmann equation. Phys Rev E. 1997;56:6811e8. [13] Hussaini IS, Wang CY. Visualization and quantification of cathode channel flooding in PEM fuel cells. J Power Sources 2009;187:444e51. [14] Li H, Tang Y, Wang Z, Shi Z, Wu S, Song D, et al. A review of water flooding issues in the proton exchange membrane fuel cell. J Power Sources 2008;178:103e17. [15] Liu HC, Yan WM, Soong CY, Chenc F, Chud HS. Reactant gas transport and cell performance of proton exchange
[16]
[17] [18]
[19]
[20]
[21]
membrane fuel cells with tapered flow field design. J Power Sources 2006;158:78e87. Martys NS, Chen H. Simulation of multicomponent fluids in complex three-dimensional geometries by the lattice Boltzmann method. Phys Rev E. 1996;53:743e50. Mench M. Fuel cell engines. John Wiley & Sons, Inc; 2008. Metz T, Paust N, Mu¨ller C, Zengerle R, Koltay P. Passive water removal in fuel cells by capillary droplet actuation. Sens Actuators A e Phys July 2007;5 [in press], Corrected Proof. Minor G, Zhu X, Oshkai P, Sui P, Djilali N. Water transport dynamics in fuel cell micro-channels. In: Kakac S, Pramuanjaroenkij A, Vasiliev L, editors. Fundementals and applications. Dordrecht: Springer; 2008. Mini-Micro Fuel 75 Cells. Yuan XZ, Wang H, Zhang J, Wilkinson DP. Bipolar plates for PEM fuel cells-from materials to processing. J New Mater Electrochem Syst 2005;8(4):257e67. Zhu X, Liao Q, Sui PC, Djilali N. Numerical investigation of water droplet dynamics in a low-temperature fuel cell microchannel: effect of channel geometry. J Power Sources 2010;195:801e12.
Please cite this article in press as: Ben Amara MEA, Ben Nasrallah S, Numerical simulation of droplet dynamics in a proton exchange membrane (PEMFC) fuel cell micro-channel, International Journal of Hydrogen Energy (2014), http://dx.doi.org/ 10.1016/j.ijhydene.2014.09.077