CHEMICAL ENGINEERING TRANSACTIONS VOL. 76, 2019 A publication of The Italian Association of Chemical Engineering Online at www.aidic.it/cet Guest Editors: Petar S. Varbanov, Timothy G. Walmsley, Jiří J. Klemeš, Panos Seferlis Copyright © 2019, AIDIC Servizi S.r.l. ISBN 978-88-95608-73-0; ISSN 2283-9216 Modelling the Motion of a Single Solid Bead in a Newtonian Fluid by Two-Phase CFD Methods Lívia Gyurika,*, Attila Egedyb, Zsolt Ulberta, Kevin Croninb, Denis Ringb aUniversity of Pannonia, Department of Process Engineering, Egyetem St. 10, 8200 Veszprém, Hungary bUniversity College Cork, Ireland (Republic) gyurikl@fmt.uni-pannon.hu Sedimentation is important in many processes for example in cement, chemical and food industries. Investigation of the sedimentation process of solid bodies showed a phenomenon as settling velocity is not the constant value as it is expected from the analytical deduction; rather it has fluctuating behaviour. Seeking an explanation, detailed two-phase hydrodynamic models are developed. In the present work, two-dimensional incompressible flow equation in conservative Euler form is solved numerically on PC by an implicit finite difference method, and the fluid-solid interaction is described by the Immersed Boundary Method. For additional model validation, Particle Tracing module of COMSOL is applied as well. The self-developed simulator based on an a priori (white box) model, and the 3D model by the commercial CFD software is used for modelling the two-phase system and the particle motion in the fluid column. Settling velocity is calculated and compared to the velocity of a particle recorded by a high-speed camera measurement system in the laboratory-scale sedimentation equipment. The results are promising, and after improvements, the simulator can be a useful tool in supporting energetically optimal industrial design and operation. 1. Introduction Sedimentation is the name of the phenomena that particles of larger density than the density of the fluid settle down. The settling process is used in many fields of industry and science such as in geology (Van Santen et al., 2007), wastewater treatment (Ochowiak et al., 2017) and various areas of chemical engineering. Sedimentation can be a method of separating, and thus a part of producing processes. A separator device usually makes enhanced sedimentation by multiplying the gravity force by spinning. Solid particles with greater density than the density of the fluid in still liquid also show settling behaviour due to natural gravity. This motion is caused by the sum of forces acting on the system. Gravity acts in a positive direction, while buoyancy and drag forces act as counterforces. Describing the falling sphere is often divided into two regimes according to Reynolds number. Stokes law is valid where Re<<1, and Newton’s drag works for larger Re number flows. There are multiple equations in the literature for calculating the mean velocity, i.e. hindered settling (Kondrat’ev and Naumova, 2007), Richardson-Zaki velocity (Richardson and Zaki, 1997), Batchelor’s formula (Batchelor, 1972), however, microscale velocity is determined by hydrodynamics. Beside numerous theoretical and empirical studies, computer-aided estimations are getting more and more often used since numerical methods, and computers which can evaluate those calculations are developed. Ready to use software modules are comfortable, but in most cases, it is hard to modify. Create an own simulator for the given problem can be more adjustable and gives more control over the implementation and solving methods. Building up an a priori (white box) model besides finding and implementing the best method for the numerical solution is a challenging and interesting task for a researcher. In the following paragraphs the most common two-phase modelling approaches are introduced briefly. Two-fluid method (TFM) or Euler-Euler method treats both phases as a fluid even when one of them consists of solid particles. Two sets of momentum equations should be calculated with the addition of a volume fraction term and an interaction term to both sets of equations. Earlier (Noetinger, 1989) and recent (Patel et al., 2017) works are published by using TFM for modelling the sedimentation of a solid-liquid suspension. The method is capable of simulating a gasifier too (Tran et al., 2018), where the two phases are solid and gas. Use of TFM DOI: 10.3303/CET1976030 Paper Received: 15/03/2019; Revised: 15/04/2019; Accepted: 29/04/2019 Please cite this article as: Gyurik L., Egedy A., Ulbert Z., Cronin K., Ring D., 2019, Modelling the Motion of a Single Solid Bead in a Newtonian Fluid by Two-Phase CFD Methods, Chemical Engineering Transactions, 76, 175-180 DOI:10.3303/CET1976030 175 results in a less detailed solution, but calculations can be performed with lesser computation cost and can give us overall insight about the operation of the devices. CFD software offer opportunities to handle immiscible phases too, however calculating flow properties of solid-fluid two-phase flow, these tools are usually limited. The most common disperse methods root in the Euler-Euler method. TFM solutions could be applied to simulate the movement of a crowd of particles in a fluid. Thus its advantage the neglect of the shape of the particles, so for example eggshell particles can be modelled as well as perfect spheres (Jurena et al., 2017). Discrete element method (DEM) is a modelling approach which results in detailed information about each element of the solid phase. It computes the particle-particle and particle-wall collision, and the movement is calculated by the sum of contact forces emerge on the particle. (Zhu et al., 2007) This approach gives detailed information about the solid phase, but it needs large computational capacity. Simulation of two-phase flow and calculation of the motion of solid particles by DEM can also be performed by CFD software packages using user-defined functions. CFD-DEM coupling was carried out to simulate particle flow in many research fields such as crystallisation (Kerst et al., 2017) and fluidisation (Zhang et al., 2015). Both TFM and DEM methods use approximately one order of magnitude larger mesh than the size of a typical particle of the solid phase. Direct numerical simulation (DNS) techniques are another branch of two-phase flow modelling approaches. These use one order of magnitude smaller mesh to calculate the flow field and the replacement of the solid particles. Two main subtypes of this method are the one which uses body-fitted mesh, and the other which use regular Cartesian time-independent mesh. By use of body-fitted mesh, the particles occur as a boundary or wall, while by the use of one of the methods of the other subtype, a so-called immersed boundary takes place where the particle occurs only virtually without a physical boundary. This makes the method computationally acceptable, while we have to face other kinds of challenges by using immersed boundary-type methods, for example how to define the connection between the node points of the fluid and the solid (Ghosh and Stockie, 2015). Immersed Boundary Method (IBM) belongs to the second subtype of DNS which uses the regular steady Cartesian grid. The method was first introduced for describing blood flow in the heart where valves acted as a moving boundary (Peskin, 1972), structured calculation mesh can be applied, and the inner boundaries (e.g., particles) can be virtually defined using body forces and extra force term added to the momentum equation. Structured (or uniform, Cartesian) mesh is often called Eulerian grid while the node points of the immersed boundary referred to as Lagrangian points. These two coordinate systems are independent, but an interpolating function gives the relation between locations and flow variables. All these methods are called Computational Fluid Dynamics (CFD) techniques which are quickly spreading tools to calculate and simulate flow fields and flow related problems. Nowadays there are numerous available commercial and open-source CFD software, which offer approaches to handle two-phase flows. Level Set and Phase Field methods are applicable mainly for the fluid-fluid system to track the interface between the two immiscible fluids. Moving mesh is a body-fitted DNS technique, hence its computational cost is high. Particle Tracing is a module of COMSOL Multiphysics, which calculates the trajectory of the particles with given properties according to the previously calculated momentum balance (flow-based particle tracing). The latter method seems proper to use to model the motion of a single solid bead in a Newtonian fluid. By modelling a solid-liquid settling system, a constant terminal velocity is expected according to rough physics, however, in real systems perturbation around this settling velocity is observed. Measurements show that the settling velocity of the falling body fluctuates in both direction and magnitude. Our system consists of 18 °C water in a 0.1 m diameter, 1 m high cylinder, and a 5.9 mm diameter nylon sphere. The measured terminal velocity of the bead (143 mm/s) is in agreement with the analytical solution for high Reynolds number flows; the prevailing particle Reynolds number was around 840. In the present study use of a 2D model is chosen as a good and computationally feasible approximation of the real 3D system of a settling cylinder with a single spherical bead. The model is calculated by a numerical solution of the incompressible, isotherm, inviscid flow equations (continuity equation and momentum equations) by a central finite difference method. The proper pressure field and the correct velocity field is calculated by using a pressure correction method based on Patankar’s SIMPLE algorithm (Patankar, 1980), and fluid-structure interaction is handled by the Immersed Boundary method. Another modelling approach is studied as well, where the geometry is defined in 3D, and the model is solved by the Particle Tracing module of COMSOL Multiphysics. Calculations are performed on a Dell Optiplex 790 PC with 16 GB memory. The results are compared to measurement data to validate the model. Once we have a proper model, it can be used to different other materials, moreover other solid-liquid problems according to interest. 2. Materials and methods 2.1 Experimental setup and video processing A laboratory-scale (diameter 0.1 m and height 1 m) settling cylinder is given for experimental validation (Figure 1a). The tube has a lid with a particle-size hole in the centre to help by the possible most accurate 176 releasing. The column is filled with water, and the whole equipment is placed in an approximately 18 °C (291.15 K) room. The spherical particles are from 1,114 kg/m3 density nylon plastic, 5.9 mm in diameter and weight of 0.12 g. The bead is produced by plastic injection technique, which may cause tiny spikes by the junction of the two hemispheres. Figure 1: a) The laboratory-scale cylinder for sedimentation, b) Trajectory of the bead in water (horizontal scale magnified compared to vertical scale) Position is recorded by the camera in every 33.33 ms The movement of the particle was capture, and we gained the position data in two dimensions (Figure 1b) by data processing. The oscillating behaviour is well observable, although note that the scale strongly differs. By repeated experiments the bead runs different trajectory, but the oscillating vertical velocity is observable in all measured cases. The video is recorded in a 40 cm section of the cylinder, where the vertical velocity has already achieved the settling velocity, and the bottom effect is still negligible. 2.2 Modelling approach Since liquids are incompressible and we suppose in the current case that the system is isotherm, the governing equations consist of the incompressible continuity equation and momentum equations in x and z directions. We use them in conservative form for convenience. As the flow is driven by a pressure gradient, proper calculation of the pressure field is necessary, and thus we can get the right velocity field which fulfils the continuity equation. One widely used method for this task is the iterative SIMPLE algorithm (Patankar, 1980). Particle Tracing of COMSOL Multiphysics is a quasi-two-phase modelling approach. A single spherical particle with given properties (density and diameter) can be introduced to the liquid column at the top boundary, and the bottom can be defined as the outlet. Gravity and drag forces are calculated according to physics, and an additional adjustable force can be added to balance the forces and gain a terminal velocity. Flow equations In this study, Eq(1) - Eq(2b) are discretised by a semi-implicit method called SIMPLE for gaining a numerical solution of the flow field in two dimensions. 𝜕𝑢 𝜕𝑥 + 𝜕𝑣 𝜕𝑧 = 0 (1) 𝜌 𝜕𝑢 𝜕𝑡 + 𝜌𝑢 𝜕𝑢 𝜕𝑥 − 𝜌𝑣 𝜕𝑢 𝜕𝑧 = − 𝜕𝑝 𝜕𝑥 (2a) 𝜌 𝜕𝑣 𝜕𝑡 + 𝜌𝑢 𝜕𝑣 𝜕𝑥 + 𝜌𝑣 𝜕𝑣 𝜕𝑧 = − 𝜕𝑝 𝜕𝑧 (2b) where u and v are velocity in x and z-direction, x and z are space directions, ρ is density, t is time, p is pressure. These unsteady partial differential equations exhibit a mixed elliptic-parabolic behaviour. At the vertical walls, a no-slip condition is applied. Inflow velocity is 0 m/s, as only the motion of the particle forces the liquid to move. 177 Fluid-solid interaction Coupling between liquid-solid two-phase flows is modelled by the Immersed Boundary Method. This method requires the introduction of an extra force into the momentum equation to modify the fluid flow by the solid body surface Eq(3a) and Eq(3b). 𝜌 𝜕𝑢 𝜕𝑡 + 𝜌𝑢 𝜕𝑢 𝜕𝑥 − 𝜌𝑣 𝜕𝑢 𝜕𝑧 + 𝜕𝑝 𝜕𝑥 + 𝐴 ∙ 𝐹𝑥 = 0 (3a) 𝜌 𝜕𝑣 𝜕𝑡 + 𝜌𝑢 𝜕𝑣 𝜕𝑥 + 𝜌𝑣 𝜕𝑣 𝜕𝑧 + 𝜕𝑝 𝜕𝑧 + 𝐴 ∙ 𝐹𝑧 = 0 (3b) where A is the stencil to switch the effect of the force, Fx and Fz are the extra forces in x and z-direction respectively. The applied forcing scheme that is used to calculate the body force follows the direct forcing method as described in (Fadlun et al., 2000). The main idea of this method is to set the value and direction of F to balance momentum at those grid points where the solid immersed object takes place (Eq(4a) and Eq(4b)). The A stencil matrix has the size of the Eulerian node points and contains ones and zeros according to the presence or absence of the particle. In the present case, stepwise interpretation of the particle is introduced to A. However it will be smoothened by using a distance-based weight function to follow the circular shape more accurately. 𝐹𝑥 = 𝑢 𝜕𝑢 𝜕𝑥 + 𝑣 𝜕𝑢 𝜕𝑧 + 1 𝜌 𝜕𝑝 𝜕𝑥 − 𝑢0 𝜕𝑡 (4a) 𝐹𝑧 = 𝑢 𝜕𝑣 𝜕𝑥 + 𝑣 𝜕𝑣 𝜕𝑧 + 1 𝜌 𝜕𝑝 𝜕𝑧 − 𝑣0 𝜕𝑡 (4b) where u0 and v0 are the velocity components in the previous time step. The three-dimensional COMSOL model by Particle Tracing calculates by the fluid-solid interaction (included viscosity) too, however in a hidden way. It uses both user-defined and automatically calculated forces to act on the particle. An extra force is set according to reach the measured velocity and is introduced as a manually fitted parameter by identification. 3. Results and discussion According to the measurement of the mean velocity of the falling bead, the Reynolds number of the flow is around 844, see Eq(6). Stokes’ drag does not fit to the case, rather Newton’s drag. 𝑅𝑒 = 𝜌𝑣𝐿 𝜇 (6) where ρ is the density of the fluid (998 kg/m3), v is the measured settling velocity (143 mm/s), L is the characteristic length (diameter of the bead, 5.9 mm), μ is the dynamic viscosity of the fluid (10-3 Pa s). Newton’s second law of motion with Newton’s drag force Eq(7) gives analytically a similar velocity value for settling. 𝐹𝐷 = 1 2 𝜌𝑓𝑙𝑢𝑖𝑑 𝑣 2𝑐𝐷 𝑑2 4 𝜋 (7) where cD is the drag coefficient, which is close to 0.44 for spheres (at high Reynolds numbers). Settling velocity calculated analytically is 141 mm/s which is in good agreement with the measured 143 mm/s velocity value. The model using the Immersed Boundary Method is calculated in a reduced 2D computational domain (in a narrower and shorter cylinder), but the main characteristics are gained. In Figure 2a and 2b the horizontal (x) component of the fluid velocity is shown in the first and the last computed time-step. While the sphere falls, the fluid goes away and fills the place again, where the particle has already gone. Direct forcing does not let the fluid flow into the solid phase, so fluid velocity is not calculated inside the virtual boundary, however the bead’s velocity is depicted (Figure 2a - 2d), because that is defended by the direct forcing. The vertical velocity pattern shows increasing values by the particle according to our expectation because of the narrowed diameter caused by the bead (Figure 2d). The other part of the coupling is the fluid phase’s effect on the solid phase. The drag force is calculated from the pressure field of the fluid around the circle. Running the code for 20,000 time-steps with dt = 2.5·10-5 s, the terminal velocity is reached (Figure 2e). 178 Figure 2: x (horizontal) component of the fluid velocity changes to pass by the particle at the first time-step (a) and the last computed time-step (b). z (vertical) component of the fluid velocity at the first time-step (c) is still, but at the last computed time-step (d) it increases as the diameter narrows where the particle takes place. After approximately 14,000 time-steps, the terminal velocity has reached (e). Our Particle Tracing based 3D CFD model resulted in six different settling velocities according to the extra counterforce. The value of the added force is identified to 8.28·10-5 N to balance the other forces and give the required settling velocity (Figure 3). Figure 3: a) Position of the particle in time, in case of different strength of the drag force, b) 8.28·10-5 N force would result in the settling velocity as measured (1.43·10-1 m/s) The results gained by the models are summarised in Table 1, which also helps to compare with the measured values. The IBM model is under development, but it gained an already proximate settling velocity, while the Particle Tracing model uses a fitted force parameter (making it some kind of grey box model), so we got right the expected value. Including the fluid’s effect onto the force calculating is a further step which can whiten the scheme. Table 1: Comparing the terminal velocities of the models and the experimental data Method Dimensionality Computational time Settling velocity Error IBM modelling 2D 2 h 0.15 m/s 5% Particle Tracing 3D 21 h 0.143 m/s 0% Video-based measurement 3D - 0.143 m/s 4. Conclusions Modelling the settling of a single solid bead in a Newtonian fluid is an important and challenging task. In this study, two modelling approaches are implemented, and the results are discussed. Particle Tracing is a built-in module of COMSOL Multiphysics, and it offers small editing. Achievement of the current study is that drag force is fitted according to material properties. One of the advantages of using commercial CFD software is the 179 containment of the viscosity to the modelling equations, while in our self-developed simulator it has not been implemented yet. Small perturbations and transient phenomena can already be examined by the model based on the Immersed Boundary Method, which gives a big added value to it, compared to empirical equations and simpler models to determine the terminal settling velocity. Mesh independence study would be a good extension to the current results, and it is among the near future plans. The simulator should be computationally more economic to make it able to calulate on larger computational domain. The accurately detailed model will help answer design and operational questions. Chemical and food industrial processes can be more plannable by understanding the exact processes during settling and the emission can be reduced, and energetically more efficient, greener technologies can be designed. Acknowledgements The research was supported by EFOP-3.6.1-16-2016-00015 Smart Specialization Strategy (S3) - Comprehensive Institutional Development Program at the University of Pannonia “Pannonia to Promote Sensible Individual Education and Career Choices” project. The financial support of MOL Nyrt. and Peregrinatio I. Foundations is gratefully acknowledged. References Batchelor, G.K., 1972. Sedimentation in a dilute dispersion of spheres. Journal of Fluid Mechanics 52, 245. doi: 10.1017/S0022112072001399 Fadlun, E.A., Verzicco, R., Orlandi, P., Mohd-Yusof, J., 2000. Combined Immersed-Boundary Finite-Difference Methods for Three-Dimensional Complex Flow Simulations. Journal of Computational Physics 161, 35–60. doi:10.1006/jcph.2000.6484 Ghosh, S., Stockie, J.M., 2015. Numerical Simulations of Particle Sedimentation Using the Immersed Boundary Method. Communications in Computational Physics 18, 380–416. doi: 10.4208/cicp.061113.050115a Juřena, T., Dohnal, M., Hájek, J., 2017. Simulation of multiphase flow in stirred tank of nonstandard design. Chemical Engineering Transactions 625–630. doi: 10.3303/CET1761102 Kerst, K., Roloff, C., Medeiros de Souza, L.G., Bartz, A., Seidel-Morgenstern, A., Thévenin, D., Janiga, G., 2017. CFD-DEM simulations of a fluidized bed crystallizer. Chemical Engineering Science 165, 1–13. doi: 10.1016/j.ces.2017.01.068 Kondrat’ev, A.S., Naumova, E.A., 2007. Calculation of the hindered settling velocity of a bimodal mixture of solid particles of arbitrary shape in a Newtonian fluid. Theoretical Foundations of Chemical Engineering 41, 216– 220. doi: 10.1134/S0040579507020170 Noetinger, B., 1989. A two fluid model for sedimentation phenomena. Physica A: Statistical Mechanics and its Applications 157, 1139–1179. dpi: 10.1016/0378-4371(89)90037-X Ochowiak, M., Matuszak, M., Włodarczak, S., Ancukiewicz, M., Krupińska, A., 2017. The modified swirl sedimentation tanks for water purification. Journal of Environmental Management 189, 22–28. Patankar, S.V., 1980. Numerical heat transfer and fluid flow, Series in computational methods in mechanics and thermal sciences. Hemisphere Publ. Co, New York. Patel, R.G., Desjardins, O., Kong, B., Capecelatro, J., Fox, R.O., 2017. Verification of Eulerian-Eulerian and Eulerian-Lagrangian simulations for turbulent fluid-particle flows. AIChE Journal 63, 5396–5412. doi: /10.1002/aic.15949 Peskin, C.S., 1972. Flow patterns around heart valves: A numerical method. Journal of Computational Physics 10, 252–271. doi: 10.1016/0021-9991(72)90065-4 Richardson, J.F., Zaki, W.N., 1997. Sedimentation and fluidisation: Part I. Chemical Engineering Research and Design 75, S82–S100. doi: 10.1016/S0263-8762(97)80006-8 Tran, K.A., Le, P.T.K., Pham, V.V., Nguyen, T.L.M., Nguyen, T.T., Tran, T.N., Le, K.A., 2018. Experimental and computational fluid dynamics investigation of rice husk updraft gasifier with various gasification agents. Chemical Engineering Transactions 223–228. doi: 10.3303/CET1863038 Van Santen, P., Augustinus, P.G.E.F., Janssen-Stelder, B.M., Quartel, S., Tri, N.H., 2007. Sedimentation in an estuarine mangrove system. Journal of Asian Earth Sciences 29, 566–575. Zhang, Y., Ye, M., Zhao, Y., Gu, T., Xiao, R., Liu, Z., 2015. Emulsion phase expansion of Geldart a particles in bubbling fluidized bed methanation reactors: A CFD–DEM study. Powder Technology 275, 199–210. Zhu, H.P., Zhou, Z.Y., Yang, R.Y., Yu, A.B., 2007. Discrete particle simulation of particulate systems: Theoretical developments. Chemical Engineering Science 62, 3378–3396. doi: 10.1016/j.ces.2006.12.089 180