Resolving Multiphase Flow through Packed Bed of Solid Particles

Sep 28, 2017 - Exploring a Multiphysics Resolution Approach for Additive Manufacturing. Alvaro Antonio Estupinan Donoso , Bernhard Peters. JOM 2018 70...
0 downloads 0 Views 2MB Size
Subscriber access provided by UNIV OF ESSEX

Article

Resolving Multiphase Flow through Packed Bed of Solid Particles Using eXtended Discrete Element Method with Porosity Calculation Maryam Baniasadi, and Bernhard Peters Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.7b02903 • Publication Date (Web): 28 Sep 2017 Downloaded from http://pubs.acs.org on October 1, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Industrial & Engineering Chemistry Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Resolving Multiphase Flow through Packed Bed of Solid Particles Using eXtended Discrete Element Method with Porosity Calculation Maryam Baniasadi ∗ and Bernhard Peters Université du Luxembourg. Faculté des Sciences, de la Technologie et de la Communication, Campus Belval, 2, Avenue de l'Universite L-4365, Esch-sur-Alzette, Luxembourg E-mail: [email protected];[email protected]

Phone: +532 4666445505. Fax: +352 46664435505

Abstract Multiphase ow reactors such as trickle bed reactors are frequently used reactors in many industries. Understanding the uid dynamics of these kind of reactors is necessary to design and optimize them. The pressure drop and liquid saturation are the most important hydrodynamic parameters in these reactors which depend highly on the porosity distribution inside the bed. The eXtended Discrete Element Method (XDEM) was applied as a numerical approach to model multiphase ow through packed beds of solid particles. This method has the ability to be coupled with Computational Fluid Dynamics (CFD) through interphase momentum transfer which makes it suitable for many Eulerian-Lagrangian systems. The XDEM also calculates the porosity distribution along the bed which not only eliminates the empirical correlations but also makes it possible to investigate the malddistribution of liquid saturation inside the bed. The

1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

results for the hydrodynamics parameters were compared with experimental data and satisfactory agreement was achieved.

Introduction Innovative and novel methods are necessary in order to reach a higher level of understanding of some critical issues. Multiphase ows through a packed bed of solid particles is one issue and exists in a broad spectrum of engineering disciplines such as chemical industries, petroleum engineering, wastewater treatment. One frequently used reactor of this type is a trickle bed reactor that usually contains particulate phase of which the interstitial space is lled with uid phases. Based on the direction of the uid ow they can be classied as cocurrent downow and counter-current trickle bed reactors. 1 The ow behavior in trickle bed reactors is very complicated and depends not only on hydrodynamics but also on mass and heat transfer. 2 The most important and eective parameters on the ow in trickle bed reactors can be classied as the physical properties of solid particles and uids, the geometry of bed and the operating conditions which make the ow more complex. 24 Moreover, the local interaction between all phases such as uid-uid, uid-solid and solid-solid interactions must be taken into account in order to better predict the performance of these reactors. Nowadays, numerical methods are widely applied to model trickle bed reactors. Some frequently used numerical methods will be briey discussed. Numerical approaches to model multiphase ow phenomena including droplets or solid particles may be classied into two categories: Continuum and Discrete. In the rst approach, all phases are treated as a continuum on macroscopic level; two-uid and Volume of Fluid (VOF) are the most well-known representatives. Almost all studies on trickle bed reactors use the continuum approach where the liquid is distributed uniformly and solid particles are completely (or partially) wetted. 2,5,6 Computational Fluid Dynamics (CFD) is widely used particularly based on the porous media concept to discretize the governing 2

ACS Paragon Plus Environment

Page 2 of 40

Page 3 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Nomenclature

b c d D eψ E1 , E2 Fij g K L m ˙ ij n N p p0 R I t T U V VREV x

band width speed of sound in liquid (ms−1 ) particle diameter (m) reactor diameter (m) mean absolute relative error Ergun's constant momentum exchange between phase i and j (J) gravitational accelerating (ms−2 ) Kernel function length (m) mass ow rate between phase i and j (Kgm−3 s−1 ) number of uid phases, number of neighbour cells number of experimental points, number of particles in a CFD cell pressure (P a) reference pressure (P a) reactor radius (m) and universal gas constant (Jmol−1 K −1 ) identity matrix time (s) temperature (K) velocity (ms−1 ) volume (m3 ) representative elementary volume (m3 ) coordinate of CFD cell or porticle (m)

Greek Symbols

∂ ∇ ∆ ρ ρ0 ν  τ µ βT βP φ α η ψ

dierential operator nabla operator dierence density (Kgm−3 ) reference density (Kgm−3 ) velocity (ms−1 ) volume fraction stress strain tensor (Kgm−1 s−2 ) dynamic viscosity (Kgm−1 s−1 ) thermal expansion rate (K −1 ) compressibility rate (P a−1 ) porosity volumetric saturation weight of a particle variable

Superscripts and Subscripts

f g i, j l p T c

uid gas phase liquid particle transpose CFD cell

3

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

equations. 3,4,712 The two phase Eulerian model describes the ow of gas and liquid through a packed bed on the basis of volumetric phase average concept. 10 In the second approach, discrete techniques are applied to the discontinuous entities e. g. liquid droplets and/or solid particles. Therefore, they represent the discrete structure of liquid ow 2,11 or solid particles in packed bed reactors. In the latter case, CFD is used to solve the governing equations of the uid phases and Discrete Element Method (DEM) calculates the position, orientation and velocity of the each particle. The literature on treating trickle bed reactors by the CFD-DEM method are scarce as was stated by Wang et al., 13 while several studies were conducted on the packed bed reactors including either the gas or the liquid phase. 1416 The combined CFD-DEM method is more ecient where the role of solid particles becomes more dominant. For instance when there is a chemical reaction or heat and mass transfer in the interface of uid phases and solid particles. On the other hand, the discrete element method provides a useful framework to quantitatively understand the inuence of particle size and shape on characteristics of the randomly packed beds. 17 Despite dierent mathematical models for trickle bed reactors, there is still no universal agreement on the right model for this phenomena. We know that predicting the ow path of dierent uid phases in packed beds plays a vital role in improving productivity as well as achieving a stable operation. 18 Understanding the mechanism of ow maldistribution on one side and solid particles on the other side is important for the improvement of the performance of industrial reactors. The critical question is how to eectively produce reliable data to provide a model for a trickle bed reactor. In order to answer this question and develop a strategy to model several liquid ows through packed beds in chemical reactors, a coupling method between continuous and discrete approached may be applied. To the best knowledge of the authors, the unresolved CFD-DEM method has not been developed yet for trickle bed reactors. In this paper, therefore, we present the application of the eXtended Discrete Element Method (XDEM) to trickle bed reactors. The XDEM developed by Peters 19 is an extension of classical DEM that can be coupled with OpenFOAM. Although the XDEM 4

ACS Paragon Plus Environment

Page 4 of 40

Page 5 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

platform has already solved the heat and mass transfer in a packed bed reactor for a multi species gas phase through solid particles, 20 it has not been yet developed for several uid phases. Therefore, the XDEM method in this study was extended to multiphase ow to predict pressure drop and liquid saturation for packed bed reactors without considering heat and mass transfer for the moment. This novel multiphysics, multi-scale model can be applied to many industrial processes where the role of solid particles become more dominant. This method has the potential to include all transport phenomena such as heat and mass transfer between all phases which are the motivations for future works. In this particular study, the hydrodynamics of trickle bed reactors was investigated using the XDEM. Here, the focus is mainly on the understanding of gas and liquid ow distribution, spatially and temporally, through a packed bed of particles using an Eulerian-Lagrangian approach. One challenging parameter in trickle bed reactors is the porosity distribution (axial and radial) which can have a considerable eect on the calculated hydrodynamic parameters specially where the ratio of reactor diameter over particle diameter ( dDp < 20) is low. 21 Although other studies benet from dierent correlations to take this eect into consideration, 3,4,7,21 the XDEM does not require these correlations because it can calculate the porosity distribution inside the bed by knowing the number and position of the particles in each CFD cell. The empirical correlations are mainly based on averaged properties with homogeneous assumption, 16 while the XDEM provides a realistic, random packing structure for the packed bed. This advantage of the XDEM makes it suitable for any kind of reactor including particles. The most important hydrodynamic parameters of a trickle bed reactor are pressure drop and liquid/gas hold up which highly depend on bed geometry, interaction between phases, physical properties and kinetic characteristics in case of reaction. Therefore, the main aim of this study is to develop a discrete-continuous model for a system of gas and liquid through a packed bed of particles that can be extended to a system with any number of phases including heat and mass transfer between uid phases and/or uid and solid phases.

5

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 40

Mathematical model and numerical solution The performance of multiphase reactors can be modelled by solving the conservation equation for mass, momentum and energy in combination with species transfer, phase change and chemical reaction. Due to the highly coupled and nonlinear nature of the governing equations, a complete solution depends on a reliable physical model, an advanced numerical algorithm and sucient computer power. 4 In this study we consider cocurrent gas-liquid ow through a packed bed of solid particles. The solid particles are immobile and nondeforming like those in packed and trickle bed reactors. These phases are categorized according to their dierent physical quantities. The uid phases consists of dierent components, but in this current study the components are neglected and only the phases are considered. In order to model a system with the above mentioned assumptions, a coupled discrete-continuous method based on the XDEM and CFD was applied to model trickle bed reactors. The XDEM is based on classical DEM method 22 in which additional properties such as thermodynamic state, stress/strain or electro-magnetic elds for each particle are also modelled. This discrete method has the ability to be coupled with continuous numerical methods such as CFD and FEM to address many challenging engineering problems. 20,2326 In the application of the XDEM to trickle bed reactors, gas and liquid phases are considered as continuous media where solid particles are tracked in the Lagrangian framework. Therefore, an Euler-Lagrange model is proposed to describe the cocurrent ow of gas and liquid in both trickle bed and packed bed reactors. The entire process including uid phases and solid particles may be represented in the following way: Entire Process =

X

Particle processes + CFD

(1)

In the XDEM, the solid particles which can have dierent type of shapes 20 are considered as a porous medium, where uid phases are either trapped or passed in the porous structure and allowed to exchange heat and mass transfer with solid particles. In the packed bed reactors 6

ACS Paragon Plus Environment

Page 7 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

the position and orientation of each particle is specied using the XDEM. The conservation of mass, momentum and energy, as well as species transfer within the pore volume of porous particles, is described by transient and one-dimensional dierential conservation equations. For more details about the XDEM method, the reader is referred to. 19,20,27 In the uid phase, the governing equations, such as conservation of mass and momentum, are characterized as a continuum phase through porous media. It is not always possible to track the local properties at the interfaces in multiphase systems. Therefore, an appropriate averaging is required to obtain macroscopic governing equations. It must be pointed out that the macroscopic aspect of multiphase ow is more important to the design and operation of a multiphase system. 28 In order to achieve this goal, the multi-uid model of volume average method, that performs spatial averaging for each individual phase was applied.

Figure 1: Representative Element Volume(REV) for volume averaging method Since Eulerian volumetric average is the most important and widely-used method of averaging, the conservation equations of mass, momentum, energy and species of the current project are derived based on this method. The Eulerian volumetric average is usually performed over a Representative Element Volume (REV) in the ow as shown in Figure 1. For multiphase ow the volume fraction of ith phase, i is dened as follows: i =

Vi VREV

(2)

where Vi is the volume of each phase and VREV is the representative elementary volume. 7

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 40

The volume fraction of all phases must sum to one. If the denition of volumetric uid phases saturation (αi ) for each uid phase dened as Eq. (3) the volume fraction of each uid phase can be rewritten in terms of porosity ( φ) as Eq. (4). Vi αi = Pn

i=1

(3)

Vi

Here, n is number of uid phases and porosity is φ =

Pn

V i=1 i VREV

.

i = φαi

(4)

Therefore, the mass conservation equation can be derived for each uid phase separately: Conservation of mass n X ∂ (φαi ρi ) + ∇ · (φαi ρi v~i ) = m ˙ ij ∂t j=1

(5)

where ρi and v~i are the phase density and the phase intrinsic velocity 28 respectively. The right hand side of the continuity Eq. (5) represents mass transfer per unit volume to other phases due to phase change, which is often non-zero when one phase is transferred to another phase. The momentum equation of uid phases consists of unsteady, convection, pressure gradient, diusion and the source term is as follows: Momentum Equation n X ∂ (φαi ρi v~i ) + ∇ · (φαi ρi v~i v~i ) = −φαi ∇p + φαi ρi g + ∇ · τ¯i + Fij ∂t j=1

(6)

where p is the pressure, g is gravitational acceleration and τ¯i is the stress strain tensor based on Newton's law of viscosity for phase i and dened as Eq. (7). The last term in Eq. 8

ACS Paragon Plus Environment

Page 9 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(6) is the interphase momentum exchange between every two phases, which is the only force considered in this modelling. The most well-known and used drag force model for the uid phase through a packed bed of solid particles is Ergun equations and it is derived as Eq. (8). 2 τ¯i = φαi µi [∇~ vi + ∇~ vi T ] − φαi µi ∇ · v~i I 3

"

Fij

= φαi

E2 (1 − φαj )ρi E1 (1 − φαj )2 µi (~ vi − v~j ) + |~ vi − v~j |(~ vi − v~j ) 2 2 (φαi ) dj φαi dj

(7)

#

(8)

The closure of the system of equations can be achieved by providing a sucient number of equations to determine all unknowns. There isn't still an independent equation to calculate pressure eld in the system. With the assumption that all phases share the same pressure, the pressure equation can be derived by summation over mass conservation equations like Eq. (5) of all uid phases, and the densities should be replaced by pressure according to their relationship based on equations of state. Total mass balance for a two phase system (liquid and gas)

∂ (φαg ρg + φαl ρl ) + ∇ · (φαg ρg v~g ) + ∇ · (φαl ρl v~l ) = 0 ∂t

(9)

Pressure and density are related by an equation of state, which in this study perfect gas is considered for gas phase and constant speed of sound for liquid phase: • for ideal gas

ρg =

pg RT

where R is the global gas constant. • liquid phase

9

ACS Paragon Plus Environment

(10)

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ρl = ρ0,l (1 + βT T + βP (pl − p0,l ))

Page 10 of 40

(11)

In this equation, βT is the Thermal expansion rate which is equal to zero in isothermal cases and βP is Compressibility rate of liquid, which is dened as: βP =

1 1 dρl = ρl dpl ρ l c2

(12)

where c is the speed of sound in liquid phase. Then, the equation of state for isotherm liquid phase could be written as: ρl = ρ0,l +

1 (pl − p0,l ) c2

(13)

where ρ0,l and p0,l are reference density and pressure respectively. In this particular application of the XDEM, the coupling between uid phases and solid particles is due to the drag force that solid particles are inserting on uid phases as Eq. (8). In order to calculate these forces, porosity and particle diameters for each CFD cell are calculated by the XDEM.

Porosity calculation The void fraction in xed bed reactors is not uniform and varies from the vicinity of the walls toward the bulk of the bed, 3 and it needs to be incorporated into the model to also consider the wall eects. The XDEM has the ability to provide porosity distribution inside the ow eld by considering the position and number of particles in each CFD cell. 30 This advantage eliminates the need for correlations to calculate porosity distribution and in consequence, it predicts the uid properties such as velocity and volume fraction in a more accurate and precise way. The porosity calculation algorithm is based on the method proposed by Xiao et al. 31 which considers that the volume of each particle is distributed to all the neighbour cells according to its distance to the cell centers. This model is modied in a way to distribute 10

ACS Paragon Plus Environment

Page 11 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

the particle volume among the cells, which particle touches. Thus, the volumetric average porosity for each CFD cell is calculated based on the number of particles in the CFD cell with their corresponding weights.

φc = 1 −

N 1 X ηj,c Vj Vc j=1

K(|xj − xc |/b) m=1 K(|xj − xc,m |/b)

ηj,c = Pn

(14)

(15)

Where ηj,c is the weight of the j th particle to the given cth CFD cell. N is the number particles whose centers are located within the cth CFD cell. n is the number of touching neighbour cells. A kernel-based interpolation procedure is used to normalize the weight function, which is noted as K , discussed in more details by Xiao et al. 31 x is the coordinate of particles and CFD cells and b is the band with, which is suggested to be chosen between 1.5∆ and 2.5∆,where ∆ is the average cell dimension.

For instance, in the upper left corner of Figure 2, two particles ( P1 , P2 ) are shown, which belong to the cell 7th . The volume of particle P1 is going to be distributed among four touching neighbour cells ( C6 , C7 , C10 , C11 ) while the volume of particle P2 is going to distributed among C3 , C7 . Therefore, the calculated porosity for C7 can be written as:

φC7 = 1 −

1 (ηP ,c VP + ηP2 ,c7 VP2 ) Vc7 1 7 1

(16)

In this gure the calculated porosity by the XDEM for a system of two particles and sixteen CFD cell are shown (the right side).

Solution procedure The computational domain according to the experimental study of Gunjal et al. 32 on trickle bed reactors is depicted in Figure 3 with the lenght of L = 1m and diameter of D = 0.11m. Particles with specied diameters are introduced into the domain and settled due to the 11

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 40

Figure 2: Illustration of the averaging procedure for porosity calculation by the XDEM for a system of two particles and sixteen CFD cells in both 2D and 3D view.

Gas

Liquid

(a)

(b)

(c)

Figure 3: (a) Geometry of the case based on experimental study of Gunjal et al. 32 (L = 1m, D = 0.11m), (b) dp = 0.003m , (c) dp = 0.006m. 12

ACS Paragon Plus Environment

Page 13 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

gravitational force to reach the desire height, which is one meter in this study using the DEM method. Complete description of the DEM module of the XDEM is described in, 33,34 which is not the focus of this study. This part is the most expensive part of the whole simulation in terms of computational time. For instance, the case with 3mm particle diameter requires 2160 min CPU time with 10 cores in parallel using Gaia cluster of the University of Luxembourg. The packing structure can be classied as poured random packing as dened by Klerk. 35 The gas and liquid phases are introduced from the top (liquid entering through 25 percent of the inlet area) with uniform inlet velocities. No slip boundary condition was used for the walls and xed atmospheric pressure was specied for the outlet boundary condition. The same discretization schemes were used for all cases where the rst order implicit Euler scheme was applied to temporal derivatives and Gauss linear or limitedlinear (rst order) for convective and Laplacian terms. In order to solve the governing equations, we assumed that: • there is no inter-phase mass transfer. • the system is isotherm. • compressible, where the liquid phase is weakly compressible in order to treat uid

phases in the same way. • all uid phases share the same pressure since prewetted bed is assumed. 32 In this

study the capillarity pressure can be neglected since in the trickling zone there is a low interaction between the gas and liquid phase. 3 • the turbulence eects are also negligible 12 for both liquid hold up and pressure drop

as the liquid and gas ow rates studied in this work are following the trickling zone. This assumption was also made by other authors. 3,4 In this numerical model, the size of the CFD cell must be larger than the particle diameter known as unresolved model. 36 This limitation made this method not to be suitable for 13

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

very low

D dp

ratio but for some pilot scale and industrial scale applications. The segregated

projection algorithm, based on the pressure implicit of operators (PISO) algorithm implemented in OpenFOAM is used to obtain variable elds at dierent time steps. 37 The PISO algorithm is based on time stepping procedure that solves the momentum equation using pressure from the previous time step. The predicted velocity eld should satisfy the total mass equation as well. This is achieved by substituting the uxes derived from the predicted velocity eld into the total mass equation. The solution of the resulting equation called pressure correction equation, would be predicted pressure eld. In the next step, velocities and uxes are corrected according to the new pressure eld. The solution algorithm for the coupled XDEM solver can be summarized as Table 1. Table 1: Solution algorithm. Calculate porosity and mean particle diameter for each CFD cell by XDEM For each time step, do For each outer corrector loop, do 1. Compute momentum exchange force Eq.(8) 2. Solve momentum equation using pressure from previous time step Eq.(6) 3. Solve pressure equation Eq.(9) 4. Update velocities and fluxes based on new pressure field 5. Compute volume fraction Eq.(5) 6. Update densities Eqs.(10, 11) End outer corrector loop End time step

The suggested number of outer corrector loop 37 is 4 to 6 loops, where in this study, four outer corrector loops were used to reach satisfactory residuals. Preliminary numerical cases with dierent cell numbers (see Figure 6) were set up to study grid independence. In Figure 3, the calculated porosity distribution by XDEM for each grid resolution is also shown for a slice in z direction. The average axial and radial porosity calculated by XDEM for dierent number of CFD grids are shown in Figures 4,5. 14

ACS Paragon Plus Environment

Page 14 of 40

0.9

6400 12500

0.75

Porosity

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

21600 29400

0.6

0.45 0.3 0

0.2

0.4

0.6

0.8

1

length

Figure 4: Average axial variation in porosity calculated by the XDEM for dierent number of CFD cells (dp = 0.006m).

1 6400

0.9

porosity

Page 15 of 40

12500

0.8

21600

0.7

29400

0.6 0.5 0.4 0.3 0

0.01

0.02

0.03

0.04

0.05

radius

Figure 5: Average radial variation in porosity calculated by the XDEM for dierent number of CFD cells (dp = 0.006m).

15

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The calculated pressure drop and liquid hold up for dierent number of cells are depicted in Figure 7 and 8, where the optimal grid size of 21600 cells were used for further evaluations since no signicant changes were observed by rening the mesh. Adjustable time step was used to satisfy Courant number of less than 0.5 and residuals of less than 104 were accepted.

Figure 6: Geometry of the case with dierent cell numbers and the calculated porosity for a slice in bed

Results and discussion In order to investigate the eect of porosity distribution on the hydrodynamic parameters of trickle bed reactors, the XDEM was used to compute the porosity eld for the computational domain. Initially, the particles were placed in the domain (see Figure 3), then porosity was calculated based on the position of the particles in each CFD grid. The calculated porosity eld for two dierent particle sizes ( 0.003m, 0.006m) are depicted in Figure 9, which repre16

ACS Paragon Plus Environment

Page 16 of 40

10 8 6 4 2 0 0

10000

20000

30000

40000

Number of cells

Figure 7: Pressure drop for dierent number of cells

0.19

Liquid hold up

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Pressure drop [KPa/m]

Page 17 of 40

0.185

0.18 0.175 0.17 0

10000

20000

30000

40000

Number of cells

Figure 8: Liquid hold up for dierent number of cells

17

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

sents higher porosities at the vicinity of walls and lower porosities in the bulk of particles. In this gure, the porosity variation in terms of three components (axial, radial and azimuthal) is shown, which reveals the benet of this model by covering more detailed information compare to the porosity calculation correlations specially in the azimuthal direction. The local porosity calculated by the XDEM at dierent angles on a slice at Z = 0.5m for 0.006m, are depicted in Figure 10. Although, The local radial porosity shows an oscillation behaviour, the average radial porosity calculated by the XDEM obey an exponential behaviour.

Figure 9: Calculated porosity distribution by the XDEM, (a) dp = 0.003m , (b) dp = 0.006m. Therefore, the XDEM results for the porosity calculation is plotted against two exponential correlations 38,39 in radial direction (see Figure 11). Almost all exponential correlations predict constant porosity in the bulk of the particles while the XDEM porosity calculation does not, as shown in Figure 11. In order to compare the XDEM results with the experi18

ACS Paragon Plus Environment

Page 18 of 40

Page 19 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 10: Local porosity calculated by the XDEM over the dimensionless radius at three dierent angles (dp = 0.006m).

19

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

mental observations, a case with

D dp

= 16.77 was constructed as shown in Figure 12 and the

XDEM results were compared with the experimental data reported by Goodling et al. 40 The main reason for inability of the XDEM to capture the complete porosity oscillation near the wall, is the restriction of the CFD cells to the particle diameter in the unresolved method. The total average porosity calculated by the XDEM is 0.432, while the experimental data is 0.415. The error could be due to the dierence in the packing structure. Despite this shortness, the XDEM method is still a reliable methods since it can also be applied to the cases with dierent particle sizes 41 and also can treat non-spherical particles such as cylinders. 42 For the simplicity, the authors have started with spherical particles with constant particle diameters based on the experimental studies and considering other particle congurations are the aim of future contributions.

1.05 Sun et al. [38] Hunt and Tien [39]

0.9

XDEM prediction

0.75

porosity

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 40

0.6

0.45 0.3 0

1

2

3

4

5

6

7

8

9

(R-r)/dp

Figure 11: Comparison between the XDEM prediction for average radial porosity variation and exponential correlations ( dp = 0.006m). Accurate prediction of porosity distribution is an asset because it aects the simulation results highly, which will be discussed later in this study. 20

ACS Paragon Plus Environment

Page 21 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 12: Comparison between the XDEM prediction for average radial porosity variation and experimental data 40 for dDp = 16.77. In the rst step, the model was validated for a single phase of gas and liquid through a packed bed of spherical particles using experimental work of Li and Ma. 43 The experiment was carried out in a pipe with the inner diameter of 90mm and the height of 635mm, lled with spherical particles. Table 2: Properties of gas and liquid phases. Gas Phase:

ρ = 1.18 µ = 0.000017

Air kg/m3 Density kg/(ms) Viscosity

Liquid Phase:

ρ = 1000 µ = 0.001 c = 1482

Water

kg/m3 Density kg/(ms) Viscosity m/s speed of sound in 20◦ C

The liquid and gas properties based on water and air properties at room temperature and atmospheric pressure are listed in Table 2. Since the pressure drop in packed bed reactors is highly sensitive to bed porosity and particle diameter, the model was validated for three dierent particle diameters. In Figure 13 the calculated pressure drop against experimental 21

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

25

Pressure drop [Kpa/m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 40

dp = 0.006 m [43] dp = 0.003 m [43] dp = 0.0015 m [43] XDEM prediction

20

15

10

5

0 0

0.2

0.4

0.6

0.8

1

1.2

Gas velocity [m/s]

Figure 13: Pressure drop of gas ow in packed bed for three dierent particle diameters data for single gas phase with air properties is shown and the same data for single liquid phase with water properties is depicted in Figure 14. The XDEM predictions holds in fairly good agreement with the experimental data. As was shown in Eq. (8), the semi empirical model of the Ergun equation has two terms, the rst term is known as the Viscous loss and the second term is known as Inertial loss. In dierent studies, dierent sets of Ergun constants were used 3,7,21,44 in order to nd the optimum values for these constants. For instance, Illiuta et al. 44 used E1 = 312 and E2 = 1.98 for the bed porosity of  = 0.38 and E1 = 178 and E2 = 1.2 for the bed porosity of  = 0.35, while Gunjal et al. 7 applied E1 = 215 and E2 = 1.8 for a 3mm particle and higher values of E1 = 500 and E2 = 2.4 for 6mm particles. Janecki et al. 21 studied and compared the results for three dierent sets of Ergun constants. In another study, Atta et al. 3 considered E1 = 180 and E2 = 1.8 for all their cases. Dierent studies show a wide acceptable range of

Ergun constants; about 150 − 500 for E1 and 1.2 − 5 for E2 . In this modeling, dierent sets 22

ACS Paragon Plus Environment

Page 23 of 40

14 dp = 0.006 m [43] dp = 0.003 m [43]

12

Pressure drop [Kpa/m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

dp = 0.0015 m [43]

10

XDEM prediction

8 6 4 2 0 0

0.005

0.01

0.015

0.02

Liquid velocity [m/s]

Figure 14: Pressure drop of liquid ow in packed bed for three dierent particle diameters of Ergun constants were examined for inlet liquid velocity of 0.002m/s and 0.01m/s shown in Figures 15 and 16 against experimental data of Gunjal et al. 32 These gures show the most reliable values in this modelling are close to the ones used by Gunjal et al.; 32 E1 = 215 and E2 = 2 for 3mm particles and E1 = 500 and E2 = 3.5 for 6mm particles since the same experimental data were used. Changes in operational conditions such as gas and liquid velocity, could highly aect the pressure drop and liquid saturation in packed bed reactors, where both these factors are of high importance to improve the productivity and eciency of the reactor. Therefore, several cases were set up based on experimental data 32 to nd the evolution of pressure drop and liquid saturation in packed bed reactors. The model predictions are compared to experimental data in Figure 17 for pressure drop and Figure 18 for liquid holdup for constant inlet gas velocity of 0.22m/s. The mean absolute relative error < eψ > with the denition expressed in Eq. (17) 45 gives 0.1047 and 0.0396 for pressure drop and liquid holdup 23

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

14

Pressure drop [KPa/m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 40

dp = 0.003 m [32] E1 = 200, E2 = 2 E1 = 200, E2 = 2.5 E1 = 215, E2 = 2

12 10 8

6 4 2 0 0

0.002

0.004

0.006

0.008

0.01

0.012

Liquid velocity [m/s]

Figure 15: Pressure drop for dierent sets of Ergun constant for dp = 0.003m respectively.

< eψ >=

N 1 X |ψexperimental − ψmodel | N 1 ψexperimental

(17)

where ψexperimental represents the value of the variables such as pressure drop or liquid saturation from experimental data and ψmodel is the calculated value from the XDEM results. The prediction accuracy of the XDEM model has been found through low absolute relative errors for both pressure drop and liquid hold up. Fairly good agreement between the XDEM results and experimental data in Figures 17, 18 also shows this matter. The importance of porosity distribution on calculated results is investigated by comparing the constant porosity results versus the XDEM prediction results (see Figures 19 and 20). The average balk porosity calculated by the XDEM method is 0.412 for dp = 0.003m and 0.392 for dp = 0.006m, which were considered for these gures. The mean absolute relative

24

ACS Paragon Plus Environment

Page 25 of 40

5

Pressure drop [KPa/m]

dp = 0.006 m [32] E1 = 500, E2 = 3 E1 = 500, E2 = 3.5 E1 = 300, E2 = 4

4 3 2

1 0 0

0.002

0.004

0.006

0.008

0.01

Liquid velocity [m/s]

Figure 16: Pressure drop for dierent sets of Ergun constant for dp = 0.006m

12

Pressure drop [KPa/m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

dp = 0.006 m [32]

10

dp = 0.003 m [32]

8

XDEM prediction

6 4

2 0 0

0.002

0.004

0.006

0.008

0.01

0.012

Liquid velocity [m/s]

Figure 17: Experimental data for pressure drop against the XDEM prediction 25

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

120.2

[KPa/m] drop Pressure up hold Liquid

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 40

dp = 0.006 m [32]

10 0.15 8

dp = 0.003 m [32] XDEM prediction

60.1 4 0.05 2 00 00

dp = 0.006 m [32] dp = 0.003 m [32] XDEM prediction

0.002 0.004 0.004 0.006 0.006 0.008 0.008 0.002

0.01 0.01

0.012 0.012

Liquidvelocity velocity[m/s] [m/s] Liquid

Figure 18: Experimental data for liquid hold up against the XDEM prediction error of 0.6707 for constant porosity results indicates large dispersion of computed results with respect to experimental values. As shown in Figures 19 and 20, the constant porosity model is over predicting the results since the porosity magnitude is not accurate at walls and also in the bulk of solid particles. In order to perform further analysis of the dierence between constant porosity results and XDEM prediction results, the porosity distribution, liquid saturation and intrinsic liquid velocity for a single case is shown and compared in Figure 21. For a part of the bed indicated in Figure 21, the three top gures are calculated results obtained by XDEM, while the other cases are results based on constant porosity. The velocity prole and liquid saturation for the constant porosity bed seems to be higher than the XDEM prediction. This comparison reects the importance of porosity distribution inside the bed and its eect on the calculated parameters. Therefore, it is crucial to consider the porosity prole for trickle bed reactors, specially at the wall region, to predict the hydrodynamic parameters more accurately. Introducing the porosity distribution for each CFD cell by the XDEM, provides the 26

ACS Paragon Plus Environment

Page 27 of 40

0 0

Pressure drop [KPa/m]

20

0.002

0.004

0.006

0.008

0.01

Liquid velocity [m/s]

dp = 0.003 m [32] XDEM prediction

16

porosity = 0.412

12 8

4 0 0

0.002

0.004

0.006

0.008

0.01

0.012

Liquid velocity [m/s]

Figure 19: Comparison between pressure drop computed by the XDEM and constant porosity against experimental data of dp = 0.003m.

dp = 0.003 m [32]

Liquid hold up

12

dp = 0.006 m [32] XDEM prediction 10 XDEM prediction porosity = 0.392

Pressure drop [KPa/m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

8 6

0.15 0.1

0.05

4

0

2 0.002

0.004

0.006

0.008

0.01

0.012

Liquid velocity [m/s]

0

0

0.002

0.004

0.006

0.008

0.01

0.012

Liquid velocity [m/s]

Figure 20: Comparison between pressure drop computed by the XDEM and constant porosity against experimental data of dp = 0.006m. 27

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 21: Comparison between the eect of the XDEM porosity distribution and constant porosity of 0.382 on liquid saturation and intrinsic liquid velocity.

28

ACS Paragon Plus Environment

Page 28 of 40

Page 29 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

opportunity to predict the liquid phase maldistribution inside the packed bed. Figure 22 reects the eect of the liquid inlet velocity on liquid saturation and liquid maldistribution. This gure shows that higher liquid saturation occurs with an increase of inlet liquid velocity. The increase in inlet liquid velocity enhances the interaction between liquid and gas phase in addition to interaction between liquid and solid phase. These interactions greatly inuence the predicted liquid hold ups.

Figure 22: Eect of inlet liquid velocity on liquid saturation. The intrinsic liquid velocity vectors, as well as intrinsic gas velocities for the same case, are demonstrated in Figures 23 and 24 respectively. These gures show the eect of porosity distribution on velocities which are not as uniform as constant porosity. An increase of inlet liquid velocity leads to a higher pressure drop (see Figure 25) due to the increase in the shear stress among phases. 44 In order to show the evolution of liquid phase maldistribution inside the bed, liquid saturation contours and intrinsic gas velocity vectors for dierent times are shown in Figure 29

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 23: The eect of inlet liquid velocity on intrinsic liquid velocity.

Figure 24: The eect of inlet liquid velocity on intrinsic gas velocity. 30

ACS Paragon Plus Environment

Page 30 of 40

Page 31 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 25: The eect of inlet liquid velocity on pressure drop along the bed. 26 to reach steady state condition. More detailed information can be derived from this gure, for instance it shows the interaction between the liquid and gas phases as well as solid particles and uid phases. The experimental data of Gunjal et al. 32 has previously been used by several researchers to validate their simulation results. 3,21,32,46 In Figures 27 and 28, the XDEM prediction results for the pressure drop and liquid hold up have been compared to the simulated results of Gunjal et al. 32 and Atta et al. 3 The XDEM results agreed well with the experimental data in comparison with other studies. The reason for better prediction of the proposed model is the correct and accurate calculation of the porosity distribution while other studies mostly used available correlations for porosity calculation. The XDEM porosity calculation describes the hydrodynamics parameters in a more realistic way. This Eulerian-Lagrangian characteristic of the XDEM method make it suitable not only for predicting porosity distribution but also for including heat and mass transfer between 31

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

Time [s]

0

5

10

15

Figure 26: Contours of liquid saturation and intrinsic gas velocity vectors for dierent times.

Pressure drop [KPa/m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 40

5

dp = 0.006 m [32] Gunjal et al . [32] Atta et al. [3] XDEM prediction

4 3 2 1 0 0

0.002

0.004

0.006

0.008

0.01

Liquid velocity [m/s]

Figure 27: Comparison of XDEM results of pressure drop with data reported in the literature. 32

ACS Paragon Plus Environment

Page 33 of 40

0.2

Liquid hold up

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

dp = 0.006 m [32] Gunjal et al . [32] Atta et al. [3] XDEM prediction

0.15

0.1

0.05 0

0.002

0.004

0.006

0.008

0.01

Liquid velocity [m/s]

Figure 28: Comparison of XDEM results of liquid hold up with data reported in the literature. solid particles and uid phases which distinguish this model from classical models. The mentioned advantage of the XDEM makes it suitable for many chemical reactors where the interaction between the uid phases and the particulate phase is important which is going to be outlined in future studies.

Conclusions In this contribution, a numerical model for multiphase ow through a packed bed of particles using continuous-discrete approach was investigated. The eXtended Discrete Element Method was applied as a numerical method coupled with OpenFoam to treat the uid-solid phase interactions using the Ergun semi empirical model. The XDEM provides the distribution of porosity and the mean diameter of spherical particles for each CFD cell to calculate the drag force on the uid phases and resulting pressure drop and liquid hold up calculation along the bed. This Eulerian-Lagrangian model make it possible to solve the governing equa33

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

tions for the continuous phase using CFD and the discrete phase using the XDEM where the grid size of the CFD cells must be larger than the discrete phase size. The model eliminates the need to use correlations to describe porosity distribution. The XDEM results were validated for single phases of gas and liquid as well as a two phase system of liquid-gas through a packed bed of spherical particles using experimental data of Li and Ma 43 and Gunjal et al. 32 The uids were water and air at atmospheric pressure, which are the most extensively tested uids. The eect of liquid and gas velocities on the hydrodynamic parameters for instance pressure drop and liquid saturation, were investigated. The comparison of the XDEM prediction and experimental results shows very good agreement which shows the prediction capability of the model for single and two phase ows. The results show that the higher uid velocities result in higher pressure drops along the bed. It was also found that the liquid saturation is higher at high liquid velocities. Furthermore, the porosity distribution provided by the XDEM made it possible to predict liquid spatial distribution in packed beds. The eect of porosity distribution on liquid and gas velocities as well as on liquid saturation were also discussed and compared to constant porosity results. This comparison shows the importance of porosity distribution and the accuracy of the XDEM results. In other words, the results prove that the XDEM model is able to predict the pressure drop and liquid saturation for multiphase ow through packed bed of solid particles accurately.

Acknowledgements The authors would like to thank and express their appreciation to the Luxembourg National Research Fund (FNR) and High Performance Computers (HPC) team at the University of Luxembourg (https://hpc.uni.lu/systems/gaia/ ) for supporting this research. They also would like to thank Florian Homann and Xavier Besseron for their scientic advice. 34

ACS Paragon Plus Environment

Page 34 of 40

Page 35 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

References (1) Mederos, F. S.; Ancheyta, J. Review on criteria to ensure ideal behaviors in trickle bed reactors. Appl. Catal., A

2009,

355, 1.

(2) Schwidder, S.; Schnitzlein, K. A new model for the design and analysis of trickle bed reactors. Chem. Eng. J.

2012,

207-208, 758.

(3) Atta, A.; Roy, S.; Nigam, K. D. Investigation of liquid maldistribution in trickle-bed reactors using porous media concept in CFD. Chem. Eng. Sci.

2007,

62, 7033.

(4) Jiang, Y.; Khadilkar, M.; Al-Dahhan, M.; Dudukovic, M. CFD of multiphase ow in packed-bed reactors: I. k-Fluid modeling issues. AlChE J.

2002,

48, 701.

(5) Bazmi, M.; Hashemabadi, S.; Bayat, M. CFD simulation and experimental study of liquid ow mal-distribution through the randomly trickle bed reactors. Int. Commun.

Heat Mass Transfer

2012,

39, 736.

(6) Hamidipour, M.; Chen, J.; Larachi, F. CFD study and experimental validation of trickle bed hydrodynamics under gas, liquid and gas/liquid alternating cyclic operations .

Chem. Eng. Sci.

2013,

89, 158.

(7) Gunjal, P. R.; Ranade, V. V.; Chaudhari, R. V. Liquid distribution and RTD in trickle bed reactors: experiments and CFD simulations. Can. J. Chem. Eng.

2003,

81, 821.

(8) Nikou, M. K.; Ehsani, M. Turbulence models application on CFD simulation of hydrodynamics, heat and mass transfer in a structured packing. Int. Commun. Heat Mass

Transfer

2008,

35, 1211.

(9) Jiang, Y.; Khadilkar, M.; Al-Dahhan, M.; Dudukovic, M. CFD of multiphase ow in packed-bed reactors: II. Results and applications. AlChE J.

2002,

48, 716.

(10) Beni, A. H.; Khosravi-Nikou, M. R. Modeling Hydrodynamics of Trickle-bed Reactors at High and Low Pressure Using CFD Method. Pet. Sci. Technol. 35

ACS Paragon Plus Environment

2015

, 33, 1770.

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 40

(11) Jiang, Y.; Guo, J.; Al-Dahhan, M. H. Multiphase Flow Packed-Bed Reactor Modeling: Combining CFD and Cell Network Model. Ind. Eng. Chem. Res.

2005,

44, 4940.

(12) Lopes, R. J. G.; Quinta-Ferreira, R. M. Volume-of-uid-based model for multiphase ow in high-pressure trickle-bed reactor: Optimization of numerical parameters. AIChE J. 2009,

55, 2920.

(13) Wang, Y.; Chen, J.; Larachi, F. Modelling and simulation of trickle-bed reactors using computational uid dynamics: A state-of-the-art review. Can. J. Chem. Eng.

2013,

91, 136. (14) Tabib, M. V.; Johansen, S. T.; Amini, S. A 3D CFD-DEM Methodology for Simulating Industrial Scale Packed Bed Chemical Looping Combustion Reactors. Ind. Eng. Chem.

Res.

2013,

52, 12041.

(15) Blais, B.; Lassaigne, M.; Goniva, C.; Fradette, L.; Bertrand, F.Development of an unresolved CFD-DEM model for the ow of viscous suspensions and its application to solid-liquid mixing. J. of Comp. Phy.

2016

, 318, 201.

(16) Bai, H.; Theuerkauf, J.; Gillis, P. A.; Witt, P. M. A Coupled DEM and CFD Simulation of Flow Field and Pressure Drop in Fixed Bed Reactor with Randomly Packed Catalyst Particles. Ind. Eng. Chem. Res..

2009

, 48, 4060.

(17) Donohue, T.; Wensrich, C. A numerical investigation of the void structure of brous materials. Powder Technol.

2008,

186, 72.

(18) Yagi, J.-I. Mathematical modeling of the ow of four uids in a packed bed. ISIJ Int. 1993

, 33, 619.

(19) Peters, B. Thermal conversion of solid fuels; WIT press, 2003. (20) Mahmoudi, A. H.; Homann, F.; Peters, B. Application of XDEM as a novel approach to predict drying of a packed bed. Int. J. of Therm. Sci. 36

ACS Paragon Plus Environment

2014,

75, 65.

Page 37 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(21) Janecki, D.; Burghardt, A.; Bartelmus, G. Inuence of the porosity prole and sets of Ergun constants on the main hydrodynamic parameters in the trickle-bed reactors.

Chem. Eng. J.

2014,

237, 176.

(22) Cundall, P. A.; Strack, O. D. L. A discrete numerical model for granular assemblies.

Geotechnique

1979,

29, 47.

(23) Peters, B.; Homann, F.; Senk, D.; Babich, A. Experimental and numerical investigation into iron ore reduction in packed beds . Chem. Eng. Sci.

2016,

140, 189.

(24) Michael, M.; Vogel, F.; Peters, B. DEM-FEM coupling simulations of the interactions between a tire tread and granular terrain. Computer Methods in Applied Mechanics and

Engineering

2015,

289, 227.

(25) Donoso, A. A. E.; Peters, B. {XDEM} employed to predict reduction of tungsten oxide in a dry hydrogen atmosphere . Int. J. Refract. Hard Met.

2015

, 49, 88.

(26) Mahmoudi, A. H.; Markovic, M.; Peters, B.; Brem, G. An experimental and numerical study of wood combustion in a xed bed using Euler-Lagrange approach (XDEM). Fuel 2015,

150, 573.

(27) Doiugys, A.; Peters, B. An approach to simulate the motion of spherical and nonspherical fuel particles in combustion chambers. Granular Matter

2001,

3, 231.

(28) Faghri, A.; Zhang, Y. Transport phenomena in multiphase systems; Academic Press, 2006. (29) Homann, F. Modeling heterogeneous reactions in packed beds and its application to the upper shaft of a blast furnace model heat-up. Ph.D. thesis, University of Luxem-

bourg, 2014.

37

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 40

(30) Peters, B.; Besseron, X.; Estupinan, A.; Mahmoudi, A.; Mohseni, M. A Discrete/Continuous Numerical Approach to Multi-physics. IFAC- PapersOnLine

2015,

48, 645. (31) Xiao, H.; Sun, J. Algorithms in a robust hybrid CFD-DEM solver for particle-laden ows. Com. in Comp. Phy.

2011,

9, 297.

(32) Gunjal, P. R.; Kashid, M. N.; Ranade, V. V.; Chaudhari, R. V. Hydrodynamics of Trickle-Bed Reactors: Experiments and CFD Modeling. Ind. Eng. Chem. Res.

2005,

44, 6278. (33) Samiei, K. Assessment of Implicit and Explicit Algorithms in Numerical Simulaion of Granular Material. Ph.D. thesis, University of Luxembourg , 2012. (34) Michael, M. A Discrete Approach to Describe the Kinematics between Snow and a Tire Tread. Ph.D. thesis, University of Luxembourg , 2014. (35) de Klerk, A.Voidage variation in packed beds at small column to particle diameter ratio. AIChE J.

2003

, 49, 2022.

(36) Mondal, S.; Wu, C.-H.; Sharma, M. M. Coupled CFD-DEM simulation of hydrodynamic bridging at constrictions . Int. J. of Multiphase Flow

2016,

84, 245.

(37) Miller, S.; Jasak, H.; Boger, D.; Paterson, E.; Nedungadi, A. A pressure-based, compressible, two-phase ow nite volume method for underwater explosions. Comp. and

Fluids

2013,

87, 132.

(38) Sun, C. G.; Yin, F. H.; Afacan, A.; Nandakumar, K.; Chuang, K. T. Modelling and Simulation of Flow Maldistribution in Random Packed Columns with Gas-Liquid Countercurrent Flow. Chem. Eng. Res. and Design

2000

, 78, 378.

(39) Hunt, M.; Tien, C. Non-Darcian ow, heat and mass transfer in catalytic packed-bed reactors. Chem. Eng. Sci.

1990,

45, 55. 38

ACS Paragon Plus Environment

Page 39 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(40) Goodling, J.; Vachon, R.; Stelpfug, W.; Ying, S.; Khader, M. Radial porosity distribution in cylindrical beds packed with spheres. Powder Techn.

1983,

35, 23.

(41) Mahmoudi, A. H.; Homann, F.; Peters, B.; Besseron, X. Numerical study of the inuence of particle size and packing on pyrolysis products using XDEM. Int. Commun.

Heat Mass Transfer

2016

, 71, 20.

(42) Mohseni, M. Drying and Conversion Analysis of Biomass by a DEM-CFD Coupling Approch. Ph.D. thesis, University of Luxembourg , 2017. (43) Li, L.; Ma, W. Experimental study on the eective particle diameter of a packed bed with non-spherical particles. Transp. Porous Media

2011,

89, 35.

(44) Iliuta, I.; Fourar, M.; Larachi, F. Hydrodynamic Model for Horizontal Two-phase Flow Through Porous Media. Can. J. Chem. Eng.

2003

, 81, 957.

(45) Narasimhan, C.; Verma, R.; Kundu, A.; Nigam, K. Modeling hydrodynamics of tricklebed reactors at high pressure. AIChE J.

2002

, 48, 2459.

(46) Abdolkarimi, V. Hydrodynamics modeling of a bench scale trickle bed reactor using CFD. Pet. Coal

2013,

55, 188.

39

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Graphical TOC Entry

40

ACS Paragon Plus Environment

Page 40 of 40