Coarse-Grained-Particle Method for Simulation of LiquidâSolids

Subscriber access provided by Northern Illinois University

Article

Coarse grained particle method for simulation of liquid-solids reacting flows Liqiang Lu, Kisoo Yoo, and Sofiane Benyahia Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.6b02688 • Publication Date (Web): 08 Sep 2016 Downloaded from http://pubs.acs.org on September 12, 2016

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 45

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

Coarse grained particle method for simulation of liquidsolids reacting flows Liqiang Lu, Kisoo Yoo, Sofiane Benyahia*

National Energy Technology Laboratory, Morgantown, WV 26507, USA

* Corresponding author: [email protected]

Abstract A coarse-grained particle method (CGPM) that tracks fewer particles and their collisions by lumping several real particles into computational parcels is developed. We demonstrate that the accuracy of CGPM can be improved by properly selecting interpolation schemes used to compute flow variables at particle and computational grid locations. The computational speed of CGPM was shown to be similar to that of a widely used and less accurate particle in cell method (PIC). Finally, a chemical reaction mechanism based on rare earth elements (REE) leaching from coal byproducts was implemented with CGPM to study the effects of several flow and design parameters on the leaching process. A counter-current reactor was studied and optimized to maximize the mass fraction of REE in the liquid solution. The CGPM can now be employed to cheaply and accurately solve industrial-scale problems containing millions of computational parcels. Key words: Computational fluid dynamics; discrete particle method; coarse grained particle methods; liquid- solids fluidized bed; rare-earth elements

ACS Paragon Plus 1Environment


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

1.

Page 2 of 45

Introduction

Computational fluid dynamics tools are widely used to improve our understanding of multi-phase flows and chemical reactor designs. The two fluid model (TFM)1 plays a dominent role in the simulation of engineering problems due to its relatively lower computation cost and easy access through commercial and open-source software. However, since the particles are considered to be a continuum phase, it is difficult to avoid numerical diffusion2, as well as accurately take into account the particle shape3, 4, particle-scale chemical reaction models5 and other types of particle-particle forces (such as cohesion)6. The computational fluid dynamic-discrete element method (CFD-DEM)7-9 avoids these complications by following the trajectory of each particle using Newton’s law of motion. However, the number of particles in real systems is usually much larger than current computational capabilities even by taking advantage of the intrinsic parallelism of processing the particle motion and particle-particle interactions and using multi-scale discrete supercomputing10 with modern accelerators such as graphics processing units (GPUs)11 and many integrated cores (MICs)12. To overcome the computation gap, many coarse grained particle methods (CGPM) were developed by tracking coarse grained particles (CGP) instead of real particles.13-19 The multiphase particle-in-cell method (MP-PIC) even abandoned the computation of collision forces in favor of a simpler continuum solids pressure model in order to increase computational speed.20, 21 The currently published CGPMs were mainly applied for non-reacting gas-solids flows. However, gas-solids flow systems show heterogeneous structures due to the large gas-solids density ratio22, 23 and so require correction to the homogeneous inter-phase drag model24-27 and restitution coefficient16, 17. Such corrections are not required in this study as it only involves liquid-solids flow systems that tend be more homogeneous due to smaller density ratio compared with gas-solids flow.28, 29 The calculation


Page 3 of 45

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


of fluid-particle interaction terms, such as drag force, can be significantly affected by the particle-togrid interpolation schemes. In traditional CFD-DEM, the computational fluid dynamics (CFD) grid is much larger than particle diameter and the simple center based method is commonly used. However, in CGPM, the diameter of coarse grained particle may be close to the size of CFD grid. In this case, the centroid method can introduce large statistical errors

30

. Thus, different interpolation schemes will be

discussed in this study. In this study, we propose a CGPM for liquid-solids reacting flows and show its computational speed to be closely comparable to that of a MP-PIC model. Instead of tracking CGPs whose physical properties are difficult to derive, the parcel in this method shares the same diameter and other physical properties as the real particles except for a statistical weight which represents the number of real particles in a parcel. Thus, the drag, gravity, and pressure forces are calculated at real particle scale. The collision force is, however, calculated using the collision diameter of the computation parcel. The proposed CGPM was verified by simulating a liquid-solids fluidized bed, and the results compared well with those obtained by both experiments and traditional CFD-DEM. This validated model was then applied to a hypothetical countercurrent fluidized bed reactor used for leaching valuable rare-earth elements (REE) from coal-byproducts. This study is detailed in the following sections: the physical model is addressed in Sec.2 and the simulation parameters are listed in Sec. 3. The simulation results are discussed in Sec. 4 and the main conclusions are drawn with prospects on future work in Sec.5.

2.

Description of the coarse grained particle method



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 4 of 45

Figure 1 shows that the main idea behind CGPM is the lumping of several real particles into a coarsegrained or computational parcel (CGP). In this figure, we take a computation cell in original system as an example to explain the assumptions and calculation method for drag force (b), chemical reactions (c) and collision forces (d). In Eulerian-Lagrangian method, the fluid phase is resolved at grid scale while the solids phase is resolved at particle scale. All the particles in a computation cell share the same (or an interpolated value) fluid velocity and fluid species fractions. Also the same drag force is applied to all the particles in a CGP due to their same collective motion (recall that the size of CGP is smaller than the cell-size). For cases involving chemical reactions, the reaction rate is influenced by species concentrations, temperature, particle size, and relative velocity. The chemical reaction rates are the same for all particles within a CGP as they share the same (or an interpolated value) fluid species fractions. Finally, the collision force is calculated using the collision diameter of the sphere representing the CGP. Then, the calculated force is divided equally on the particles constituting the CGP. The only difference with previously published CGPM is that instead of tracking CGPs, real particles with a statistic weight W (representing the number of particles per parcel) are tracked directly. It should be noted that these two approaches to CGPM are identical. The main motivation for using our current approach is to avoid modifying, and potentially introducing bugs to, numerical code for chemical reactions and heat transfer that was already written and validated for discrete particle methods. The detailed governing equations are given in the following parts.


Page 5 of 45

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 1. Schematic showing assumptions in CGPM (for explanation of equations see section 2.2)

2.1

Statistic weight and properties of the tracked particles

In CGPM, the number of real particles inside a CGP is represented by statistic weight W. It has a major influence on the simulation accuracy and directly related to computation cost. If W equals one, this method will converge to CFD-DEM. However, in most applications only the collective behavior of the particles, rather than the trajectories of all individual particles, are of primary interest. Statistically speaking, the CGPM is only valid when there is large number of particles in a system so that the collective behavior of these particles can be reproduced with fewer ones. In gas-solids flows, the collective motion of particles form heterogeneous structures called clusters, which reduce the drag force and have a significant influence on the flow field. The size of the clusters is influenced by the flow conditions and physical properties of fluid and particles. For this type of flows, the parameter W must be able to resolve the evolution (including deformation, aggregation and breakage) of clusters. Thus, the diameter of the coarse grained particles (CGPs) (dCGP) should be smaller than most cluster sizes, which can be determined by empirical corrections or using energy minimization multi-scale (EMMS) model as shown in the literature16. However, liquid-solids flows are more homogenous and clusters are rarely observed. We consider it only constrained by the geometry of the reactor. But, due to



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 45

the different methods used for fluid phase (Eulerian) and solid phases (Lagrangian), we need to map data from Eulerian grid to Lagrangian point and vice-verse. Thus, the diameter of the CGPs should be smaller than the CFD cell to avoid unphysical solid volume fraction in these cells. Also, the use of interpolation schemes can reduce the statistical errors of estimating the mean fields (such as void fraction) using Lagrangian data.26, 27 By lumping particles together, we lose information by reducing the frequency of particle collisions which directly influence the collision energy dissipation. As a quick remedy, some researchers suggested that the simulation accuracy can be improved by using a modified restitution coefficient in CGPM to ensure the energy dissipation due to particle collision is kept at the same level as the original systems.16, 17 This conclusion is drawn from gas-solids systems where particles speed and their collision frequency are high. But for liquid-solids systems, the particles speed and their collision frequency are much lower than in gas-solids systems. Thus, in our simulation the physical properties of CGPs are kept same as original real particles. The collision force is then calculated using the collision diameter which assuming all the solid mass in a CGP are lumped into a sphere. Then, the calculated force is divided equally into its member particles.

2.2

Particle governing equations

The CGPM is similar to other traditional CFD-DEM methods based on the Lagrangian tracking of particles. The conservation of linear momentum of each particle takes the familiar form of Newton’s second law of motion as:

mp

dv p dt

= Fb + Fgp + Fdrag + Fc


(1)

Page 7 of 45

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


Where mp is the mass of a real particle. Fb is the gravity force, Fgp is the pressure gradient force and Fdrag is the drag force. These forces are calculated following the same process as in traditional CFDDEM. As shown in Figure 2, the contact force is calculated by the soft-sphere DEM using dCGP as collision diameter, the force is then divided by statistic weight to represent the force on a real particle. Since our model was implemented in the MFIX-DEM code, the details of the DEM model can be found in the official MFIX- DEM documentation31. Here, we just make a simple description of this model as it applies to coarse-grained particles.

Figure 2. Collision force is calculated using diameter of CGP In DEM , the spring-dashpot model32 is used to model the particle-particle and particle-wall collisions. The normal contact force of particle i due to its collision with that of particle j is calculated as: •

Fijn = ( − k nδ n + υ n δ n )n ij

(2)

Where kn is the spring stiffness, δn is the thickness of their overlapped region and the dot representing time derivative, nij is a unit vector along the line of centers from j to i, and νn is the damping coefficient

υn = 2kn mCGP

ln(1/ eCGP )

π + [ln(1/ eCGP )]2 2

,


(3)


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 45

Where eCGP denotes the restitution coefficient for the CGPs. Similar to normal contact force, the tangential contact force can be calculated as: •

Fijt = ( − ktδ t + υt δ t )t ij

(4)

Where δt is the tangential displacement. The tangential force is also limited by the finite Coulomb friction between colliding particles with a sliding friction coefficient of µ. Thus, it is defined as:

Fijt = min(Fijt , µFijn )

(5)

Finally, the collision force on the tracked particle is:

Fc =

N

∑

(Fijn + Fijt ) / W

(6)

j =1, j ≠i

As for particle rotation, the moment of inertia I is calculated by mpdp2/10, it is W

5/3

times smaller

than the inertia based on coarse grained particle. Thus, the torque on each tracked particle should be calculated by:

T=

N

∑

( Ln × Fijt ) / W 5/3

(7)

j =1, j ≠i

Where L is the distance from particle center to contact point.

2.3

Fluid governing equations

As in traditional CFD, the fluid phase is described by the volume-averaged Navier-Stokes equations33, f ∂(ε f ρf ) + (∇ ⋅ ε f ρf uf ) = ∑ R fn ∂t n =1

N

D(ε f ρf u f ) = ∇ ⋅ S f + ε f ρf g − I Dt


(8)

(9)

Page 9 of 45

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


Where D/Dt is the substantial derivative, εf is the volume fraction of fluid, ρf is the density of fluid,

uf is the velocity of fluid, Nf is the number of fluid species, Rfn is the mass source term due to chemical reaction and I is the drag source term. S f is the fluid phase stress tensor given by

S f = − Pf I + τ f

(10)

Where Pf is the fluid phase pressure. Also, τ f is the fluid phase shear stress tensor defined as:

τ f = 2µf Df + λf ∇ ⋅ tr ( Df ) I

(11)

1 Df = [∇u f + (∇u f )T ] 2

(12)

Where D f is the strain rate tensor, and µf and λf are the dynamic and second coefficient of viscosity of the fluid phase. The mass balance for species i is

∂(ε f ρf X i ) + ∇ ⋅ (ε f ρf uf X i ) = ∇ ⋅ ( Di∇X i ) + S Xi ∂t

(13)

Where Xi is the mass fraction of fluid species i, Di is the diffusion coefficient and SXi is the source term due to chemical reaction.

2.4

Interphase drag forces

The drag source term on fluid cell k can be calculated as

1 I = vk k

β i 16 π d p3 W p ( v f (xi ) − v ip )K (xi , x k ) ∑ 1− εf i =1 Np

(14)

Where vk is the volume of cell k, Np is number of particles influence cell k, Wp is statistic weight of particle p, βι is drag coefficient of particle i in cell k, vf(xi) is the fluid velocity interpolated at particle i and K is the interpolation weight of particle i to cell k. In homogeneous liquid-solids flow, the Wen-Yu drag model is used and β is calculated as:



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

β=

3 ε f (1 − ε f ) ρf v f − v p CD0ε f −2.65 4 dp

Page 10 of 45

(15)

Where CD0 is the standard drag coefficient for a single particle

CD0

 24 0.687  Re (1 + 0.15Rep ) = p 0.44 

Rep < 1000

(16)

Rep ≥ 1000

Where

Rep =

ε f ρf v f − v p d p

(17)

µf

In an Eulerian-Lagrangian simulation, statistical errors are introduced into the system from the estimation of mean fields using Lagrangian data and vice-versa. Traditionally, the centroid method is used where the particle influences only the cell in which its center resides. When the particle diameter is close to cell size or the particle holds a large statistic weight, it will lead to a large statistic error. An analytic method to precisely calculate the particle volume in each cell could be used but is computationally expensive. Thus, approximate methods such as divide particle volume method (DPVM) and linear interpolation method (LIM) are tested in our simulation. The detailed implementation of these interpolation methods can be found in open source software MFIX31. 2.5

Chemical reaction and mass transfer models

With Eulerian-Lagrangian simulation, heterogeneous chemical reactions can be simulated with particle-scale chemical reaction kinetic models, such as shrinking core model, where the reaction rate is directly related to particle diameter. In this research, the leaching reactions of rare earth ore with magnesium sulfate are simulated. The chemical reaction mechanism is as follows34:

{ Al [ Si O ] ( OH ) } 4

4

10

−

8 m

⋅ 2 nRE(3s+) + 3nMg (2f+) → { Al4 [ Si4O10 ] ( OH )8 } ⋅ 3nMg (2f+) + 2 nRE(3s+) −

m

(18)

ACS Paragon Plus10 Environment

Page 11 of 45

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


Where s and f are solid phase and fluid phase, respectively. It is a typical heterogeneous lumped reaction in liquid-solids systems and the inner diffusion controlled shrinking core model was demonstrated to best describe the reaction mechanism by the following equations 34

2 2/3 1 − α − (1 − α ) = kt 3 r  α = 1−  c  R   p kexp (1/ min ) =

(19)

3

0.011 R p ( mm )

1.217

(20)

 9480  exp  −   RT 

(21)

Where α is the leached fraction, kexp is the chemical reaction rate from experimental measurements, rc is the radius of the shrinking core, Rp is the radius of particle and R is the universal gas constant, T is

reaction absolute temperature. In this model, the influence of Mg2+ concentration need to be considered and the reaction rate is calculated as: k (1 / min ) =

CMg CMg ,exp

kexp

(22)

Where CMg,exp is Mg2+ concentration in experiment. Since the experiment is carried in a fast stirred 500 mL container with enough MgSO4 solution, the system is considered perfectly mixed and the Mg2+ concentration is assumed constant. This reaction kinetic has already been implemented in MFIX and verified by zero-dimensional (0D) simulations.35 The particle’s mass variation due to chemical reaction is given by36:

dmp dt

Ns

= ∑ Rsn

(23)

n =1

Where Rsn is the rate of production (or consumption) of species n in this particle. The mass fraction of species n (Xpn) in this particle can be expressed as:



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

Ns 1 = ( Rsn − X pn ∑ Rsn ) dt mp n =1

dX pn

Page 12 of 45

(24)

The mass transfer of species n between particles and fluid phase can be calculated as:

R fn = −

1 vk

Nk

∑W R i

i fn

i

K (x , x k )

(25)

i =1

Where vk is the volume of cell k, Nk is number of particles influencing cell k, Wi is statistic weight of particle i, Ri fn is the rate of fluid phase species production (or consumption) due to a particle-fluid reaction between particle i and the fluid phase, and K is the interpolation weight of particle i to cell k. Since the leaching process was carried at a constant temperature, the heat transfer model is not included in our simulation. However, heat transfer model was implemented and is available in MFIX CGPM.

3. Simulation parameters

To validate the method proposed in Sect. 2, a water-glass beads fluidized bed37 was simulated. The geometry is shown in Figure 3, the fluid flows into the reactor from the bottom uniformly with a constant velocity of 7.62 cm/s, and at the outlet of the reactor the fluid gage pressure is set to a constant value of zero. The physical properties and simulation parameters are listed in Table 1.


Page 13 of 45

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 3. Geometry and boundary conditions of the simulated water and glass beads fluidized bed



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 14 of 45

Table 1. Simulation parameters of water and glass beads fluidized bed

Particle

Liquid

Bed

Cases Number of CGP, NCGP Statistic weight, W dCFD/dCGP Interpolation scheme Diameter, dp (cm) Density, ρp (g/cm3) Spring constant, kn (dyne/cm) Restitution coefficient, e Friction coefficient, µ Time step, ∆t (s) Density, ρl (g/cm3) Viscosity, µ l (g/cm/s) Superficial velocity, Uin (cm/s) Time step, ∆t (s) Bed diameter, D (cm) Bed height, H (cm) CFD Grid, dx x dy x dz (cm3)

CFD-DEM 40 960 1.0 4.425 Centroid\DPVM

CGPM-8 5 120 8.0 2.212 Centroid\ DPVM\LIM 0.113 1.0

CGPM-27 1 517 27.0 1.475 Centroid\ DPVM\LIM

1.0 x 105 0.7 0.1 1.0 x 10-4 1.0 0.01236 7.62 Varied 4.0 18.0 0.5 x 0.5 x 0.5

To compare the computation speed of CGPM with MP-PIC, four different scales of a hypothetical liquid-solids system is simulated with different number of computation particles. The simulation parameters are listed in Table 2.

ACS Paragon Plus Environment

14

Page 15 of 45

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


Table 2. Simulation parameters of scaled liquid-solids fluidized bed

Bed

Particle

Fluid

Cases Bed size, Lx x Ly x Lz (cm3) CFD Grid number, Nx x Ny x Nz Number of computation particles, NC Statistic weight, W Computation cores Parallel method Diameter, dp (cm) Density, ρp (g/cm3) Spring constant, kn (dyne/cm) Restitution coefficient, e Friction coefficient, µ MP-PIC time step, ∆tMPPIC (s) CGPM time step, , ∆tCGPM (s) Density, ρf (g/cm3) Viscosity, µf (g/cm/s) Superficial velocity, Uin (cm/s) Time step, ∆tMP-PIC (s)

Scale-A

Scale-B

Scale-C 19.2 x 38.4 x 19.2 40 x 80 x 40

Scale-D 38.4 x 76.8 x 38.4

4.8 x 9.6 x 4.8

4.8 x 9.6 x 4.8

20 x 40 x 20

10 x 20 x 10

160 000

20 000

1 280 000

10 240 000

8 16 OpenMP

64 16 OpenMP

64 16 OpenMP

64 32,64,128,256 MPI

80 x 160 x 80

0.04 1.34 1.0x104 0.7 0.1 1.0 x 10-3 1.0 x 10-4 1.0 0.01 0.4534 1.0 x 10-3

To validate the chemical reaction model in this method, the rare earth elements (REE) extraction process is simulated in a small bubbling fluidized bed with dimensions 8 cm x 32 cm x 8 cm and uniform CFD grid with size of 1 cm x 1 cm x 1cm. The physical properties of particles and fluid are listed in Table 3. In the experiment34, the particle size is in the range of 0.045 cm to 0.060 cm and 0.060 cm to 0.090 cm. Since the detailed particle size distribution is not given, a single particle size is used in each simulation.


15


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 16 of 45

Table 3. Simulation parameters of fluidized leaching bed Properties

value

Diameter, dp (cm)

0.045, 0.06, 0.09

Density, ρp (g/cm3)

1.34

Particle components mass fraction: REE, NH4, 0.00102, 0.0, 0.99898 Ash Particle

Liquid

Spring constant, kn (dyne/cm)

1.0x105

Restitution coefficient, e

0.7

Friction coefficient, µ

0.1

Time step, ∆t (s)

Varied

Density, ρl (g/cm3)

1.0

Viscosity, µ l (g/cm/s)

0.01

Initial liquid components mass fraction: REE, 0.0,0.002730,0.007346,0.989924 NH4, SO4, H2O Time step, ∆t (s)

Varied

According to a recent study38, column leaching provides a more efficient alternative to the batch process and constitutes an important step towards simulating the heap or in situ leaching. However, the capability of this system is limited by liquid flow rate through the densely packed bed. Kwauk

39, 40

investigated leaching process in a counter-current fluidized bed. The particles

are added from the top of the reactor and fall counter-currently against a rising stream of liquid into the fluidized leaching section. Fresh leaching liquid is enters the system through special sprayers located near the bottom of the reactor. Below the liquid inlet region, the solids are discharged from the bottom with an automatically controlled valve. In our simulation, the inlet flow is simulated with point source and the closing and opening of the valve is controlled by the total number of particles in the reactor. A counter-current fluidized bed capable of leaching REE


16

Page 17 of 45

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


shown in Figure 4 was simulated using CGPM with different operating conditions and reactor dimensions. The geometry and operating conditions are listed in Table 4. Other physical parameters are listed in Table 3.

Figure 4. Schematic of counter-current fluidized bed for leaching process Table 4. Geometry and operating conditions used in counter-current leaching-process simulations Case

Recycle50

Recycle80

RT100

SIH58

L36

L36S36

Liquid inlet flow 72 rate, L (kg/h)

72

72

72

72

36

36

Solid inlet rate, S (kg/h)

18.41

18.41

18.41

18.41

18.41

36.82

50

80

20

20

20

20

Liquid

Base

flow

18.41

recycle 20


17


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 18 of 45

fraction, (%) Residence time, (s)

80

80

80

100

100

80

80

Diameter, D1 (cm)

8

8

8

8

8

8

8

Diameter, D2 (cm)

4

4

4

4

4

4

4

Height, H1 (cm)

4

4

4

4

4

4

4

Height, H2 (cm)

8

8

8

8

8

8

8

Height, H3 (cm)

10

10

10

10

10

10

10

Height, H4 (cm)

48

48

48

48

58

48

48

Height, H5 (cm)

64

64

64

64

64

64

64

4.

Results and discussion

In this section, the simulated overall bed voidage (liquid volume fraction) is quantitatively compared with the experimental result. And the influence of statistic weight and interpolation schemes on simulation results are also discussed. After that, the computation speed of CGPM is directly compared with MP-PIC for different cases. Then, the chemical reaction model is verified through the simulation of REE leaching. At last, a countercurrent fluidized bed used for leaching valuable REE from coal-byproducts is simulated. The influence of reactor dimensions and operating conditions are discussed.

4.1

Model verification

Figure 5 shows instantaneous particle position and axial velocity distributions simulated using the CFD-DEM with centroid interpolation scheme. Initially, the particles are uniformly distributed at the bottom of the reactor. As the liquid flows into the reactor through the bottom inlet, the total bed expands uniformly reaching a steady height at about 6.0 s. A homogeneous


18

Page 19 of 45

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


flow is then established as evident by the lack of clearly defined bubbles as in other gas-fluidized systems.

Figure 5. Instantaneous particle position and velocity distributions using CFD-DEM with centroid interpolation method Figure 6 shows the instantaneous particle distribution from CFD-DEM and CGPM with centroid interpolation scheme at time 8.0 s. The total bed height of CFD-DEM and CGPM-8 (using a statistical weight W = 8) is similar while CGPM-27 (W = 27) is much higher. The particle distributions of interpolation schemes such as divide particle volume method (DPVM) and linear interpolation method (LIM) are also shown in Figure 7. These results will be qualitatively analyzed and explained in details shortly.


19


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 45

Figure 6. Particle distribution of CFD-DEM and CGPM with centroid interpolation at time 8.0 s


20

Page 21 of 45

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 7. Particle distribution of CFD-DEM and CGPM with DPVM and LIM interpolation at time 8.0 s

As shown in Figure 8, the simulated overall bed voidages using CFD-DEM and CGPM with different number of statistic weights are compared with experimental results. The centroid method is used in CFD-DEM and CGPM to compute interphase exchange terms. Failing to interpolate particle data and flow fields, such as void fraction, is shown in this figure to produce larger errors at larger statistical weight.

The results from CFD-DEM compared well with

experimental data with an error of about 1%. For case CGPM-8, the difference from experimental data is still small (2.0%) because a small statistic weight is used and the flow field of liquid-solids flow is homogeneous as shown in Figure 5. But for CGPM-27, the error increased to about 6% (Note that for this case some particles flow out of the simulation domain) because of a large statistic weight is used. For this case, the ratio of CFD grid size dCFD to particle size dCGP is only about 1.475, which introduces a large statistic error when centroid method is used. According to a previous study30, the critical ratio is about 3.82 for gas solid bubbling fluidized bed. Here, in liquid-solids flow, the critical ratio can be as small as 2.2 (case CGPM-8). When the ratio of CFD grid size to particle size is smaller than this value, more precise interpolation schemes should be used to reduce the statistic error, as will be shown next.


21


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 45

Figure 8. Overall bed voidages from experiment and CFD-DEM, CGPM simulations As discussed before, the CGPM allows a significant reduction of simulated particles, however, more accurate methods to interpolate particles and fluid data is necessary. Figure 9 compares the results with different statistic weight using SQUARE_DPVM and LIM interpolation schemes implemented in MFIX. It is clear that a more precise interpolation scheme improves the simulation results as the error for case CGPM-27 is reduced from about 6% to only 0.15%. The LIM interpolation scheme can also improve the results similarly to SQUARE_DPVM interpolation scheme.

Figure 9. Overall bed voidages with DPVM and LIM interpolation schemes


22

Page 23 of 45

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


4.2

Comparison of CGPM and MP-PIC Computational Speed

The CGPM enables the reduction of computational particles but it still needs the calculation of particle-particle collision forces which is time consuming due to the additional time spent in collision pair searching method. MP-PIC, on the other hand, does not compute particle-particle collision forces directly, and instead, a continuum particle pressure term is introduced to avoid unphysical large solids volume fractions in computation cells. Thus, the solids time step is no longer limited by collision process and a much larger solid time step can be used. So, the computation speed of MP-PIC can be faster than CGPM. However, the collisional algorithm in CGPM can benefit from modern parallel supercomputers using openMP and MPI as well as some recent development in hardware such as GPUs and MICs. In this section, the computational speed of CGPM and MP-PIC both implemented in MFIX are compared. In order to make a fair comparison, the time steps for liquid phase of both MP-PIC and CGPM are set to the range of 1.0 x 10-3 to 1.0 x 10-6 (note that all times and time-steps are in seconds). After few initial time-steps, both are automatically adjusted to 1.0 x 10-3 by the computation program. The solids phase time step of MP-PIC is also automatically adjusted to 1.0 x 10-3, as the maximum solids time-step cannot exceed that of the fluid, while the time step of CGPM is fixed at 1.0 x 10-4. Other parameters are listed in Table 2. The instantaneous particle distribution of MP-PIC and CGPM for case Scale-A (see Table 2) are compared in Figure 10 (a). The final bed height of CGPM is about 6.0 cm while MP-PIC is about 5.9 cm. Thus, the bed height computed by CGPM is about 1.7% smaller than MP-PIC. The radial fluid volume fraction and solid velocity distribution at 2.04 cm is also compared in Figure 10 (b). The results from CGPM and MP-PIC are similar due to the uniform nature of the simulated liquid-particle flows. The advantage of CGPM is that inter-parcel collisions are


23


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 45

directly resolved. In the MP-PIC model the collisions are not directly resolved. It is based on mean fields (such as solid volume fraction) and their spatial gradients, which are directly computed from the parcels location. However, the estimation of mean fields from Lagrangian data is fraught with noise stemming from statistical errors caused by finite sample size. As a result, it is not trivial to extend the closures used in two-fluid theory (based on mean fields and their spatial gradients) to model collisions in the MP-PIC model. Another advantage of CGPM is the fact that it will converge to CFD-DEM as W tends toward a value of one. This convergence is not guaranteed in MP-PIC.

(a) Instantaneous particle distribution


24

Page 25 of 45

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


(b) Radial fluid volume fraction and velocity distribution at 2.04 cm (averaged from 30 s to 40 s)

Figure 10. Simulation results of MP-PIC and CGPM for case Scale-A Figure 11 compares the computation speed of MP-PIC and CGPM for different number of computation particles (from Scale-A to Scale-C). Noting that all the cases are computed using 16 CPU cores on a single machine with openMP multi-threads. For this case, the computation speed of MP-PIC is only about 10%-20% faster than CGPM. This is unexpected as the MP-PIC technique does not compute particle-particle collisions or conduct neighbor search, and uses a larger solids time-step. If this result is reproduced by other cases in future research, it will be of great significance because the CGPM is more accurate than the MP-PIC technique as it relies on fewer assumptions. Unlike MP-PIC, it is straightforward to show the CGPM will converge to the more accurate CFD-DEM as the parcel diameter is reduced. Therefore, we may consider the CGPM as a potential candidate for solving large-scale industrial fluid-solids problems even if its computational speed is less than that of MP-PIC as large computer systems are becoming more ubiquitous.

Figure 11. Computation speed of MP-PIC and CGPM for case Scale-A


25


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 45

For CGPM, the computation can also be accelerated further by using distributed memory parallel computation with MPI. Figure 12 compares the computation speed of CGPM for case Scale-D (with 10.24 million computation particles) with different number of MPI processes. It shows that the computation speed for 10.24 million computation particles is about 0.58 physical seconds per hour (the corresponding total strong scale efficiency of 256 processes is about 40%). As mentioned earlier, the CGPM offers us an opportunity to more accurately investigate large industrial reactors.

Figure 12. Computation speed of CGPM for case Scale-D

4.3

Simulation of REE leaching in a liquid-fluidized bed

The proposed chemical reaction model in CGPM ensures the same result as with CFD-DEM if the particle-scale reaction rate is not related to local voidage or fluid field variables. If the particle scale reaction rate is analytically related to local voidage or other fluid field variables, the result will depend on the ability of CGPM to reproduce the local hydrodynamics which accuracy has already been tested in Sec. 4.1.


26

Page 27 of 45

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


To validate the implementation of chemical reaction kinetics and mass transfer in CGPM, an REE extraction process is simulated in a small bubbling fluidized bed. Figure 13 compares the REE leaching fraction between experiment and CGPMs with different statistic weight and particle diameters. Although the experimental data is from a continuously stirred three necked flask, the simulation results compare well with experimental data due to the fact that a fluidized bed is a fully mixed reactor, especially for liquid-solids flow. The results are also insensitive to statistic weight.

Figure 13. Comparison of REE leaching fraction obtained experimentally and with CGPM Figure 14 shows the instantaneous particle distribution colored by particle velocity in the vertical direction. It shows the fluidization behavior is generally uniform, i.e. not bubbling vigorously as typically seen in gas-fluidized systems (according to Geldart diagram41, the simulated particle belong to Geldart B). Similar to our previous findings in Sec. 4.1, the simulated hydrodynamics obtained with different statistic weights show similar results. Hence, these observations prove that the proposed chemical reaction model in CGPM is sufficiently accurate to further investigate other REE leaching processes of interest to our research.


27


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 28 of 45

Time = 1.0 s

Time = 10.0 s

Figure 14. Behavior of a liquid-fluidized bed of 0.045 cm diameter particles simulated with different statistical weights (left W=1000, right W=500) 4.4

Simulation results of REE leaching in a counter current fluidized bed


28

Page 29 of 45

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


To investigate the influence of operating conditions and reactor geometry, the results of several cases are analyzed in this section. For the base-case, the liquid inlet flow rate L is 72kg/h, the solids inlet flow rate S is 18.41 kg/h, the liquid recycle fraction RF is 20%, the residence time of particles in the reactor RT is 80 s and the height of solid inlet HIS is 48 cm. The influence of these design and operation parameters on REE leaching rate are investigated, and the values of these parameters are listed in Table 4. It is noteworthy to mention that several three-dimensional (3D) transient simulations were done for this REE leaching multiphase flow reactor. It became possible to conduct thousands of operating seconds because of the availability of larger and faster computers, and more importantly due to simplifying assumptions such as the use of computational parcels as in CGPM. Therefore, this final section demonstrates the ability of CFD to be used as a process simulation tool by conducting parametric flow studies in a reasonable time.

4.4.1

Simulation results of the base-case

For continuous reactors, the stability of operation is of critical importance to the safety and product quality. Thus, the inventory of solids mass in the reactor and liquid mass flow rate at top and bottom outlet are investigated. As shown in Figure 15, the particles are added into the reactor from time 0 s and reach a statistically-steady state after about 80 s. The small fluctuations, caused by the outflow of particles from the bottom, is an indication of a stable operation. Figure 16 shows transient variations of the liquid mass flow rate at top and bottom outlet. At the

beginning of the operation, no particles flow out from the bottom outlet, and the flow rates at both top and bottom outlets are constant and sum up exactly to the inlet flow rate of 20 g/s. After around 80 s, small fluctuations are also observed due to the outflow of particles. But generally, these fluctuations are small compared to the total inlet flow rate. This indicates that the reactor


29


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 30 of 45

can be operated smoothly with such liquid and solids boundary conditions. Finally, it should be indicated that the results obtained with different statistical weights are almost identical demonstrating again the applicability of CGPM to these reacting flow systems.

Figure 15. Inventory of solids in the reactor


30

Page 31 of 45

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 16. Liquid mass flow rate at top and bottom outlet For a similar counter-current system, Kwauk 39, 40 observed, a dilute phase region surmounting a dense phase region with a somewhat well-defined interface between the two. Figure 17 shows the time-averaged liquid volume fraction (also called voidage in this study) distribution in the reactor. It is clear that there is a dense phase surrounding the liquid inlet and above the dense phase is the dilute phase. The interface is well defined similar to experimental observations. The REE mass fraction gradient is well established in both dilute region and dense region as shown in Figure 18.


31


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 45

Figure 17. Time-averaged voidage distribution in the reactor. The line-plots are radiallyaveraged values of the colored surface-plots.


32

Page 33 of 45

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 18. REE mass fraction in liquid with zoomed bottom region of reactor detailing the solids flow field The counter-current interactions of solids and liquid result in a relatively complex flow field in the reactor than previously investigated co-current fluidized bed. As shown in Figure 19, the average superficial fluid velocity should be about 0.4 cm/s with a liquid flow rate of 72kg/h. At the dense phase near the bottom liquid inlet, since the fluid volume fraction is about 0.5 to 0.6 as shown in Figure 17 , thus, the actual fluid velocity can be larger than 0.8 cm/s. There are some downward flows below both the liquid inlet and the solids inlet area due to the feeding and discharge of solid particles.


33


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 34 of 45

Figure 19. Detailed liquid flow field in the reactor To investigate the product quality, the REE mass fraction in the liquid-phase at top and bottom outlets is shown in Figure 20. In this simulation, 20% of the total outlet liquid is mixed with original liquid and recycled into the reactor. There is significant increase in REE leaching rate at first 200 s, then after about 800 s the REE mass fraction at outlet reaches 0.028% and stabilized.


34

Page 35 of 45

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


Because of recycled liquid, the REE leaching fraction at outlet is reduced from 80% to 75% at first 200 s and then stabilized as shown in Figure 21.

Figure 20. REE mass fraction at top outlet, bottom outlet and mixed inlet


35


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 45

Figure 21. REE leaching fraction based on solids REE content 4.4.2

Influence of operating parameters and reactor geometry

In the base-case, the REE mass fraction in the liquid solution at the outlet is only about 0.028%. Such a small value will waste a lot of energy in the downstream concentration process, and thus, it is important to increase it during the leaching process. It is obvious that the REE mass fraction at the outlet can be increased by increasing the fraction of recycled fluid. Figure 22 compares the REE mass fraction at the top outlet and bottom outlet with different liquid recycle fractions. It is clear that the REE mass fraction is increased from 0.028% to 0.11% when the recycle fraction is increased from 20% to 80%. Theoretically, the outlet REE mass fraction can be increased 4 times if the recycle fraction is increased 4 times. The computed mass fraction is very close to this value, albeit a little smaller due to a slight decrease of reactant concentration as indicated by the decrease in leaching fraction.

Figure 22. REE mass fraction at top outlet and bottom outlet with different liquid recycle fraction


36

Page 37 of 45

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 REE mass fraction in the liquid solution at the outlet can also be increased if we can increase the leaching fraction in the solids. To do so, the solids residence time is increased from 80 s to 100 s in case RT100 (by increasing the number of particles in reactor when valve should open) and the solids inlet position is adjusted from a height of 48 cm to 56 cm in case SIH56. Figure 23 compares the REE leaching fraction at the bottom outlet with different residence time and solids inlet position. It is noticeable that the leaching fraction is increased from 75% to 81.4% as the residence time is increased from 80 s to 100 s. Also, the leaching fraction was increased further by elevating the solids inlet position. With this increased leaching fraction, the REE mass fraction in the liquid solution at the outlet is also increased as shown in Figure 24.

Figure 23. REE leaching fraction variation with increased residence time and solid inlet position


37


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 45

Figure 24. Temporal REE mass fraction variations at top and bottom outlets with increased residence time and elevated solids inlet position The REE mass fraction can further be increased by decreasing the liquid to solids ratio. In case L36, the liquid flow rate is decreased two times to 36 kg/h. Thus, the liquid to solids ratio is reduced from 3.910 to 1.955. In case L36S36, the solids flow rate is increased two times to 36.82 kg/h and the liquid-solids ratio is reduced from 1.955 to 0.978. As shown in Figure 25, the REE mass fraction in liquid solution at the outlet is increased to 0.06% in case L36 and 0.12% in case L36S36.


38

Page 39 of 45

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 25. REE mass fraction at top outlet and bottom outlet with reduced liquid-solids ratio

5.

Conclusions

A coarse-grained particle model (CGPM) was developed for simulating liquid-solids reacting flows. In this method, the drag, gravity, and pressure forces are calculated at real particle scale while the collisional forces are calculated using the collision diameter of the computational particle (also called parcel). A parcel lumps several particles, which greatly reduces the computational effort of a CGPM simulation. For chemical reactions, the particle scale inner diffusion-controlled shrinking core model is adopted to describe reactant migration at real particle scale. The proposed method was verified by simulating a liquid-solids fluidized bed, and the results compared well with both experiment and traditional CFD-DEM. The computational speed of CGPM was similar to that of the less accurate MP-PIC approach. In fact, the computation speed of CGPM was only about 15% slower than MP-PIC in MFIX software with parcel numbers varying from thousands to millions. This makes the CGPM a valuable tool to tackle industrial-


39


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 40 of 45

scale problems because its accuracy and speed of execution can be adjusted by simply adjusting the number of real particles contained in a computational parcel. This verified model was then applied to a counter-current fluidized bed reactor used for leaching valuable rare-earth elements (REE) from coal-byproducts. The simulation results show that the REE concentration in liquid solution can be increased by recycling leachate, increasing solid residence time and reducing liquid to solids mass ratio. The scale-up of this liquid-solids reactor for REE leaching as well as other reactor designs will be pursued in future research.

Acknowledgment This research was supported in part by an appointment to the National Energy Technology Laboratory Research Participation Program, sponsored by the U.S. Department of Energy and administered by the Oak Ridge Institute for Science and Education.

Notation dp

diameter of particle (m)

dCGP

diameter of CGP (m)

dx

CFD grid size in x direction (m)

dy

CFD grid size in y direction (m)

dz

CFD grid size in z direction (m)

dCFD

smallest CFD grid size in x, y and z direction (m)

CD0

standard drag coefficient for a CGP

CMg,exp

Mg2+ concentration in experiment

D

Reactor diameter (m)


40

Page 41 of 45

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


eCGP

restitution coefficient of CGP

g

gravity (m/s2)

Fb

body force (N)

Fgp

pressure gradient force (N)

Fc

contact force (N)

Fdrag

drag force (N)

Fijn

normal contact force (N)

Fijs

tangential contact force (N)

H

Reactor height (m)

kn

normal spring stiffness (N/m)

kt

tangential spring stiffness (N/m)

kexp

chemical reaction rate from experiment (1/min)

K

interpolation weight

mp

mass of particle (kg)

Nf

number of fluid species

Np

number of particles

NCGP

number of coarse grained particles

Nx

number of CFD grid size in x direction

Ny

number of CFD grid size in y direction

Nz

number of CFD grid size in z direction

nij

normal unit vector between the center of particle i and particle j

rc

radius of reaction core (m)

R

universal gas constant (J/mole/K)


41


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

Rp

radius of particle (m)

sij

tangential unit vector between CGP i and CGP j

T

reaction absolute temperature (K)

uf

velocity of fluid (m/s)

Uin

superficial inlet velocity of fluid (m/s)

vs

velocity of particle (m/s)

vk

volume of cell k (m3)

Vsys

volume of simulation domain (m3)

W

statistic weight

Xi

mass fraction of fluid species i

Page 42 of 45

Greek letters β

drag coefficient with structure in a control volume (kg/m3 s)

εf

voidage of the dilute phase

δn

thickness of the overlapped region (m)

δs

tangential displacement (m)

νn

normal damping coefficient

νs

tangential damping coefficient

µ

friction coefficient

µg

gas viscosity (Pa s)

∆t

time step (s)

Subscripts


42

Page 43 of 45

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


f

fluid phase

s

solid phase

p

particle

CGP

coarse-grained particle

References

(1) Anderson, T. B.; Jackson, R., Fluid Mechanical Description of Fluidized Beds. Equations of Motion. Ind. Eng. Chem. Fundam. 1967, 6, 527-539. (2) Guenther, C.; Syamlal, M., The effect of numerical diffusion on simulation of isolated bubbles in a gas–solid fluidized bed. Powder Technol. 2001, 116, 142-154. (3) Peirano, E.; Delloume, V.; Leckner, B., Two- or three-dimensional simulations of turbulent gas–solid flows applied to fluidization. Chem. Eng. Sci. 2001, 56, 4787-4799. (4) Hua, L.; Zhao, H.; Li, J.; Wang, J.; Zhu, Q., Eulerian–Eulerian simulation of irregular particles in dense gas–solid fluidized beds. Powder Technol. 2015, 284, 299-311. (5) Blamey, J.; Zhao, M.; Manovic, V.; Anthony, E. J.; Dugwell, D. R.; Fennell, P. S., A shrinking core model for steam hydration of CaO-based sorbents cycled for CO2 capture. Chem. Eng. J. 2016, 291, 298-305. (6) Gu, Y.; Ozel, A.; Sundaresan, S., A modified cohesion model for CFD–DEM simulations of fluidization. Powder Technol. 2016, 296, 17-28. (7) Tsuji, Y.; Kawaguchi, T.; Tanaka, T., Discrete particle simulation of two-dimensional fluidized bed. Powder Technol. 1993, 77, 79-87. (8) Hoomans, B. P. B.; Kuipers, J. A. M.; Briels, W. J.; van Swaaij, W. P. M., Discrete particle simulation of bubble and slug formation in a two-dimensional gas-fluidised bed: A hardsphere approach. Chem. Eng. Sci. 1996, 51, 99-118. (9) Zhong, W.; Yu, A.; Zhou, G.; Xie, J.; Zhang, H., CFD simulation of dense particulate reaction system: Approaches, recent advances and applications. Chem. Eng. Sci. 2016, 140, 1643. (10) Ge, W.; Lu, L.; Liu, S.; Xu, J.; Chen, F.; Li, J., Multiscale Discrete Supercomputing – A Game Changer for Process Simulation? Chemical Engineering & Technology 2015, 38, 575-584. (11) Xu, J.; Qi, H.; Fang, X.; Lu, L.; Ge, W.; Wang, X.; Xu, M.; Chen, F.; He, X.; Li, J., Quasi-real-time simulation of rotating drum using discrete element method with parallel GPU computing. Particuology 2011, 9, 446-450. (12) Li, B.; Zhou, G.; Ge, W.; Wang, L.; Wang, X.; Guo, L.; Li, J., A multi-scale architecture for multi-scale simulation and its application to gas–solid flows. Particuology 2014, 15, 160-169.


43


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 44 of 45

(13) Nasato, D. S.; Goniva, C.; Pirker, S.; Kloss, C., Coarse Graining for Large-scale DEM Simulations of Particle Flow – An Investigation on Contact and Cohesion Models. Procedia Engineering 2015, 102, 1484-1490. (14) Sakai, M.; Koshizuka, S., Large-scale discrete element modeling in pneumatic conveying. Chem. Eng. Sci. 2009, 64, 533-539. (15) Hilton, J. E.; Cleary, P. W., Comparison of non-cohesive resolved and coarse grain DEM models for gas flow through particle beds. Appl. Math. Model. 2014, 38, 4197-4214. (16) Lu, L.; Xu, J.; Ge, W.; Yue, Y.; Liu, X.; Li, J., EMMS-based discrete particle method (EMMS–DPM) for simulation of gas–solid flows. Chem. Eng. Sci. 2014, 120, 67-87. (17) Benyahia, S.; Galvin, J. E., Estimation of Numerical Errors Related to Some Basic Assumptions in Discrete Particle Methods. Ind. Eng. Chem. Res. 2010, 49, 10588-10605. (18) Chu, K.; Chen, J.; Yu, A., Applicability of a coarse-grained CFD–DEM model on dense medium cyclone. Miner. Eng. 2016, 90, 43-54. (19) Lu, L.; Xu, J.; Ge, W.; Gao, G.; Jiang, Y.; Zhao, M.; Liu, X.; Li, J., Computer virtual experiment on fluidized beds using a coarse-grained discrete particle method—EMMS-DPM. Chem. Eng. Sci. 2016, 155, 314-337. (20) Andrews, M. J.; ORourke, P. J., The multiphase particle-in-cell (MP-PIC) method for dense particulate flows. Int. J. Multiphase Flow 1996, 22, 379-402. (21) Snider, D. M., An Incompressible Three-Dimensional Multiphase Particle-in-Cell Model for Dense Particle Flows. J. Comput. Phys. 2001, 170, 523-549. (22) Ge, W.; Li, J., Physical mapping of fluidization regimes—the EMMS approach. Chem. Eng. Sci. 2002, 57, 3993-4004. (23) Lim, K. S.; Zhu, J. X.; Grace, J. R., Hydrodynamics of gas-solid fluidization. Int. J. Multiphase Flow 1995, 21, Supplement, 141-193. (24) Benyahia, S.; Sundaresan, S., Do we need sub-grid scale corrections for both continuum and discrete gas-particle flow models? Powder Technol. 2012, 220, 2-6. (25) Yang, N.; Wang, W.; Ge, W.; Li, J., CFD simulation of concurrent-up gas–solid flow in circulating fluidized beds with structure-dependent drag coefficient. Chem. Eng. J. 2003, 96, 7180. (26) Igci, Y.; Andrews, A. T.; Sundaresan, S.; Pannala, S.; O'Brien, T., Filtered two-fluid models for fluidized gas-particle suspensions. AlChE J. 2008, 54, 1431-1448. (27) Benyahia, S., On the Effect of Subgrid Drag Closures. Ind. Eng. Chem. Res. 2010, 49, 5122-5131. (28) Zheng, Y.; Zhu, J.-X.; Marwaha, N. S.; Bassi, A. S., Radial solids flow structure in a liquid–solids circulating fluidized bed. Chem. Eng. J. 2002, 88, 141-150. (29) Richardson, J. F.; Zaki, W. N., Sedimentation and Fluidization. Trans. Instn. Chem. Engrs. 1954, 35-53. (30) Peng, Z.; Doroodchi, E.; Luo, C.; Moghtaderi, B., Influence of void fraction calculation on fidelity of CFD-DEM simulation of gas-solid bubbling fluidized beds. AlChE J. 2014, 60, 2000-2018. (31) Garg, R.; Galvin, J.; Li, T.; Pannala, S. Documentation of open-source MFIX-DEM software for gas-solids flows. (32) Cundall, P. A.; Strack, O. D. L., A discrete numerical model for granular assemblies. Géotechnique 1979, 29, 47-65. (33) Syamlal, M.; Rogers, W.; O’Brien, T. J., MFIX documentation: Theory guide. National Energy Technology Laboratory, Department of Energy, 1993.


44

Page 45 of 45

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


(34) Xiao, Y.; Chen, Y.; Feng, Z.; Huang, X.; Huang, L.; Long, Z.; Cui, D., Leaching characteristics of ion-adsorption type rare earths ore with magnesium sulfate. Transactions of Nonferrous Metals Society of China 2015, 25, 3784-3790. (35) Yoo, K.; Lu, L.; Benyahia, S., Modeling and simulation of rare earth element extraction by ion exchange reaction in batch reactor. submitted to Journal of Hydrometeorology 2016. (36) Musser, J. Modeling of Heat Transfer and Reactive Chemistry for Particles in Gas– solid Flow Utilizing Continuum-discrete Methodology (CDM). Ph.D., West Virginia University, Morgantown, WV, US, 2011. (37) Cornelissen, J. T.; Taghipour, F.; Escudié, R.; Ellis, N.; Grace, J. R., CFD modelling of a liquid–solid fluidized bed. Chem. Eng. Sci. 2007, 62, 6334-6348. (38) Moldoveanu, G. A.; Papangelakis, V. G., An overview of rare-earth recovery by ionexchange leaching from ion-adsorption clays of various origins. Mineralogical Magazine 2016, 80, 63. (39) Kwauk, M., Fluidized leaching and washing (in chinese). Science Press: Beijing, 1979. (40) Kwauk, M., Fluidization — idealized and bubbleless, with applications. Ellis Horwood: New York, 1992; p 277. (41) Geldart, D., Types of gas fluidization. Powder Technol. 1973, 7, 285-292.

“For Table of Contents Only (3.25 X 1.75)


45

Coarse-Grained-Particle Method for Simulation of LiquidâSolids

Recommend Documents

Coarse-Grained-Particle Method for Simulation of LiquidâSolids