Baleta, Jakov

Duić, Neven

Pachler, Klaus

Vujanović, Milan

FAMENA Repository

11th World Congress on Computational Mechanics (WCCM XI) 5th European Conference on Computational Mechanics (ECCM V) 6th European Conference on Computational Fluid Dynamics (ECFD VI)  (Eds)    IMPLEMENTATION AND VALIDATION OF THE MATHEMATICAL MODEL OF SURFACE TENSION INTO CFD WALL FILM MODULE JAKOV BALETA*, MILAN VUJANOVIĆ*, KLAUS PACHLER†, NEVEN DUIĆ* * University of Zagreb,  Ivana Lučića 5, 10002 Zagreb, Croatia,  e-mail: jakov.baleta@fsb.hr  †Advanced Simulation Technologies,  AVL List GmbH,  Alte-Poststraße 152, A-8020 Graz, Austria,  e-mail: klaus.pachler@avl.com  Key Words: wall film, computational fluid dynamics, surface tension, droplet spreading, capillary force  Summary. Horizontal liquid film flow sheared by an external air flow field is encountered in many engineering applications. Taking into account its physical complexity and high cost of experimental investigation, numerical simulations are nowadays considered as a valuable alternative. They are becoming a useful tool for detailed understanding of complex flow characteristics and transport phenomena, especially in situations where experimental measurements are infeasible or too expensive. The focus of this paper is the implementation and validation of a mathematical model of surface tension effects within the existing numerical framework in order to achieve a more accurate description of the liquid wall film phenomena. After literature review, optimum mathematical model has been chosen and implemented in commercial computational fluid dynamics (CFD) code. Validation was carried out using a well established case of isothermal droplet spreading for which there are analytical expressions. Comparison of simulation results with non-dimensional droplet profile shows excellent agreement with analytical results and gives confidence for commercial application of implemented model.  1 INTRODUCTION Liquid film flow sheared by an external air flow field is a physical phenomenon encountered in many engineering applications such as: burners, rocket nozzles, heat exchangers, steam turbine blades and especially internal combustion (IC) engines due to present and upcoming stringent environmental regulations [1].  In IC engines it has been observed that unburnt fuel that goes directly to the manifold causes an increase in the emissions of unburnt hydrocarbons in the petrol engines and the larger product of soot in the compression-ignited engine. Also, injection of precursor substance into exhaust gases before the catalyst leads to the formation of liquid wall film due Jakov Baleta, Milan Vujanović, Klaus Pachler, Neven Duić   2 to unsteady engine working conditions [2]. Above mentioned examples show great importance of the correct prediction of wall film behaviour. Assumptions of thin liquid film lead to the implementation of the wall film model as a 2D finite volume method on the air flow wall boundaries. Current wall film model in commercial CFD code Fire contains mathematical description of physical phenomena in the form of three conservation laws: conservation of mass, momentum and energy. However, surface tension effects which are of great importance at the late film spreading stages when the inertial forces are negligible [3] and also in cases when shear-driven film comes at a sharp expanding corner [4] are not taken into account. The film surface tension effects that are present at the boundary edges of the liquid phase play a crucial role in slowing down the film progression. These effects are not currently accounted and the approximation of constant film height over the control volume does not allow reconstructing the characteristics of the interface. That way scattering of liquid film causes high numerical diffusion which does not correspond to physical reality. The threshold value can be set by the user, and only the cells where the film is higher than this value are considered active. This approach is not acceptable and therefore the object of this work is improvement of current model by incorporating proper mathematical description of surface tension effects.  The continuum surface force (CSF) method of Brackbill has been employed extensively over the last 13 years to model surface tension in various fixed (Eulerian) mesh formulations for interfacial flows, in particular in the volume-of-fluid (VOF), level-set (LS) and front tracking (FT) interface representation techniques [5]. The main drawback of mentioned approaches is that there are unsuitable for finite area Eulerian approach because of computational demands. They are usually used on small scale interfaces such as impinging droplet or few droplet collisions. Another approach of describing surface tension was made by Bai [6] who developed the solution procedure where capillary forces are taken into account through a capillary pressure term. It is suitable to incorporate this pressure term into momentum equation of wall film and details of the procedure are given in the following section.  2 MATHEMATICAL MODEL The equations of continuum mechanics are based on the conservation laws for mass, momentum and energy. The general form of the time averaged conservation equation for any dependent variable φ of the continuous phase in the differential form is:    (1) where ρ is the density, uj Cartesian velocity, Γφ diffusion coefficient, and Sφ is the source term of the dependent variable φ. The source term Sφ is used for the coupling of the liquid and the gaseous phases [7]. Fundamental assumptions and simplifications of the wall film model incorporated in commercial CFD code Fire are listed below:  - Gas and wall film flow are treated as separate single phases; the coupling of the two phases is achieved by a modified set of boundary conditions based on a semi-( ) ( ) ( ) ,jj j ju St x x x           Jakov Baleta, Milan Vujanović, Klaus Pachler, Neven Duić   3 empirical relations; - The film thickness is very small in relation to the mean diameter of the gas flow, so no adaptation of the volume grid to the film surface is necessary; - Due to the thin film and its small velocity, wall friction and interfacial shear stress dominate the film behaviour – a momentum equation could be dropped for a steady state; - Wall temperature is below the Leidenfrost point; - The wavy surface of the film is modelled as a mean film thickness with a superimposed film roughness; - Mean film surface is assumed to be parallel to the solid wall. The above assumptions lead to the implementation of the wall film model as a 2D finite volume method on the air flow wall boundaries.  The film thickness equation is the basic governing equation for the wall film flow. It represents a slightly modified formulation of the continuity equation where, instead of mass, the wall film thickness is conserved property. The Cartesian formulation of the film thickness equation is: 1 21 21( )mD mVu uS St x x A        (2) where δ is the film thickness, ρ is the film density, u1 and u2 are film velocity components, SmD and SmV are source terms and A is the surface of the film. If we assume that the source terms are provided, equation (2) can be solved explicitly if the velocity components are known. Film momentum equation describes dynamics of liquid film interaction with its environment - wall, air stream above film, impinging droplets, etc. Equation (3) gives mathematical formulation of wall film momentum conservation law:  MiiiLiLjiii SmgdLnpdLnVuudtdM  ˆˆ)(   (3) where Mi is film momentum, ρ is the film density, ui is film velocity, Vj is wall velocity, inˆ  is normal to the face cell facing outwards, L is length of the face cell boundary, δ is the film thickness, p is film pressure, m is film mass, gi is gravity vector, Γi is the term that takes into account all shear stresses and Sm present various source and sinks terms such as film entrainment, spray droplets impingement and film evaporation.   Inclusion of surface tension effects into pressure term of momentum equation was made by calculating capillary pressure. This force drives the surface towards a minimal energy state characterized by a configuration of minimum surface area and is represented by the combined action of liquid surface tension σ and film surface curvature C [3]: Cp    (4) As shown in [8, 9], the mean curvature of the free-surface can be approximated using the following expression: Jakov Baleta, Milan Vujanović, Klaus Pachler, Neven Duić   4 22212111RRRRC  (5) On quality computational meshes surface area patches are so small than first two terms of equation (5) could be neglected. So, film surface curvature could be approximated with Laplacian of film thickness.   Figure 1: Schematic representation of wall film cells Using schematic representation of wall film cells given on the Figure 1, Laplacian of film thickness could be numerically approximated with following equation: edgesedges niiABiniii LABnSLnS 1111  , (6) where S is surface patch area, Li is length of the neighbour edge, ni is unit normal on common edge facing outwards and AB  is length between cell centers A and B. Equation (6) was incorporated into pressure term of momentum equation of existing numerical framework.  3 ANALYTICAL CASE OF DROPLET SPREADING In order to test the numerical scheme and to check if the capillary pressure effect can be predicted properly, there is need for a simple case, where only capillary force effects are relevant and possibly for which there is analytical solution. Thus, the spreading of an isothermal droplet on a solid surface is simulated. The droplet is driven with Laplace (capillary) pressure as dominant force and other forces being absent or negligible. Similar tests were obtained by Diez et al. [10] (proposed analytical solution using lubrication theory which neglects convective terms in momentum equation and obtained experiments). The full momentum equations are solved together with the interface mass continuity equations, despite the fact that the convective term in the momentum equations may be negligible. The problem can be observed as two-dimensional because of axial symmetry. Diez et al. [10] showed that the normalized film thickness h/h0 can be expressed as a single function of the scaled radial position r/rf, irrespective of the time level.  Jakov Baleta, Milan Vujanović, Klaus Pachler, Neven Duić   5 4 NUMERICAL SIMULATIONS The simulation domain with relevant boundary conditions is shown in Figure 3. A three dimensional computational mesh with dimensions of 5 x 5 x 0.1 mm with 20 000 orthogonal hexahedron cells was used for the simulation. The wall boundary condition was defined at the bottom of the domain, whilst a static pressure outlet was imposed at all other sides. The domain pressure was 1 bar and temperature was equal to the droplet temperature. Pressure velocity coupling of momentum and continuity equation was obtained using the SIMPLE/PISO algorithm. The central difference discretisation scheme was used for the convective term in the continuity equation with a blending factor of 1, whilst a MINMOD Relaxed with a blending factor of 0.5 was used for the convective terms in the momentum equations. Turbulence model was deactivated since quiescent air is necessary for comparison with analytical expression. Energy equation was also deactivated in order that the only relevant force stays capillary pressure which tends to spread given droplet shape.   Figure 2: Computational domain with boundary conditions and initial shape of droplet  The simulation starts from a drop which initially takes the shape of a rotational paraboloid, shown on the right part of Figure 2, with volume of 0.12 mm3 and thickness at the center of 0.08 mm. The spatial and time discretization increment were 0.05 mm and 10-5 s. Simulation time was 0.5 s since spreading of the droplet practically remained unchanged around and after 0.4 s.  5 RESULTS Qualitative depiction of droplet evolution is presented on Figure 3. As can be seen, droplet shape changes only after 1 ms. Reason behind this phenomena is that only relevant force is capillary pressure which produces fairly small spreading velocities whose effects on droplet spreading is visible only after 1ms. Further time evolution shows that droplet shape doesn’t change much between 0.25 and 0.4 s, since the wall film thickness gradients on the edge of the droplet tend to zero. On the same picture also could be noticed that the behaviour of the droplet during the whole period of the simulation remains symmetrical which is qualitatively in agreement with analytical results of Diez [10].  Jakov Baleta, Milan Vujanović, Klaus Pachler, Neven Duić   6  Figure 3: Simulation of droplet spreading  Quantitative comparison given on Figure 4 shows that the present numerical predictions are in excellent agreement with the experimental data, indicating that the numerical modelling of the capillary pressure is reasonable in terms of optimum between accuracy of results and computational demands. However, the model does not include correction of the film curvature at the film front due to the wettability effects represented by (dynamic) contact angle. These effects are omitted from modelling due to the relatively large spatial discretisation which prevails in most of the practical applications where wall film phenomena are important.    Figure 4: Comparison of simulated results with theoretical droplet profile Jakov Baleta, Milan Vujanović, Klaus Pachler, Neven Duić   7 6 CONCLUSIONS After literature review, appropriate surface tension model was chosen taking into account an optimum between computational demands and accuracy of the solution. Thus, effects resulting from capillary forces are modelled as a part of a pressure term of the wall film momentum equation. Simulation results conducted with the new model were compared to analytical non-dimensional profile of isothermal droplet spread on flat surface. Agreement is excellent and gives confidence for commercial application of implemented model.     In this paper was shown that incorporation of surface tension effects into existing wall film model can increase its capabilities in prediction of film propagation with a high accuracy, providing this way a powerful tool for the design and optimization of devices where liquid film phenomena play an important role. The continuation of this work would entail implementation of proper wall film rupturing model in which surface tension has important role in film stabilization opposing inertial forces.  REFERENCES [1] European Union, Regulation (EC) No 715/2007, 2007. [2] F. Birkhold, U. Meingast, P. Wassermann, O. Deutschmann, Modeling and simulation of the injection of urea-water-solution for automotive SCR DeNOx-systems. Appl. Catal. B: Environ., No. 70, pp. 119-127, 2007.  [3] Kristijan Horvat. Computational Modelling of Spray Impingement Accounting for the Wall Film Formation, Forschungsberichte Stromungslehre Und Aerodynamik, Shaker Verlag, Gemany, 2007. [4] Friedrich, Lan, Wegener, Drallmeier, Armaly: A Separation Criterion With Experimental Validation for Shear-Driven Films in Separated Flows, J. Fluids Eng., Volume 130, 051301, 2008. [5] M. M. Francois, S. J. Cummins, E. D. Dendy, D. B. Kothe, J. M. Sicilian, M. W. Williams, A balanced-force algorithm for continuous and sharp interfacial surface tension models within a volume tracking framework, Journal of Computational Physics 213, pp. 141–173 , 2006. [6] C. Bai, Modelling of spray impingement processes, Ph.D. thesis, Imperial College, University of London (1996). [7] H. Mikulčić, E. von Berg, M. Vujanović, P. Priesching, R. Tatschl, N. Duić, Numerical analysis of cement calciner fuel efficiency and pollutant emissions. Clean. Techn. Environ. Policy. 15, pp. 489–499, 2013. [8] S. B. G. O’Brien, L. W. Schwartz, Thin film flows: theory and modeling, Encyclopedia of surface and colloid science: second edition (2006) 6304–6317. [9] L. W. Schwartz, D. E. Weidner, Modeling of coating flows on curved surfaces, Journal of engineering mathematics 29 (1995) 91–103 [10] Diez, J.A., Gratton, R., Thomas, L.P., Marino, B., (1994), Laplace pressure driven drop spreading, Phys. Fluids, 6, 24-33  

Implementation and validation of the mathematical model of surface tension into CFD wall film module

Abstract

Similar works

Full text

Available Versions

FAMENA Repository