Two-Fluid-Model-Based Full Physics Simulations of Mixing in

Jun 18, 2019 - (1) or Fries et al.,(2) do not follow such a coupled approach: these studies ... (10) proposed a mathematical model for particle motion...
0 downloads 0 Views 11MB Size
Article Cite This: Ind. Eng. Chem. Res. 2019, 58, 12323−12346

pubs.acs.org/IECR

Two-Fluid-Model-Based Full Physics Simulations of Mixing in Noncohesive Wet Fluidized Beds Maryam Askarishahi,† Mohammad-Sadegh Salehi,‡ and Stefan Radl*,‡ †

Research Center Pharmaceutical Engineering GmbH, Inffeldgasse 13/III, 8010 Graz, Austria Institute of Process and Particle Engineering, Graz University of Technology, Inffeldgasse 13/III, 8010 Graz, Austria



Downloaded via NOTTINGHAM TRENT UNIV on August 13, 2019 at 07:25:33 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Mixing in a noncohesive wet fluidized bed was studied using the two-fluid model (TFM) and a zerodimensional (0D) approach. The employed TFM was extended to simulate droplet deposition on the particles, droplet evaporation, and particle drying. To quantify the bed uniformity, the variance of temperature and the particles’ loss on drying (LoD) field were computed. Subsequently, our TFM simulation data is used to support the assumptions adopted in our 0D model. Specifically, the simulation results suggest the formation of two well-mixed zones: a spraying/ wetting zone and a drying zone. Furthermore, it was demonstrated that the 0D model can accurately predict the gas and particle temperatures, as well as the moisture content with a maximum error of 4.3% when the following criteria are met: (i) low enough temperature and LoD variance (i.e., less than 5%); (ii) deep droplet penetration into the bed; (iii) no droplet loss.

1. INTRODUCTION Fluidized beds (FBs) are characterized by high rates of heat and mass transfer. FBs with liquid injection (i.e., “wet” FBs) are frequently used to produce coated particles in the food and pharmaceutical sectors. For particle coating purposes, these FB processes utilize hot gas to supply the necessary momentum for mixing and the heat for drying. Typically, vertically downwardpointing spray is used to inject a liquid to coat the particle surface. Various phenomena occur in such wet fluidized beds (WFBs) which need to be considered in addition to that in dry FBs: (i) deposition of droplets on the particle surface due to particle− droplet collisions, (ii) evaporation of freely flowing droplets, and (iii) drying (i.e., evaporation of liquid from the particle surface stemming from the flow of the heated gas over the wet particles). As a result, WFB modeling is significantly more involved compared to dry FBs. Consequently, developing simplified models for WFBs is strongly needed to gain a better quantitative understanding of the above phenomena, at least qualitatively. The most striking characteristics of WFBs, i.e., high rates of heat and mass transfer and intense agitation of the particles, as well as the position of the spray lead to the formation of distinct zones in the bed. These zones can be approximated to be wellmixed. The size of these well-mixed zones (see section 1.1.1 for a discussion of recent literature); the exchange of gas, droplet, and particle between the zones (i.e., interzone exchange, see section 1.1.2 for details); and the residence time of the particles in each zone (see section 1.1.3 for details) are the three key factors that influence the performance of WFBs. With “performance” we refer to the overall mass and heat exchange rates within the bed. © 2019 American Chemical Society

These three factors are closely coupled such that they must be analyzed together to develop process understanding. In contrast, even most of the previous flow studies on WFBs, e.g., that of Suzzi et al.1 or Fries et al.,2 do not follow such a coupled approach: these studies assumed zones with a predefined size that was not affected by particle and/or droplet flow. There is now ample evidence in the literature suggesting that such assumptions are weak, and hence should not be used in detailed flow studies. This is also of critical importance to develop the compartment models, since such models are often based on the insight gained from flow studies. In the following section we will briefly summarize the literature focusing on the aforementioned three key factors, and how they depend on the operating parameters of an FB. 1.1. Factors Influencing Wet Fluidized Bed Performance. 1.1.1. Identification of Zones and Their Size in Wet FBs. Identifying well-mixed zones in WFBs is of high importance since it helps to formulate simplified models (e.g., compartment models) for coating processes operated on the industrial scale. Extensive research has been devoted to identifying these zones. For instance, Jiménez et al.3 and Duangkhamchan et al.4 identified three well-mixed zones (i.e., wetting, isothermal, and heat transfer zones). Later, Turchiuli et al.5 identified two wellmixed zones in an FB granulator. Specifically, (i) a cool wetting zone under the nozzle and (ii) an isothermal/heat transfer zone Received: Revised: Accepted: Published: 12323

March 11, 2019 June 13, 2019 June 18, 2019 June 18, 2019 DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

kinetic theory of granular flow (KTGF). The validity of the TFM approach for the modeling of noncohesive fluidized beds is accepted to a wide extent (e.g., see the study of Fullmer and Hrenya16 that compared direct and TFM-based simulation data, or the more recent work of Ostermeier et al.17 and Banaei et al.18). KTGF was even successful in predicting the scaledependent, anisotropic fluctuating kinetic energy in a fluidized bed. It should be noted that the models and numerical parameters used in a TFM need to be carefully chosen. Next, we describe the most recent EL- and EE-based studies that focused on WFB modeling and simulation in greater detail. 1.2.1. The Euler−Lagrange Approach. EL approaches have been adopted to gain a better understanding of the WFB process. For example, applying a CFD−DEM approach, Heine et al.19 demonstrated that the spray zone of a two-fluid nozzle in a fluid bed granulator (FBG) can be described with a conical shape. Duangkhamchan et al.4 studied FBC performance via the EL approach though they disregarded the particle−droplet interaction and particle drying. Sutkar et al.20 used a CFD− DEM approach to study a WFB. They assumed that a continuous film is formed upon droplet−particle impact. In another study, Börner et al.14 employed the CFD−DEM approach to simulate an FBC. To reduce computational effort, a scaling approach21 was utilized for coarse-graining the particles. However, the application of such scaling laws for DEM simulation is still a controversial issue.22,23 Recently, Askarishahi et al.24 extended and verified25 a CFD−DEM code for the simulation of particle−droplet−fluid interaction in a WFB. They revealed the importance of considering incomplete surface coverage by droplets when estimating the drying rate. However, due to the constraint imposed by computational time, the simulated domain was very small. This necessitates the application of more affordable tools and/or speeding up the simulations. 1.2.2. The Euler−Euler Approach. Duangkhamchan et al.26 used the EE approach for the simulation of FBCs with an airblast nozzle. However, the interactions between the particles and the droplet were not accounted for. Later, Duangkhamchan et al.27 investigated the effect of atomization air pressure on the voidage distribution in an FBC using CFD, though they did not simulate the droplets. In another study, Ronsse et al.28 examined the flow behavior of gas and particle and droplets in an FBC using CFD. Nevertheless, droplet deposition and particle drying were not considered in their simulation. Szafran and Kmiec29 used the EE approach to study heat transfer in a spouted bed dryer. They reported the same temperature for the particles and the gas, which indicates fast heat transfer rates between these phases. Their axisymmetric assumption can significantly affect the reliability of their results since particle flow across the center line cannot be correctly captured with such an assumption.30 One may argue against the application of a TFM for wet fluidized bed systems since a rigorous rheological model for wet cohesive powders is not available yet. In fact, numerous studies31−33 were devoted to developing such rheological models of wet granular materials in the recent past. For instance, Roy et al.31 used DEM simulations to obtain a better insight into the rheology of wetted granular flow. They proposed that conventional rheological models need to be modified to consider other factors such as cohesion, contact softness, etc. Also, the study of Liu et al.33 has as the future goal the prediction of cohesive powder flow in the context of a kinetic theory based model. Hence, there is hope that TFMs including advanced cohesive models can be built.

(i.e., a “nonactive” zone) were identified. This necessitates the existence of a high-temperature gradient zone next to the wetting zone.6 As reported by Turchiuli et al.5 and Börner et al.,7 the size of the wetting zone was tremendously influenced by the FB’s operating parameters. In another study, Maronga and Wnukowski8,9 claimed that the wetting zone can be divided into two subzones: (i) an active-spraying zone situated next to the nozzle and (ii) an active-drying zone located between the acitvespraying zone and a nonactive zone. According to Maronga and Wnukowski,8,9 the size of the active-spraying zone is a few percent of the bed volume. 1.1.2. Exchange Rates between Zones. The flow behavior of particles plays a key role in the bed performance since it determines the rate of the particle exchange between the zones described in the previous section. We next review findings from the literature related to these exchange rates. According to Jiménez et al.,3 a higher spray rate induces faster particle circulation and consequently reduces the drying time, as supported by experimental observations of Börner et al.7 Vanderroost et al.10 proposed a mathematical model for particle motion in a fluidized bed coater (FBC). They claimed that their model can predict the solid flow behavior reliably in the absence of atomization air injection. As it will be shown in our study, atomization air significantly influences the solid flow, and clearly needs to be considered. Experimental results of Börner et al.7 revealed that particles in the freeboard are pushed downward to the spray zone due to the momentum of atomization air and the injected droplets. The particle wetting in the spray zone is assumed perfect in their study. They observed that particles in the lower circulation loop in the drying zone enter the spray zone from the bottom. They claimed that wetting of the particle entering the spray zone from the bottom is selective and imperfect. We note in passing that Börner et al.7 concluded that, in the range of umf and 4umf (umf is the minimum fluidization velocity), similar hydrodynamic behavior is observed for both rectangular and cylindrical beds. Thus, the shape of the FB is not essential, motivating studies in rectangular beds that are easier to analyze both with experiments and with simulations. 1.1.3. Particle Residence Time Distribution in Each Zone. Another important parameter influencing the coating performance of FBCs is the particle residence time in each zone and its distribution. These distributions were already the focus of previous studies.11−13 However, these studies mainly focused on agitated or drum coaters or Wurster coater. Turning to the studies devoted to fluid bed coating, Börner et al.14 employed CFD−DEM to compute the residence time of particles in the spray and drying zones. In another study,15 they demonstrated the significant effect of spray rate on the particle residence time in the spraying zone. This importance of the spray rate in FBCs will be addressed in detail in our present study. 1.2. Simulation Approaches. A variety of approaches was adopted to analyze the bed performance, of which (i) detailed models based on computational fluid dynamics (CFD) (Euler− Euler and Euler−Lagrange approaches) and (ii) compartment models are the most important ones. For CFD analysis of FBs, typically two approaches are adopted. In the Euler−Lagrange (EL) approach, the gas phase is considered as continuum and the particles are considered as discrete entities. In the Euler−Euler (EE) approach, both gas and particles are simulated as interpenetrating continua. In the EE approach, the rheological behavior of the solid phase needs special attention as it needs to be modeled, e.g., through the 12324

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

temperature and humidity field. However, these zones should be related to the regions featuring similar heat and mass exchange characteristics. As a matter of fact, the obtained temperature and humidity are only a consequence of the contribution of these exchange phenomena. Therefore, detailed flow models should be applied to identify well-mixed zones based on exchange rate characteristics. • The degree of mixing in each zone is rarely quantified, e.g., via an analysis of the temperature and moisture variance in each zone. • Particle residence times are often considered to be fixed14 despite the fact that it is well-known that they are correlated with bed configuration40 and operating parameters. 1.4. Goals and Outline. By considering the gaps in the literature addressed in the previous section, we distilled the following main goals of our present study: 1. Employ the EE approach to identify well-mixed compartments in a wet fluidized bed by examining the spatial distribution of the involved exchange phenomena (e.g., the distribution of the particle drying rate). This is done in addition to considering classical gas and particle quantities (e.g., the temperature) in a wet fluidized bed. 2. Quantify the degree of uniformity of particle LoD (stands for the loss on drying in a pharmaceutical context, which means the amount of the particle content which is evaporable upon drying) and gas temperature in the bed. This is aimed at assessing the “well-mixed” condition for individual zones in the studied fluidized bed. 3. Evaluate and improve the validity of a simple 0D model to predict wet fluidized bed performance (i.e., the temporal evolution of vapor content of gas, temperature, and LoD) via a comparison of (i) the gas and solid phase quantities predicted by a 0D model and (ii) the corresponding domain-averaged values predicted by the EE approach. This evaluation should be performed for different spray properties (injection rate and droplet velocity), bed configurations (bed size, and bed aspect ratio, i.e., fill level), and initial bed temperatures. 4. Propose criteria for the validity of the developed 0D model against data from the EE approach. These criteria are based on the degree of uniformity, the solid flow behavior, and the droplet penetration length. To achieve these goals, the open-source code MFiX (Multiphase Flow with Interphase eXchange)41 was extended to consider all relevant exchange phenomena in a wet fluidized bed (see section 2 for details). As thoroughly explained in section 1.2.1, the rheology of the granular phase was assumed to be not affected by the liquid content in most parts of our present contribution. Likewise, we consider situations in which liquid bridges are too weak to affect granular flow, or simply do not exist due to asperities on particles, or the presence of pores in the particles. In summary, MFiX was utilized to follow the EE approach to perform simulations of a so-called “two-fluid model” (TFM) using a classical rheological model for a noncohesive powder that is wetted by a liquid spray. Ultimately, particle flow and the quality of particle mixing were also studied in the case of applying a cohesive particle rheology as described in Appendices C and D in the Supporting Information. After successful verification of the newly implemented models (see Appendix A in the Supporting Information) and validation of the developed TFM approach for an experimental lab-scale fluidized bed42 (see Appendix C in the Supporting Information), a set of simulations was performed to investigate the bed performance at various operating conditions. Afterward, a 0D

Unfortunately, there is currently no consensus on a rigorous rheological model for cohesive wet granular flow. Apart from this, the liquid content of the particles is often low such that (i) the particle surface is only partially covered, (ii) the characteristic particle Stokes number is high (i.e., surface tension and liquid viscosity are comparably small, as expressed by Bond and capillary numbers, see Boyce et al.34), (iii) liquid is immediately soaked up in pores, or (iv) liquid is “hidden” between asperities as analyzed by Halsey and Levine.35 Consequently, WFBs may operate in a regime in which cohesive forces are weak, or even absent.34 As a result, one may speculate that it is reasonable to adopt the TFM approach with classical (noncohesive) rheological closures even for WFBs. Indeed, in our present work, we will demonstrate via a validation study using classical kinetic theory and a modification for cohesive particles (see Appendix C in the Supporting Information) that this is indeed the case, at least for a comparably “dry” operating point of a WFB. Specifically, to examine the fluidization behavior of cohesive particles, we implemented the rheological model recently developed by Gu et al.36 for van der Waals forces in the TFM (see Appendix D in the Supporting Information), and adopt it to model cohesive wet powder. 1.2.3. Zero-Dimensional (0D) Approaches: Compartmental Models. In what follows, we consider so-called zerodimensional models, i.e., models that do not resolve the spatial dependency of flow variables, but only consider compartments with a certain arrangement in space. A pioneer in the description of FBCs via compartment models was Sherony,37 who developed a two-compartment model based on the stochastic model for surface renewal initially presented by Hulburt and Katz.38 Later, Maronga and Wnukowski8,9 investigated a topspray FB considering a three-compartment model. However, the size of the spraying zone was a priori fixed to 10% without further justification. Adopting the same strategy, Hussain et al.39 varied the circulation time and compartment size. Nonetheless, these quantities should be determined based on detailed flow simulations. This emphasizes the necessity of using detailed flow modeling (such as CFD) for a more reliable prediction of parameters when following a compartment approach. Therefore, linking these approaches enables us to take the advantage of a high-fidelity model (via detailed simulation) and fast modeling (using a reduced-order model). 1.3. Summary of Gaps in the Literature. From the above literature it can be concluded that previous studies on wet fluidized beds suffer from several deficiencies. Specifically, the following phenomena have been often neglected in the open literature: • Droplet−particle interaction was often neglected when identifying the spraying zone.4,26−28 • When calculating the rate of drying, surface coverage (with droplets) of the particle was disregarded. This simplification can lead to huge overpredictions of the rate of drying.2,24 • The rate of droplet evaporation was marginalized in most previous studies.2,4,14,19,20,26−28 • Compartment models are often not compared against comprehensive flow models (e.g., CFD or CFD−DEM). Apart from this, the zone size and exchange rate of particles between these zones (or the residence time in these zones) are often only estimated.8,9,39 A notable exception is the work of Freireich et al.,12 which combined a DEM flow model with a sophisticated compartment model considering residence time distributions. • As reported in the open literature,3,5,8,9 the zones used in compartment models are often demarcated based on the gas 12325

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

2.2. Mass and Heat Transfer. Due to the exchange of heat, vapor, and liquid between gas, liquid, and solid phases, it is required to solve the heat and mass balance equations as presented in Table 4. The deposited liquid and evaporated vapor from the particle surface were considered as sink/source terms in the corresponding equations. In addition, the change in enthalpy due to phase change was also considered in the heat balance equation. Our modifications to the MFiX code41 are described in the following paragraphs in greater detail. 2.2.1. Closure for the Rates of Droplet Evaporation and Particle Drying. In this study, the rate of droplet evaporation on the particle surface was calculated based on the driving force to transfer water vapor between the wetted particle surface and the gas phase (for droplet evaporation, see below). This rate was computed based on the saturation density of water vapor at the respective temperature, i.e., ρsat w as

model was developed and validated against the experimental data of Heinrich and Mörl42 as well as TFM data to identify the range of reliability of this compartmental model. In section 2, the description of the governing and constitutive equations for the hydrodynamics, as well as the heat and mass transfer, will be presented.

2. MATHEMATICAL MODEL 2.1. Hydrodynamics. As described in section 1.2, the TFM approach is reliable enough to simulate WFBs characterized by a low liquid content, and/or low surface tension and viscosity of the liquid. Consequently, the Euler−Euler approach was applied to simulate the gas−solid flow in the present study. In this approach, different phases are treated mathematically as interpenetrating continua. Conservation equations were derived for each phase and are linked by correlations for interphase transfer rates of momentum, heat, and mass. For the calculation of solid phase rheological properties, the kinetic theory of granular flow (KTGF) in connection with a simple frictional stress model proposed by Schaeffer43 was used. It should be added that we assume the rheological behavior of the solid phase is not influenced by the particles’ moisture content in most parts of our study. A complete list of the governing equations is available in our previous study,44 and is summarized in Tables 1−4. The modification of standard KTGF to account for the rheology of cohesive particles is summarized in Appendix D in the Supporting Information.

̇ = (ρsat − ρ μ )adpβ εs Sdry w g wv p

This rate was added as a source and sink term for water vapor and solid liquid content, respectively, as reported in Table 4. In eq 1, βp is the mass transfer coefficient which can be calculated based on the Sherwood number defined as Sh = (βpdp)/Dvap. This coefficient has been calculated in analogy to the heat transfer coefficient correlation developed by Deen et al.45 as reported in eqs T4.11 and T4.12 in Table 4. In the case of zero gas−particle slip velocity, this correlation approaches the correct limit (i.e., Nu = Sh = 2). This assumption is realistic for the computation of the free-flowing droplet evaporation rate due to low droplet volumetric concentration in the spray region and the low droplet−gas slip velocity. In eq 1, ρsat w can be estimated based on the ideal gas law and the equilibrium vapor pressure of water; adp is the surface area of the particle available for liquid evaporation and is computed based on the correlation developed by Kariuki et al.,46 in which the surface coverage (i.e., the ratio of adp and the total particle surface area) is calculated as

Table 1. Momentum Equations Used in TFM Simulation Momentum Conservation Equations for Gas and Solid Phases

∂(εgρg Ug) ∂t ∂(εsρs Us) ∂t

+ ∇·(εgρg UgUg) = − εg∇Pg + ∇·τg̿ − Igs + εgρg g

(T1.1)

+ ∇·(εsρs UU s s) = − εs∇Pg + ∇· τs̿ + Igs + εsρs g − ∇Ps

(T1.2)

Interphase Momentum Transfer

Igs = βgs(ug − us)

ψliq = 1 − [1 − f ]Φp / f

(T1.3)

momentum exchange coefficient

βgs = 18ρg |ug − us|εg(1 − εg) F(εg , Re) = 10

1 − εg εg2

Ap

(T1.4)

coating number given by Nd f, where Nd is the number of droplets deposited on the particle surface, given by LoD (dp/ dd).3 The rate of evaporation of freely flowing droplets is computed based on the same methodology as for the evaporation from the particle surface. Specifically, we use

+ εg 2(1 + 1.5 1 − εg ) 1

+ 3εg(1 − εg) + 8.4Re−0.343 0.413Re εg + 24εg 2 1 + 103(1 − εg)Re−(1/2)(1 + 4(1 − εg))

Rep =

(T1.5)

ρg |ug − us|d p μg

̇ = (ρsat − ρ μ )εlw adβ Sevap w g wv d

(T1.6)

solid phase stress tensor ÅÄÅi ÑÉÑ 2 y Å Ñ τs = − εsÅÅÅÅjjjξs − μs zzz(∇·us)I ̿ + μs (∇us + (∇us)T )ÑÑÑÑ ÅÅÇk ÑÑÖ 3 {

(T1.7)

d

noted that droplets are not simulated as a separate phase. Instead, droplets are considered as one of the components of the gas phase; i.e., they share the same speed as the gas. Besides, it was assumed that the droplet mass loading is too low to change the density of the gas phase. This is because the spray rate is much smaller than the fluidization gas flow rate in our present study.

(T1.9)

gas phase stress tensor ÄÅ ÉÑ ÅÅi ÑÑ 2 y τg = − εgÅÅÅÅjjjξg − μg zzz(∇·ug)I ̿ + μg (∇ug + (∇ug)T )ÑÑÑÑ ÅÅÇk ÑÑÖ 3 {

(3)

where ρsat w is the saturation density of water vapor in the gas phase at the gas temperature, εlw is the volume fraction of liquid water (i.e., droplets) in the gas phase, and ad is the specific 6 surface area of a single droplet given by ad = d . It should be

solid shear viscosity μs = μs,KTGF + μs,f (T1.8) solid pressure Ps = Ps,KTGF + Ps,f

(2)

Here the parameter f is the fraction of the particle surface A coated by a single droplet (defined as d,projected ). Φp is the particle

F(εg , Re) ds 2

(1)

(T1.10) 12326

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research Table 2. Transport Equations for Granular Kinetic Energy

yz 3 ijj ∂(εsρs Θ) jj + ∇·(εsρs usΘ)zzzz = (− PI s ̿ + τs̿ ): ∇us − ∇· q − Π s + Js 2 jk ∂t {

pseudothermal conductivity of solid phase ij κ * yzÅÄÅÅi y yi 12 12 2 κs = jjjj s zzzzÅÅÅÅjjj1 + ηεsg0zzzjjj1 + η (4η − 3)εsg0zzz j g zÅÅÇk 5 5 { k { k 0{ ÉÑ Ñ 64 Ñ + (41 − 33η)η2(εsg0)2 ÑÑÑÑ ÑÑÖ 25π

κs* =

κ=

(T2.1)

(T2.2)

εsρs g0Θκ εsρs g0Θ +

6βκ 5εsρs

(T2.3)

75ρs d p π Θ 48η(41 − 33η)

(T2.4)

pseudothermal energy generation due to interaction of solid phase with gas

Πs = − 3β Θ +

81εsμg 2 |ug − us|2 g0ρs d p3 π Θ

(T2.5)

pseudothermal energy dissipation due to inelastic particle−particle collision εsg 48 Js = η(1 − η) 0 Θ3/2 π dp (T2.6)

2.2.2. Closure for the Droplet Deposition Rate. In order to simulate the deposition of droplets on the particle surface, a clean-bed filter model was utilized as suggested in the work of Kolakaluri.47 In his model, the droplet deposition rate can be calculated from ̇ = −λ|ud − u p|μ εgρ Sdep lw g

expanded bed height predicted by the TFM simulation, as well as the cross-sectional area of the FB. A schematic illustration of the developed 0D model is available in Figure 1b. The exchange rates in the bed (i.e., the last two terms in eq 6) have been computed as follows: (i) The rate of drying is calculated based on the mean voidage in the dense bed (see eq 8). (ii) The amount of evaporation (from the droplets) is estimated based on the mean travel time of the droplets in the spraying zone. In eq 6, the rates of evaporation, ṁ̇ evap, and drying, ṁ̇ dry, are computed based on the driving force available for evaporation, and are given by

(4)

where |ud − up| is the slip velocity between the fluid phase and the particle. As reported in Table 5, the filtration coefficient is a function of (i) the flow Reynolds number and (ii) the droplet Stokes number. It should be noted that the droplet Stokes number in eq T5.5 is computed based on the true gas−particle slip velocity, while Rem in eq T5.6 is calculated using the superficial velocity. After successful implementation and verification of the abovementioned models (see Appendix A in the Supporting Information for details), a set of simulations was performed to examine the effect of operation conditions on uniformity of the LoD and temperature distribution in the bed. 2.3. Compartment Model. Mass and heat balance equations of a (0D) compartment model were derived and numerically solved to track the temporal evolution of gas and particle properties (e.g., temperature, moisture, and vapor content) under the assumption of perfect mixing (see section 3.3 for a justification). Parameters of our 0D model were adopted from the results of TFM simulations. The mass balance equation for air in the gas phase can be derived as dry dry dry ṁ air = ṁ air + ṁ spray out in

(7)

ṁ dry = (ρwsat (Tp) − ρg μwv )βpadpVp

(8)

tot

Here adVdtot denotes the total surface area of droplets in the spraying zone. The total volume of droplets, Vdtot, is obtained from the spray rate and the evaporative loss of droplet volume during the droplets’ travel time in the spray region. The specific surface area of droplets can be estimated by 6/dd. In eq 8, adpVptot represents the total surface area available for evaporation in the drying zone. Particularly, the specific surface area term (i.e., adp with units of m2 per m3 particle volume) considers the particle surface coverage ψliq, and is given by 6/ dpψliq. Moreover, the total volume of particles, Vptot, is calculated based on the known initial voidage and the total bed volume. The mass transfer coefficients for the particle drying in eq 8, βp, and evaporation in eq 7, βd, can be calculated using the same method as presented in section 2.2.1. It should be added that for calculation of βp, the mean solid volume fraction in the dense bed (i.e., the drying zone) was used. The rate of droplet deposition is calculated based on the assumption of no droplet loss but considering droplet loss due to evaporation. Consequently, the mass balance equation for the particle liquid adhering to the particles can be derived as

(5)

For water vapor, the mass balance equation, using the water vapor mass loading μwv as the state variable, is given by ∂ dry dry (ρ Vgμ ) = ṁ air μ − ṁ air μ + ṁ evap + ṁ dry in wvin out wv ∂t g wv

ṁ evap = (ρwsat (Tg) − ρg μwv )βd adVd tot

(6)

Here Vg denotes the volume of the gas in the drying zone, which is calculated based on the mean voidage of the dense bed and the volume of the drying zone. The latter is estimated using the 12327

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research Table 3. Constitutive Equations for the Calculation of the Solid Stress Tensor Solid Viscosity

solid shear viscosity ÅÄÅ μs* y yi ij i 2 + α yzÅÅÅ 8 8 jj1 + g0ηεszzzjjj1 + η(3η − 2)g0εszzz zzÅÅÅ μs,KTGF = jjj 5 5 { {k k 3 {ÅÅÅÇ g0η(2 − η) k ÉÑ ÑÑ 3 Ñ + ημ b ÑÑÑÑ ÑÑ 5 ÑÖ

μs* =

(T3.1)

εsρs Θg0μ εsρs Θg0 +

5 μ= ρ dp π Θ 96 s

2βμ εsρs

(T3.2) (T3.3)

256 μb = μεsεsg0 (T3.4) 5π solid bulk viscosity

ξs =

4 Θ εsρ d pg (1 + e) 3 s 0 π

solid pressure Ps = εsρs Θ[1 + 4ηεsg0]

η=

1+e 2

(T3.5)

(T3.6)

(T3.7) Frictional Stress T

τs,f̿ = Ps,f I ̿ + + μs,f [∇us + (∇us) ]

(T3.8)

solid frictional viscosity l o yz ij P sin(φ) o o z j o o minjjjj s , μmmax zzzz, for εg < εg* o o o z j μs,f = o m ̿ { k 2 I2D o o o o o o for εg ≥ εg* o 0, n μmax = 1000 P and εg* = 0.5 in this study

(T3.9)

m

second invariant of deviatoric stress tensor ÄÅ 2 2 Å ∂us, y yz i ∂u ∂u y 1 ÅÅÅÅijj ∂us, x zz + jjj s, y − s, z zzz + − I2D ̿ = ÅÅÅjjj zz jj zz ∂y { ∂z { 6 ÅÅk ∂x k ∂y ÅÇ ÅÄÅ 2 2 ∂us, y zy i ∂u ∂u y 1 ÅÅji ∂us, x zz + jjj s, y + s, z zzz + ÅÅÅÅjjjj + z j z z j z Å ∂x { ∂y { 4 ÅÅk ∂y k ∂z ÅÇ

∂us, x yz ij ∂us, z z jj j zx − ∂x zz { k

É ÑÑ ÑÑ ÑÑ ÑÑ ÑÖ

2Ñ ÑÑ

∂us, x zy i ∂us, z zz + jjjj + ∂z z{ k ∂x

É ÑÑ ÑÑ ÑÑ ÑÑ ÑÖ

2Ñ ÑÑ

(T3.10)

frictional solid pressure

Ps,f = Fr

(εs − εs,min)n (εs,max − εs) p

Fr = 0.05, n = 2, and p = 5

mp

tot

∂ dry μ = ṁ spray μwl − (ṁ dry + ṁevap) ∂t lp

(T3.11)

∂ (m p ΔHPtot) = ṁ deposΔHliq − ṁ dry ΔH vap tot ∂t

(9)

+ hVp a p(Tg − Tp)

Here μlp represents the mass loading of the liquid in the solid phase (i.e., the particles). If LoD is considered as wet-basis mass fraction of liquid in the solid, the mass loading μlp can be related to LoD via μlp = LoD/(1 − LoD). The energy balance equation of the gas and solid phases, assuming no heat exchange between the FB and its surroundings, can be derived as

tot

We note in passing that the total energy balance of the FB can be obtained by summing the two equations above, which yields ∂ (ρ Vg ΔHg + m p ΔHPtot) tot ∂t g dry dry = ṁ air ΔHin + ṁ spray ΔHspray + (ṁ evap + ṁ depos)ΔHliq in dry − ṁ air ΔHout out

∂ dry dry (ρ Vg ΔHg) = ṁ air ΔHin + ṁ spray ΔHspray + ṁevapΔHliq in ∂t g + ṁ dry ΔH vap −

dry ṁ air ΔHout out

− hVp a p(Tg − Tp) tot

(11)

(12)

The left-hand side of eq 10 indicates the accumulation of the thermal energy in the gas phase of the FB. The first, second, and third terms on the right-hand side of the eq 10 denote the energy

(10) 12328

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research Table 4. Thermal Energy and Species Transport Equations, and Constitutive Equationsa yz ij ∂Tg εgρg CP g jjjj + ug ·∇Tg zzzz = −∇·qg + ha p(Ts − Tg) − ΔHg ∂ t { k

y i ∂T εsρs CP sjjjj s + us·∇Tszzzz = −∇·qs − ha p(Ts − Tg) − ΔHs ∂ t { k conduction heat flux in gas and solid phases qg = − εgkg∇Tg (T4.3) qs = − εsks∇Ts

(T4.1) (T4.2)

(T4.4) Heat Transfer Coefficient

h=

kgNu dp

(T4.5)

Nusselt number

Nu = (7 − 10εg + 5εg 2)(1 + 0.7Rep0.2Pr1/3) + (1.33 − 2.4εg + 1.2εg 2)Rep0.7Pr1/3

(T4.6)

Species Transport Equations ∂ ̇ ̇ − Sdepos ̇ (εgρg x wl) + ∇·(εgρg x wlug) = ∇·Dgn∇x wl + Sspray − Sevap ∂t

∂ ̇ + Sevap ̇ (εgρg x wv) + ∇·(εgρg x wvug) = ∇·Dgn∇x wv + Sdry ∂t

(T4.7)

(T4.8)

∂ ̇ ̇ (εsρs xls) + ∇·(εsρs xlsus) = ∇·Dsn∇xls + Sdepos − Sdry (T4.9) ∂t Mass Transfer Rate and Coefficient ̇ = (ρ Sdry − ρ μ )adpβ w,sat

D Sh β = wv m dp

(T4.10)

g vap

(T4.11)

Sherwood number

Shm = (7 − 10εg + 5εg 2)(1 + 0.7Rep0.2Sc1/3) + (1.33 − 2.4εg + 1.2εg 2)Rep0.7Sc1/3

(T4.12)

a

For droplet evaporation Sh was set to 2.

compared to the reference condition. For example, for the gas streams, this enthalpy can be calculated via

Table 5. Equations for the Calculation of the Droplet Deposition Rate filtration coefficient

λ = ηs

3 εs 2 ds * Steff

ΔH = Cpair(T − Tref ) + μwv [C p (T − Tref ) + ΔHlatent,ref ] vap

(T5.1)

(13)

3.2

single-collector efficiency

ηs =

effective droplet Stokes number

* = [A(εs) + 1.14Rem1/5(1 − εs)−3/2 ] Steff

A−εs function

* 4.3 + Steff

A(εs) =

3.2

It should be noted that the developed 0D model was validated against the experimental data reported by Heinrich and Mörl.42 For more details, the readers are referred to Appendix C in the Supporting Information.

(T5.2)

6 − 6εs

St 2

(T5.3)

5/3

6 − 9εs1/3 + 9εs5/3 − 6εs 2

(T5.4)

3. RESULTS AND DISCUSSION In the present work, the main focus is given to examining the validity of the 0D model. Hence, after validation of TFM and 0D model predictions against experimental data of Heinrich and Mörl,42 the TFM approach is employed to compute the temporal evolution of (i) selected sample-averaged quantities (e.g., LoD, temperature, and gas vapor content) as well as (ii) the variance (based on these samples) of LoD and temperature for a typical operating point. Afterward, the validity of the 0D model will be assessed through a comparison of the sampleaveraged quantities against the corresponding values predicted by the 0D model. Finally, the effect of operating parameters will be investigated to comprehensively evaluate the reliability of the developed 0D model.

2

droplet Stokes number particle Reynolds number

St =

|ug − us|dd ρd

Rem =

18d pρg μg εg

(T5.5)

|ug − us|(1 − εg)ρg d p μg εg

(T5.6)

input by fluidization gas, the injected spray (dry air), and the droplets (liquid) that evaporate in the spray, respectively. The fourth and fifth terms represent the energy of the outlet gas and the energy consumed by the heat exchange between gas and solid phases (i.e., Q̇ exchange in Figure 1b), respectively. In the above equations ΔH represents the specific enthalpy (in joules per kilogram dry particles or gas) of the stream 12329

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 1. Schematic illustration of the wet fluidized bed setup used for (a) TFM simulations (the black region indicates the particle bed), and (b) 0D (compartment) model simulations.

3.1. TFM Results and Qualitative Analysis. 3.1.1. Setup and Parameter Ranges. The fluidized bed setup used for simulation is schematically represented in Figure 1a. The MFiX code41 was modified to consider a conical zone as the spray injection region such that droplets are sprayed downward onto the fluidized particles. It should be noted that the droplet velocity was fixed in the spray region based on the available measurement data. The total flow rate of the injected droplets to the conical zone is equal to the injection (spray) rate from the nozzle. The fluidization gas was assumed to be uniformly distributed over the bottom inlet of the FB. The initial temperatures of gas and particle are identical in the studied cases and are summarized in Table 6. In order to examine the validity of the TFM approach, a simulation was performed for the lab-scale FB presented by Heinrich and Mörl42 for cohesive and noncohesive particles. As shown in Appendix C in the Supporting Information, the outlet gas temperature predicted by our TFM simulation is in good agreement with the observed value.42 As the bed has reached the pseudo-steady-state condition, this temperature can be considered as the wet-bulb temperature. Such an agreement demonstrated the validity of the developed TFM approach for simulations of WFBs in which particles can be assumed to be relatively dry (i.e., WFBs characterized by a low Bond number). In other words, in the studied case, the particle liquid content does not influence the flow behavior to the extent that the capillary forces defluidize the particles. It should be mentioned that a prior grid sensitivity study for mesh size of 10dp, 6dp, and 5dp proved that a grid size of 6dp is fine enough to resolve the solid and flow characteristics (i.e., voidage, temperature, and mass fractions profiles) based on the predicted time-averaged solid axial velocity, voidage, and particle temperature at various heights and grid sizes. To thoroughly evaluate the 0D model developed in our present study, a wide range of operating conditions were simulated. A summary of the studied range of various quantities can be found in Table 6. It should be noted that the spray rate has been normalized based on the drying capacity of the fluidization gas. Specifically, the normalized spray rate is norm ṁ spray =

Table 6. Simulation Conditions and Parameters Used in the TFM Approach parameter Hbed [m] Lbed [m] Hinj [m] H0/Lbed θhalf,inj [deg] ρs [kg/m3] dp [μm] ew,p epp φs Tpinit [K] LoDinit

studied range 0.75−2.25 0.23−0.69 − 0.15−0.65 − − − − − − 300−333

ρd [kg/m3] dd [μm] μg [Pa·s] ud [m/s] ṁ̇ spray [kg/s]

0.0204 Spray Properties 1000 20 1.79 × 10−5 7 2/3ṁ̇ max spray

− − − 0.1−7 0−ṁ̇ max spray

ρg [kg/m3] μg [Pa·s] Tginit [K]

Gas Phase Properties 1.188 1.79 × 10−5 333

− − 300−333

Tgi [K] u [m/s] gas solid

0.0204−0.0309



333 3.5umf Wall Boundary Condition no slip partial slip

− − −

where ṁ max spray corresponds with the spray rate at which the outlet gas would be fully saturated, (i.e., relative humidity (RH) = 100%). For more information on the procedure for calculation of ṁ max spray, the reader is referred to Appendix B in the Supporting Information. It should be mentioned that the summary of the studied operating conditions and bed configurations can be found in Table 7. To gain an insight into the transient behavior of key bed variables, exemplary contour plots for the selected gas and

ṁ spray max ṁ spray

base case Bed Geometry 0.75 0.23 0.26 0.65 15 Particle Properties 2600 400 1 1 0.5 33

(14) 12330

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Table 7. Summary of the Operating Conditions and Bed Configurations of the Cases Investigated in the Present Worka section

studied effect

3.1 3.2.1

base case 1. spray rate

3.3.1

2. initial bed temperature

3.3.2

1. bed aspect ratio 2. bed dimension 1. droplet velocity

3.3.3

2. droplet travel distance 3.3.4

atomization air

ṁ̇ norm spray

ud

0.66 0 0.5 1 0.5 1 1.5 1.5 1.5 1.5 1.5 1.5 8

7

H0/Lbed

Hbed

Lbed

Tinit

xwv

Ltrav/Ltravbase

0.65

0.75

0.23

333

0.5

1

300 300 0.65 0.15

2.25 0.75

0.69 0.96

0.1 2 3 0.1

Only the unspecified values are identical to the base case, as reported in the first row.

a

exchange rates were computed over 300 s of flow time. Averaging was started after a pseudosteady state of exchange quantities was reached. Figure 3 shows the predicted distribution for a typical case. As shown in Figure 3b,c, free droplet evaporation and deposition occur only in the spraying zone. Therefore, the evaporation of the freely flowing droplets can be assumed to be limited to a comparably small spraying zone. This finding will be used as a key assumption for computing the rates of evaporation and drying in the developed 0D model. As shown in Figure 3a, near the distributor (i.e., the bottom inlet of the FB), a higher rate of drying (with a large negative gradient in the z-direction) is predicted due to the local high driving force. In contrast, the rate of drying is almost uniform in the particle bed, which qualitatively proves the reliability of the well-mixed assumption for the studied fluidized bed. Such a uniformity of the rate of drying will be employed as a second key assumption to develop the 0D model. 3.1.3. Particle and Gas Properties. As expected, the contribution of the involved exchange phenomena influences the distribution of relevant gas and particle quantities. As depicted in Figure 4a, due to the high rate of evaporation and deposition in the spray region, droplets are totally consumed in the spraying zone. Hence, there is no droplet loss when the process parameters used to produce Figure 4 are considered. This finding will be used as a third assumption, specifically to simplify the mass balance derivation for the particle LoD in our 0D model. This high rate of evaporation results in the generation of a large amount of water vapor in the spray region as observable in Figures 3b and 4b. As is easily seen in Figure 3a,b, the rate of droplet evaporation in the spray region is more significant compared with the rate of drying in the drying zone, disregarding the zone size. Therefore, a higher mass fraction of water vapor is predicted in the spraying zone. Figure 4c reveals that the injection of droplets highly influences the distribution of the voidage in the top center of the dense bed. This is associated with the solid flow pattern affected by the droplets’ momentum acting on the particle bed: as depicted in Figure 4d, droplets push the particles downward in the central region of the bed, subsequently to the side due to high downward velocity of the droplets. According to Figure 2c, the (instantaneous) spatial distribution of the particles’ LoD is almost uniform in the drying zone (except for the region near the distributor). This is due to the efficient mixing of particles, and consequently the

particle quantities are shown in Figure 2 for t = 30 s. The instantaneous distribution of the voidage demonstrates

Figure 2. Instantaneous distribution of (a) voidage, (b) gas temperature, (c) particle LoD, and (d) water liquid mass fraction in the gas phase at t = 30 s for ṁ norm spray = 0.66 and Tinit = 330 K. The contour plot for LoD has been plotted for the cells with εs ≥ 0.1 only.

extensive bubble formation in the bed. The motion of these bubbles plays the most significant role in the overall rate of heat and mass transfer in the FB. As a result, a rather uniform distribution of temperature and particle LoD is expected in the FB, as is also discernible in Figure 2b,c. However, a lower temperature is predicted in the spraying zone due to the injection of relatively cold air and droplets, which induces a high rate of evaporation in this region. The injected droplets are totally consumed in the spraying zone due to the evaporation or deposition on the particle as is also easily discerned from the distribution of water liquid mass fraction shown in Figure 2d. 3.1.2. Exchange Rate Distribution. Aimed at assessing the “well-mixed” condition with respect to the key exchange phenomena, corresponding time-averaged values of these 12331

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 3. Distribution of key exchange rates in the wet fluidized bed: (a) drying rate, (b) evaporation rate, and (c) deposition rate for ṁ norm spray = 0.66 and Tinit = 333 K (time-averaged data).

Figure 4. Distribution of (a) water liquid mass fraction, (b) water vapor mass fraction, (c) voidage, and (d) solid flow field in a wet fluidized bed for ṁ norm spray = 0.66 and Tinit = 330 K (time-averaged data).

fairly uniform distribution of the rate of drying (see Figure 3a). It should be noted that the maximum LoD was obtained in the spraying zone because of (i) the rapid droplet deposition on the particle surface and (ii) the low driving force for drying (as noted above, the temperature is low in this region, and the air is almost saturated with vapor). 3.2. Quantification of the Degree of Mixing in the Bed. In order to quantify the uniformity of the LoD distribution in the bed, the standard deviation of this quantity was computed from TFM simulation output. To do so, the following sampling procedure was developed: 1. The computational domain is divided into sample regions (see Figure 5a). 2. Since the mass of particles in each sample must be identical, a fixed mass of solid particles is ensured in each sample as follows: the solid mass is summed up by looping over the cells in each sample region. In other words, the samples of identical solid

Figure 5. Illustration of the sampling procedure for the computation of sample-averaged LoD and its variance based on dividing the bed into sample volumes (red grid in (a); the colors indicate the voidage distribution). Illustration of the cells contributing to each sample (b) and effect of sample size on the computed LoD standard deviation (c).

volume are collected. This is necessary because the variance of the sample means depends on the sample volume. This loop continues until the desired mass of solid particles in the sample is reached. This means that, in a sample region whose voidage is higher, a larger number of cells will contribute to the sample as 12332

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 6. Influence of initial bed temperature and spray rate on temporal evolution of (a) sample-averaged LoD and (b) gas temperature.

Figure 7. Temporal evolution of LoD standard deviation in the fluidized bed for various spray rates for initial bed temperatures of (a) 333 and (b) 300 K.

shown in Figure 5b. It can also happen that the voidage in the sample is so high that the total mass of particles in all cells of that specific sample region cannot reach the desired value. Consequently, such a sample region must be skipped as indicated in Figure 5b. To evaluate the effect of number of cells in each sample, a set of simulations for various sample sizes were performed. As shown in Figure 5c, as the size of the sample increases from 4Δcell to 20Δcell in each direction, the variance of LoD decreases up to 30%. 3. The average LoD in each sample must be weighted by the solid volume fraction in the cell to get a meaningful value for the amount of liquid in each sample. Thus, we use

The same procedure was used in a slightly adopted form to also compute the temperature variance in the bed. Both procedures were successfully implemented and tested in MFiX. To test the implementation, the computational domain was divided into a certain number of cells with different LoDs and particle volume fraction, and the values computed from the TFM simulation were compared with the exact value. 3.2.1. Effect of Spray Rate and Initial Bed Temperature. After implementing the computation of LoD and temperature variance in the MFiX code,41 the effect of various operating parameters on the sample-averaged LoD and degree of uniformity were investigated. Specifically, Figure 6 presents the effect of initial bed temperature and spray rate. As easily discerned from Figure 6, the initial bed temperature has a tremendous impact on the temporal evolution of LoD. Specifically, a higher bed temperature, i.e., Tinit = 333 K, results in a fast drop of the LoD due to the larger driving force for drying. Consequently, a fast drop of the gas temperature can be observed in Figure 6b. However, the temperature stays above the temperature for the lower initial bed temperature. This defers particle liquid uptake (due to faster evaporation), followed by prolonging the drying process time, which requires more energy (for heating particles to compensate for the latent heat necessary for evaporation) and cost (for running the process until a certain LoD is obtained). The way that the spray rate influences the temporal evolution of the particles’ LoD is coupled with the initial bed temperature. Specifically, as shown in Figure 6a, at the higher initial bed temperature (e.g., Tinit = 333 K), a lower spray rate induces faster

N

LoDsample(i) =

∑ j =cell1 εscell(j , i)LoDcell (j , i) N

∑ j =cell1 εscell(j , i)

(15)

4. The sample-averaged LoD is then calculated as (since all samples have the same solid mass, no weighting is necessary here): LoDsample =

1 Nsample

Nsample



LoDsample(i)

i=1

(16)

5. The variance of the LoD can be calculated via σLoD2 =

1 Nsample − 1

Nsample



[LoDsample (i) − LoDsample ]2

i=1

(17) 12333

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

for both rheological models (i.e., the classical KTGF and the cohesive KTGF of Gu et al.36). It should be mentioned that our simulations for high Bond numbers demonstrated that, in the case of very cohesive powders, the bed will be defluidized and mixing quality is strongly reduced, as presented in Appendix D in the Supporting Information. This was also reported by Boyce et al.34 The distribution of exchange rates, as well as the gas and particle quantities in the bed, besides our quantification of the degree of mixing confirms the formation of two compartments in the studied wet fluidized bed. 1. The first compartment is a wetting/spraying zone where the following occur: (i) The droplets impact on particles. (ii) Freely flowing droplet evaporation happens. (iii) Droplets are completely consumed. (iv) The rate of particle drying is negligible compared to the rate of deposition (drying can be considered to occur only outside of this zone). 2. The second is a drying zone where the following occur: (i) Particle drying takes place uniformly in the bed. (ii) No droplet enters, and consequently the rate of droplet evaporation and deposition is negligible. (iii) The distribution of the particles’ LoD is almost uniform. (iv) It hence can be assumed as well-mixed. The communication between the two zones is established by the motion of individual particles as implied by the predicted solid flow pattern. As is easily understood from the contour plot of the exchange rates shown in Figure 3, the size of the spraying zone is relatively small compared to the drying zone (i.e., the dense bed) as is also supported by measurement data of Maronga and Wnukowski.8,9 Therefore, having considered the above-outlined characteristics of the zones, a 0D model can be employed to predict the behavior of a wet fluidized bed as follows. 1. There is a wetting/spraying compartment (called the spraying zone in what follows) in which droplets are consumed in two ways: (a) by evaporation; the rate is calculated based on the driving force for the air saturation (b) by deposition; the rate is calculated based on the mass balance for droplets (i.e., their total consumption), and the rate will be added as a source to the mass balance equation for the particles’ LoD. 2. There is also a drying/mixing compartment in which particles receive the deposited droplets from the wetting compartment and the mass and heat balance equations are derived for two phases: (a) the solid phase with the domain-averaged volume fraction and expanded bed height extracted from TFM simulation (b) the gas phase with voidage calculated based on the expanded bed height Therefore, the drying zone can be considered as the dense bed up to the expanded bed height as indicated in Figure 4c. Based on the mean voidage of the bed, the expanded bed height is given by Hexp = H0εs0 /εs . The spraying zone can be also demarcated based on the distribution of the deposition rate as depicted in Figure 3c. The length of the spray zone can be estimated as the distance between the nozzle tip and the dense bed surface, i.e., Hnozzle − Hexp, which represents the distance that droplets need

liquid loss. This is because the rate at which liquid added to the particles via deposition cannot compensate for the rate of liquid loss stemming from particle drying. However, a rise in the spray rate to ṁ norm spray = 1 leads to a situation in which the particle LoD hits a minimum, which shows the equilibrium condition in the bed. Also, reducing the spray rate increases the time when this local minimum is achieved. On the other hand, at the lower initial bed temperature (Tinit = 333 K), an increase in the spray rate to ṁ norm spray = 1 results in a faster particle liquid uptake. It is worth noting that at the lower initial bed temperature the particle LoD increases continuously as the deposition rate outweighs the drying rate. This is due to the lower driving force available for drying. Considering the temporal evolution of the temperature, according to Figure 6b, at the lower initial temperature and lower spray rate (i.e., ṁ norm spray = 0.5 and Tinit = 333 K), due to the very low rate of drying the bed temperature increases slowly over time. On the other hand, in other cases, the bed temperature decreases sharply especially for the higher initial bed temperature due to the high rate of drying. It should be added that in Figure 6b the temperature computed in the TFM simulations has been compared with the one predicted by the 0D model. This will be explained in greater detail in section 3.3. In order to evaluate the degree of uniformity with respect to LoD, the variance of LoD for these cases is computed and shown in Figure 7. As seen in Figure 7, a significant difference was predicted in the temporal evolution of LoD variance for different initial temperatures. Particularly, higher Tinit induces a higher LoD variance. This is due to particle liquid uptake in the spray region and particle liquid loss (i.e., drying) near the distributor. To be more specific, as seen in Figure 7a, at higher Tinit, δLoD increases drastically at early times, followed by a gradual decrease, and finally reaches a plateau when the rate of drying reaches its pseudo-steady-state value. In contrast to Tinit, the spray rate does not influence the temporal evolution of LoD variance significantly. Only the final value of the variance is affected by the spray rate, and to a smaller degree by the initial solid temperature. It is worth mentioning that running the simulation at various inlet gas flow rates (i.e., 1.5umf−4umf) reveals that, as long as the bed falls into the bubbling regime, the effect of the fluidization velocity on the LoD variance is not significant. This is due to the efficient particle mixing in the bubbling regime in the range of the investigated fluidization speeds. One may argue that the presence of liquid on the particle surface can influence the solid mixing behavior in the fluidized bed due to cohesion. Therefore, a set of simulations at various Bond numbers was carried out using the rheological model of Gu et al.36 This effort aimed at qualitatively evaluating the impact of the particle cohesiveness (wetness) on the fluidization behavior. As summarized in Appendix D in the Supporting Information, at a Bond number of 100 for the studied WFB, accounting for cohesive forces does not deteriorate the quality of mixing in the bed on a global level. This means that though cohesiveness decreases mixing on the particle level (as revealed by the predicted decrease of granular temperature), global mixing is not significantly impacted. Similar findings were also reported by Radl et al.,48 and more recently discussed by Luding.49 It is worth noting that, at a higher Bond number of 200, the bed is still fluidized, and a similar degree of uniformity was predicted. What is more, the validation study based on the experimental data of Heinrich and Mörl42 supports our claim: the gas temperature predicted by our TFM is in good agreement 12334

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 8. Comparison between TFM and 0D model for (a) particle LoD, (b) water vapor mass fraction in the gas phase, and (c) gas temperature at different spray rates. ṁ norm spray = 0 denotes pure particle drying.

fill level) defined as ratio of the initial bed height hinit to the bed width, (v) droplet penetration length, and (vi) atomization air flow rate. 3.3.1. Effect of Spray Rate and Initial Bed Temperature. The assumption of validity of the well-mixed condition was presented and confirmed in section 3.2 for various spray rates and initial bed temperatures. In this section, this validity will be demonstrated by a comparison of the gas and particle quantities predicted by the TFM simulation and the 0D model as reported in Figures 6 and 8. As is easily seen in Figures 6 and 8, the value of LoD, water vapor mass fraction, and gas temperature predicted by the 0D model are in very good agreement with the sample-averaged data values predicted by the TFM. As a result, it can be concluded that as long as (i) the standard deviation of LoD is small compared to the total change of LoD over the process time (LoDfinal − LoDinit) and (ii) there is no droplet loss, the bed can be assumed well-mixed and the bed performance can be predicted well by the 0D model. In other

to travel to hit the particle bed. It should be noted that in our present contribution we assume that droplet loss cannot occur. Hence, the size of the spray zone is not required for the calculation of the rate of deposition (see eq 9). As depicted in Figure 1b, rate of droplet evaporation is calculated based on the travel time of the droplet in this region and the vapor content of the gas, i.e., xwv. According to the TFM simulation results, the predicted standard deviation of xwv is less than 3% in the bed, so a uniform distribution of the water vapor can be assumed in the bed. The comparison of the bed performance using TFM and 0D model supports the validity of this assumption as described in section 3.3. 3.3. Comparison of TFM Results with 0D Model Predictions. In order to evaluate the reliability of the compartment model, the results of our TFM simulations will be next compared with the ones obtained with the 0D model for various operating parameters including (i) spray rate, (ii) initial bed temperature, (iii) bed dimension, (iv) bed aspect ratio (i.e., 12335

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 9. Comparison between TFM and 0D approaches for (a) particle LoD, (b) water vapor mass fraction fraction in gas phase, and (c) gas temperature at different bed aspect ratios.

compensated, and the average bed temperature raises, as also predicted by the TFM. This shows the deviation of the bed from well-mixed condition in the first few seconds of the process time. Since the bed temperature is lower, the driving force for drying and evaporation is lower; therefore, a lower amount of water vapor is predicted by TFM in this case. 3.3.2. Larger FB with Different Aspect Ratio. In order to examine the reliability of the developed 0D model in a larger scale, several simulations were performed for the bed three times larger in each direction, and for aspect ratios of 0.15 and 0.65 as shown in Table 7. The comparison of the results predicted by the two approaches is presented in Figure 9. As seen in Figure 9, at the early stage of the process (i.e., for t < 12 s), for both aspect ratios, the values of Tg and xwv predicted by the 0D model deviate from the TFM one (see Figure 9b,c). This can be explained by the higher residence time of the gas in the particle bed in the case of the larger bed. In other words, as the height of the bed increases, a longer time (in this case, at least H0/u = 0.75 s) is needed to reach the well-mixed condition. Besides, it takes time for the hot air to flow over all particles in the bed. Therefore, the bed deviates from the well-mixed condition, and a

words, if solid circulation is so fast that the liquid can be quickly distributed uniformly in the bed, the particle bed can be assumed well-mixed with respect to LoD. It is worth noting that a maximum droplet losses of 0.3% and 0.05% were predicted for ṁ norm spray values of 0.5 and 1, respectively. This shows that the droplet loss is negligible, and the simplification made in the droplet mass balance is valid for the 0D model. At a very high spray rate, i.e., ṁ norm spray = 3.5, a fast drop of the average bed temperature was predicted by our TFM (data not shown here). This can be associated with the high rate of evaporation stemming from the injection of a significant amount of liquid at the top of the bed, and the large droplet surface area available for evaporation. The thermal energy needed for evaporation cannot be immediately supplied from the hot gas flowing into the bed via the distributors. Therefore, the fluid flowing into the spraying zone will be cooled down. In fact, depending on the residence time of the gas in the particle bed, it can take a few seconds for the hot fluidization gas to reach the spraying zone, which causes the temperature drop in that region. As soon as the gas reaches the spraying zone, the heat loss will be 12336

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 10. Temporal evolution of the standard deviation of (a) LoD and (b) gas temperature for aspect ratio of 0.15 (ṁ norm spray = 1.5, Tinit = 333 K).

Figure 11. Comparison of time-averaged lateral and axial solid velocities at two different heights for bed aspect ratios (H0/Lbed) of (a) 0.65 and (b) 0.15.

and Wnukowski8,9 (this was done to be able to assume a uniform distribution of water vapor and temperature in the bed). In order to quantify this nonuniformity of temperature, the variance of temperature was computed in the TFM simulations and is depicted in Figure 10b. The standard deviation of around 5 K compared to total temperature change of 30 K (the latter was observed over the full flow time) proves a high degree of nonuniformity of the bed temperature with a relative standard deviation of 16.6%. (iii) Having considered the solid velocity vector field, in the case of a low bed aspect ratio, the axial convective mixing of solids (i.e., as quantified by the time-averaged axial speed) is less intense than the lateral one (see Figure 11). This hinders efficient mixing near the bed surface. Hence, the particles cannot pass the spraying zone and the injected droplets will be consumed via evaporation. In addition, the predicted solid flow pattern in Figure 12 reveals that the lateral flow of particles is more intense in the case of a lower bed aspect ratio. According to Figure 12, in the bed with the aspect ratio of 0.15 (see Figure 12c), the particles move upward in the central region of the bed, while for the higher aspect ratio (see Figure 12a), due to the deep penetration of the droplets in the dense bed, the particles are pushed downward by droplets. Hence, such a downward velocity field in the central region of the bed can be an indicator for deep droplet penetration as shown in Figure 12a. To compare the axial and lateral velocities more quantitatively, the mean velocity is calculated in the half of the bed at various heights from the distributor surface. As reported in Table

lower amount of water vapor will be generated in the bed due to the lower rate of drying as predicted by the TFM simulation. In the case of lower static bed height, i.e., considering a bed aspect ratio of 0.15, the 0D model fails to predict the wet fluidized bed performance. In fact, the highest deviations were observed for the gas temperature and the particle LoD due to the following reasons: (i) Some of the particles cannot enter the spraying zone for the liquid take-up due to deficient mixing. Consequently, the droplets will be mainly consumed for droplet evaporation, as is also predicted by the results of the TFM simulations. Such a high rate of evaporation results in the saturation of the gas in the spray region, which stems from a gas flow profile that is close to plug flow in the freeboard. Consequently, a number of droplets will leave the bed neither evaporating nor depositing. This is also supported by the rather high droplet loss of 9.6% predicted by the TFM. Besides, a lower particle liquid take-up is predicted by the TFM (i.e., a lower LoD was predicted for TFM) owing to the fact that the rate of evaporation at top of the bed in the TFM simulation is significant; (ii) Due to the lack of efficient mixing of particles, and the high rate of evaporation at top of the bed, the heat distribution cannot be assumed uniform in the bed. Specifically, two zones with comparable sizes and different temperatures will be formed in the bed. Subsequently, the variance of temperature may increase in the bed. As a matter of fact, in the developed model it was assumed that the size of the spraying zone is relatively small compared to the drying zone, as was also reported by Maronga 12337

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

from the walls and downward close to the wall. The same pattern was observed for the bed with an aspect ratio of 0.15, as shown in Figure 12c, due to the lack of droplet penetration into the dense bed. According to the voidage distribution in Figure 12a, the volume fraction of solid particles is very low in the top-central region of the dense bed which is immediately below the spraying zone. This corresponds to the inflection point observed by Börner et al.7 in the spray region. In fact, the region top of this inflection point is in the freeboard, while the region below is located in the dense bed. At this point, the spraying zone starts to shrink due to the collisions of the droplet with the particle bed. 3.3.3. Penetration Length of Droplets. Another important parameter influencing the assumption of the “well-mixed” condition in the bed is the penetration of droplets into the drying zone, governed by (i) spray position and its configuration, (ii) droplet velocity, and (iii) the distance between spray and bed surface. In this section, the effects of droplet velocity and the nozzle position on the validity of the 0D model are investigated. 3.3.3.1. Effect of Droplet Velocity. At an identical droplet injection rate of ṁ norm spray = 1.5, as stated in Table 7, the effect of droplet velocity on the bed performance was investigated. Specifically, the simulations were performed at different droplet velocities (i.e., 0.1 and 7 m/s) and compared with the 0D model. It should be noted that although the former droplet velocity is not used in classical real-world industrial applications, it was investigated to obtain the criteria for the validity of the 0D model. As depicted in Figure 13, in the case of a very low penetration length of droplet, the 0D model deviates from the TFM model due to a deficient interaction of the particles in the dense bed (i.e., drying zone) with the spraying zone. Thereby, the bed deviates from the well-mixed condition. Particularly, in the case that the droplet velocity is comparable with the fluidization velocity (in accordance with a low penetration length of droplet), an accurate computation of the droplet travel time is challenging due to the influence of upward-flowing gas on the motion of droplets. Hence, the fraction of the liquid at the top outlet must be simulated in the 0D model. In other words, significant droplet loss causes the 0D model to fail to predict the bed performance for all measured quantities. In this case, due to the deviation of the bed from the well-mixed condition, the evaporation rate and consequently water vapor mass fraction are considerably overpredicted at the very beginning of the process, and consequently a lower LoD and temperature are predicted by the 0D model. In the case of low droplet velocity, the fluidization gas will entrain the droplets. Subsequently, the gas will be saturated due to the evaporation of the entrained droplet at top of the spray region. Furthermore, this saturation will be accelerated by a gas flow profile, which is close to plug flow in the freeboard. In contrast, in the 0D model the bed is assumed well-mixed, and the gas is far from being saturated. In summary, according to Figure 13, saturation of the gas phase in the freeboard results in a higher LoD value predicted by the TFM simulation compared to the 0D model with a travel time calculated based on the droplet velocity. To support this conclusion, gas and particle quantities were computed in our 0D model for various values of droplet travel time. As shown in Figure 13, upon using different values of ttrav, it is not possible to match all quantities with identical correction factors for droplet travel time. The main reason for such a mismatch is the droplet loss in the TFM simulation which is not accounted for in our 0D

Figure 12. Predicted time-averaged voidage distribution and solid flow norm pattern for aspect ratios of (a) 0.65, ṁ norm spray = 1.5; (b) 0.65, ṁ spray = 0, ud = norm 0 (pure drying); and (c) 0.15, ṁ spray = 1.5.

8, the ratio of the lateral to axial velocity is higher for an aspect ratio of 0.15. This dominance of the lateral solid velocity leads Table 8. Comparison of the Average Upward Axial to Lateral Solid Velocity at Different Heights in the Bed for Various Bed Aspect Ratios bed aspect ratio

h/Hexp

U̅ s [m/s]

V̅ s [m/s]

U̅ s/V̅ s

0.15 0.15 0.65 0.65

0.33 0.66 0.33 0.66

0.047 0.1 0.076 0.045

0.036 0.056 0.071 0.12

1.31 1.79 1.07 0.38

the particles to having a lower chance to meet the droplets and take up liquid. Therefore, the droplets may either leave the bed (note that a nonzero water mass fraction in the outlet was predicted) or evaporate near the top of the bed. Such a droplet loss results in a higher water vapor mass fraction (see Figure 9b), as well as a lower LoD and temperature in TFM compared to the ones predicted by the 0D model as seen in Figure 9a,c. As discussed above, it can be concluded that the 0D model is appropriate if the following criteria are met: (i) no droplet loss; (ii) negligible variance of LoD and temperature compared to the total change of the corresponding values during the process time (e.g., LoDend − LoDinit); (iii) deep penetration of the droplets to the dense bed indicated by a downward pointing velocity field in the center region of the bed. Having considered the axial and lateral velocities of the solid phase in Figure 11, the injection of droplets significantly influences the motion of particles in the bed. As depicted in Figure 12a, in the case of droplet injection, the particle motion is downward in the central region of the bed in both top and bottom vortices in the dense bed. What is more, the volume fraction of particles is higher close to the wall at the top of the bed because droplets push particles toward the walls. In contrast, in the case of pure dryingas depicted in Figure 12bthe particles in the two top vortices move upward at the position far 12338

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 13. Comparison between TFM and 0D approaches for (a) particle LoD, (b) water vapor mass fraction, and (c) gas temperature for different correction factors of the droplet travel time and droplet loss.

evaporation zone is broadened at small distances from the nozzle (black solid lines in Figure 14c). This can be explained by the competition for droplets between deposition and evaporation. In fact, a higher droplet velocity gives a rise to expanded bed height. Hence, particles have a higher chance to enter the spraying zone and take up liquid (see the solid flow pattern versus water liquid mass fraction contour plot in Figure 15). Therefore, a smaller number of droplets will be available for evaporation. On the other hand, at ud = 0.1 m/s, since the particle volume fraction is very low close to the nozzle position, droplets have a very low chance to collide with particles. Thus, they either evaporate or leave the bed with the gas flow. At larger distances from the nozzle, i.e., closer to the bed surface (see solid lines in Figure 14a), a higher rate of evaporation was predicted for ud = 7 m/s compared the one at ud = 0.1 m/s. This can be justified by the fact that a larger fraction of droplets can travel downward to more distant regions from the nozzle. In other words, on the grounds that the droplets move away from the nozzle at lower ud, the evaporation rate

model. This highlights the importance of the accurate calculation of droplet travel time, as well as the droplet loss fraction, which is significant in the case of shallow droplet penetration. Therefore, the droplet loss was included via a correction factor in the 0D model to fit the particle LoD as well as the particle and gas temperatures to the corresponding values in the TFM model. As shown in Figure 13, the particle LoD and gas water vapor can be matched with one predicted by TFM with trav Closs corr = 0.83 and Ccorr = 0.025. The maximum deviation of 20% was precited for the particle temperature. Aiming at investigating the effect of droplet velocity on the involved phenomena, the rates of evaporation and deposition were obtained at different distances from the nozzle positions. As depicted in Figure 14, a higher rate of evaporation was predicted for ud = 0.1 m/s at the central region due to longer residence time of the droplets. What is more, at ud = 0.1 m/s, as plotted in Figure 14c, at higher distances from the bed surface, i.e., closer to the nozzle, the rate of evaporation outweighs the rate of deposition significantly at these distances, due to the lower probability of droplet−particle collision. Additionally, the 12339

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 14. Comparison of exchange rates for evaporation and deposition at different droplet velocities and at dimensionless vertical distances (based on spraying zone length) of (a) 0.75, (b) 0.5, and (c) 0.25 from the nozzle.

decreases and becomes smaller than the one at ud = 7 m/s, as revealed by comparing the solid lines in Figure 14. Another point discerned from Figure 14 is that at lower ud and at h/Lspray = 0.25 (see the dashed lines in Figure 14c), the rate of deposition is smaller compared to the value at higher ud. This results from the longer penetration length of droplets for higher ud, which induces particles residing close to the wall to move upward faster (for the sake of continuity). Consequently, these particles have a greater chance to enter the spraying zone and impact the droplets (see the solid lines in Figure 14c). Also, the higher rate of deposition at higher droplet speed is suggested by the filtration model that is linear in the droplet-to-particle speed. The probability for particles to enter the spraying zone highly depends on the flow behavior of particles in the bed. As depicted in Figure 15b, particles in the freeboard are pushed downward to the spraying zone due to downward velocity of droplets, and subsequently leave the bed sideways. This behavior is consistent with the observations of Börner et al.7 and Maronga and

Figure 15. Comparison of solid flow pattern (vectors colored by gas volume fraction) and water liquid mass fraction (white-to-red-to-black contours) in the gas phase for droplet velocities of (a) 0.1 and (b) 7 m/ s.

12340

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 16. Comparison of (a) LoD, (b) particle and gas temperatures, and (c) water vapor mass fraction predicted by 0D model and TFM approach for Ltrav/Ltravbase = 3 (ṁ norm spray = 1.5, Tinit = 333 K).

Figure 17. Relative deviation of 0D model prediction from TFM results (top panel) and relative standard deviation of temperature and LoD (bottom panel) for various operating parameters and bed characteristics: RLtrav = Ltrav/Ltravbase, Rṁ spr = ṁ norm spray , Rud = ud/udbase, Rfilllevel = (H0/Lbed)/(H0/Lbedbase).

particle size on the solid flow map can be found elsewhere,50,51 and supports this conclusion. 3.3.3.2. Effect of Nozzle Position. It is essential to evaluate the reliability of the developed 0D model, i.e., to investigate the situations in which the WFB has a higher degree of nonuniformity and may not fulfill the “well-mixed” condition. To do so, the nozzle position was changed in the simulation in such a way that Ltrav = 3Ltravbase; i.e., the nozzle is further away from the bed. It should be added that increasing the nozzle-tobed-distance can reduce the performance of the bed and may be far from the optimal condition for operating a WFB. Nevertheless, to evaluate the criteria proposed for the validity of the developed model, such an extreme condition, which may be out of interest for industry, was simulated.

Wnukowski.9 Despite the agreement with respect to the solid flow in the freeboard (i.e., the upper circulation loop), Börner et al.7 observed a different flow behavior for the particles in the lower circulation loop of the dense bed as the particles entered the spraying zone from the bottom. This is in contrast to the pattern predicted in the present contribution. Particularly, based on our TFM results, the particles meet the droplets at the boundary of the spraying zone due to their downward velocity. This inconsistency stems from the discrepancy in the solid flow pattern observed by Börner et al.7 and the one predicted in our present study. In fact, since the particle size used in their study falls into Geldart D classification, two vortices are observed in the dense bed. In contrast, in our current work the particles belong to Geldart B classification, and four vortices were predicted in the bed. More information regarding the effect of 12341

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

Figure 18. Comparison of (a) LoD, (b) particle and gas temperatures, and (c) water vapor mass fraction predicted by the 0D model and TFM (ṁ norm spray = 8, Tinit = 333 K).

deviation of the temperature exceeds 10%. This induces a deviation of LoD and temperature predicted by the 0D model from the one predicted by the TFM. 2. Hydrodynamics of the bed greatly affects the uniformity of the bed. To explain in more detail, downward flow of solid particles in the central region of the bed is an indication of sufficient droplet penetration. Consequently, this indicates efficient mixing of particles and droplets, i.e., complete deposition of the injected droplets. If this condition is fulfilled, no droplet loss occurs in the spraying zone, and the variance of LoD and temperature will be low enough to assume that the bed is well-mixed. Examples for an extremely insufficiently mixed bed are the cases with a low bed aspect ratio and a low droplet velocity. In the latter, the solid movement is upward in the top central region, and a remarkable deviation between 0D and TFM, for LoD and temperature, can be observed as reported in Figure 17. 3. At extremely high atomization air flow rates, the 0D model manages to predict the bed performance with the largest deviation being 4.3% for the particles’ LoD. This appears still acceptable. This deviation can be associated with the approximations made for the rate of droplet evaporation and drying in the spraying zone. In summary, to examine the validity of the developed 0D model, it is important to consider the degree of uniformity along with the droplet penetration and droplet loss. A good example is the case with low droplet velocity, as reported in Figure 17; the relative standard deviations of LoD and temperature for low droplet velocity do not exceed 6%. Nonetheless, due to insufficient penetration of droplets, as indicated by the solid flow pattern in Figure 15, the 0D model fails to predict the bed performance due to excessive droplet loss. Consequently, to use the model presented in the present work, it is essential to fulfill all proposed criteria. It should be emphasized that the extreme cases showing the high deviation are usually out of interest for industry due to depreciation of the WFB performance. What is more, the validity of our 0D model is also demonstrated based on the experimental data of Heinrich and Mö rl42 as described in Appendix C in the Supporting Information. Specifically, it was revealed that the developed model can successfully predict the outlet temperature with a maximum deviation of 8% at various spray rates.

As depicted in Figure 16 and quantified in Figure 17, in the case of Ltrav/Ltravbase = 3, the predicted values of LoD, gas temperature, and particle temperature deviate −27.9, − 15.7, and −2.3% respectively from the predicted values in the TFM. This deviation can be associated with the high degree of nonuniformity of the temperature distribution in the bed. According to the variance of the temperature presented in Figure 17b, the standard deviation of the temperature is around 10% of the total change of gas temperature over the total simulated flow time. Clearly, this nonuniformity is caused by a lack of deep penetration of droplets into the drying zone. To express this finding in more detail, in the case of a long travel distance of the droplets, the droplet kinetic energy will be dissipated, and droplets cannot penetrate the bed and push the solid particles downward. In such a case, the droplets will be entrained and will be consumed by evaporation instead of deposition on the particle surface. It should be emphasized that such a configuration that impairs the efficiency of the injected droplets is not of interest for most industrial applications, and hence a deeper discussion is skipped. 3.3.4. Effect of the Atomization Air Flow Rate. To induce a higher level of nonuniformity in the temperature distribution in the bed, both the atomization air flow rate and spray rate were extremely increased. Consequently, a high amount of gas with low enthalpy was injected into the system at ṁ norm spray = 8 and xwl = 0.1. As seen in Figures 18 and 17, the comparison between the 0D model and TFM approach reveals a relative deviation of 4.3% for particle LoD, of less than 1% for temperatures, and of 2% for xwv. This demonstrates the validity (against detailed CFD data) of the developed 0D model for the case of high atomization air flow rates. 3.4. Discussion of the Reliability of the Proposed 0D Model. Now that the dependency of bed uniformity has been quantified for various operating parameters through LoD and temperature variance, as well as a comparison between results from the 0D model and the TFM, the range of validity of the 0D model can be summarized as follows: 1. The higher the variance of LoD and temperature, the higher the degree of the deviation of the 0D model from the TFM results when predicting gas and solid quantities (i.e., LoD, vapor content, and temperature). This finding can be also supported by Figure 17. As reported in Figure 17, in the case of a low bed aspect ratio and long travel time of the droplets, the standard 12342

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

4. CONCLUSION TFM simulations of wet fluidized beds were performed to study the bed performance in terms of mixing quality and particle drying rate. The open-source code MFiX41 was extended to simulate the involved exchange phenomena including (i) droplet evaporation, (ii) droplet deposition, and (iii) particle drying in the bed. In order to examine the formation of wellmixed zones in the bed, the degree of uniformity was computed in the bed based on our TFM simulation results. It was demonstrated that the degrees of uniformity of LoD and gas temperature can be a useful indicator for the validity of the wellmixed condition in a wet fluidized bed, and consequently the validity of a 0D model. Another essential indicator to examine the “well-mixed” condition is attributed to the solid flow pattern predicted by the TFM. Specifically, if the droplets can penetrate through the dense bed and influence the flow of the particles, the bed can be assumed to be close to being well-mixed. This finding reveals the importance of considering the droplet velocity (and the spray) effect on solid flow behavior and hence bed performance for Geldart B particles. This bed−spray interaction has been ignored in most of the previous studies in which (i) the spray zone was considered as a point (and droplet trajectories as straight, massless rays1,52) or (ii) the hydrodynamics of the spray is simulated separately from the gas−solid flow in the fluidized bed to identify the spraying zone for simulation of heat and mass transfer in the bed.26,27 Based on our comparison of bed performance predicted by TFM and the 0D model, it was demonstrated that the wet fluidized bed can be assumed as “well-mixed” by a maximum deviation of around 4% for gas and particle properties for the following: (a) The degrees of nonuniformity of LoD and temperature are small in comparison to the total temporal change of the corresponding quantity during the process. To express it in a quantitative way, standard deviations of smaller than 2% for LoD and 5% for temperature (the basis for the relative temperature change is the temporal change over the process time, i.e., |Tend − Tinit|) and gas vapor content, xwv, indicate that the bed can be assumed to be well-mixed. It is worth mentioning that both criteria for the standard deviations must be fulfilled. For instance, when the bed aspect ratio was chosen to be 0.15, the standard deviation of LoD was 0.16%, while the standard deviation of the temperature exceeded 16%. Consequently, the values of LoD, temperature, and gas vapor content predicted by the 0D model deviated 10, 20, and 14% respectively (and hence appear to be unacceptable). (b) The solid flow pattern is a good indicator for sufficient droplet penetration into the bed. (c) No or only very little (i.e., less than 1%) droplet loss should occur. In the case that the “well-mixed” condition is valid and the 0D model utilizes the macroscopic bed characteristics computed from the TFM data, the 0D model can predict the bed performance accurately with a maximum deviation of 4%. In the case of excessive droplet loss from the spraying zone (e.g., as observed for a low bed aspect ratio) and low droplet penetration length (e.g., caused by a low droplet velocity), the bed deviates from the well-mixed condition and the 0D model fails to predict TFM data. Furthermore, the results of TFM simulation proved that the initial bed temperature significantly influences the temporal evolution of LoD variance, although its effect on the pseudo-

steady-state value is negligible. In contrast, the spray rate has a significant impact on the pseudo-steady-state value of LoD variance. What is more, the initial bed temperature plays a key role in the temporal evolution of the LoD. In the case of high initial bed temperature, and due to the high rate of particle drying, the LoD profile first hits a minimum, and then gradually increases. Consequently, the initial bed temperature can be optimized to reduce the process time if the goal is to increase the LoD of the particles. The 0D model presented in the current contribution can be extended to account for the entrainment of droplets. Such a development would also necessitate computing the rate of deposition based on the size of the spray zone, and the mean voidage in this deposition zone. Furthermore, an extension to a multicompartment 0D model is possible when spending additional thoughts on how to model particle exchange rates between compartments. The model can be also extended to model fluid bed granulation by including a population balance equation to describe the agglomerate size distribution. Furthermore, data from the TFM could be utilized in the calculation of agglomeration rates (e.g., via a closure for the mean particle granular temperature in the spray zone).



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.iecr.9b01344. Model verification (droplet evaporation, droplet deposition, particle drying), methodology to calculate maximum spray rate, validation of TFM and 0D model, as well as model implementation to examine the effect of powder cohesion (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel.: +43 316 873 30412. ORCID

Maryam Askarishahi: 0000-0001-7232-263X Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was partially funded by the Austrian COMET Program under the auspices of the Austrian Federal Ministry of Transport, Innovation and Technology (bmvit), the Austrian Federal Ministry of Economy, Family and Youth (bmwfj) and by the State of Styria (Styrian Funding Agency SFG). COMET is managed by the Austrian Research Promotion Agency FFG.



SYMBOLS ad = specific mass transfer surface area of droplet, 6/dd,1/m adp = specific mass transfer surface area of particle available for drying (including surface coverage of particles, (1 − εg)6/ dpψliq Cd = drag coefficient Cp = specific heat capacity, J/kg·K CPg = specific heat of fluid phase, J/kg·K CPs = specific heat of solid phase, J/kg·K dd = diameter of droplets, m dp = diameter of particles, m Dgij = rate of strain tensor, fluid phase, s−1

12343

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research Dsij = rate of strain tensor, solid phase, s−1 Dgn = diffusion coefficient of nth gas phase species, kg/m·s Dwv = molecular diffusivity of vapor in air, m2/s eeff = effective coefficient of restitution for collisions in solid phases ep,p = coefficient of restitution for collisions in solid phases ep,w = coefficient of restitution for collisions between solid particles and wall f = fraction of particle surface coated by droplet, m2/m2 g = acceleration due to gravity, m/s2 g0 = radial distribution function at contact H0 = static bed height, m Hbed = bed height, m Hexp = expanded bed height, m Hnozzle = nozzle axial distance from distributor, m ΔHg = enthalpy of phase change (evaporation) in gas phase, J/kg ΔHs = enthalpy of phase change (evaporation) in solid phase, J/kg ΔH = enthalpy, J/m3·s ΔHlatent = latent heat of vaporization, J/m3·s i, j = indices to identify vector and tensor components; summation convention is used only for these indices Igsi = interphase momentum exchange force, N/m3 Js = granular energy transfer, m2/s3 kg = fluid phase conductivity, J/m·K·s ks = solid phase conductivity, J/m·K·s Lbed = bed length, m Linj = spray zone length, m LoD = loss on drying, defined as mass fraction of liquid in particle mp = mass of particle, kg ṁ dry airin = dry-basis inlet (fluidization) air flow rate, kg/s ṁ dry airout = dry-basis outlet air flow rate, kg/s ṁ depos = rate of droplet deposition, kg/s ṁ dry = rate of particle drying, kg/s ṁ evap = rate of droplet evaporation, kg/s ṁ spray = rate of droplet injection, kg/s ṁ dry spray = dry-basis spray rate, kg/s n = constant in friction model Ncap = capillary number, μlu/σ NBo = Bond number, 6σ/ρpgdp2 Ncell = number of cells in each sample Nsample = numbers of sample in computational domain Num = Nusselt number Pg = pressure in fluid phase, Pa Ps,f = frictional pressure in solid phase, Pa Ps = solid pressure, Pa Pr = Prandtl number, Cpvfρf/kg qg = fluid phase conductive heat flux, J/m2·s qs = solid phase conductive heat flux, J/m2·s Q̇ exchange = heat exchange between gas and particle, J/s Rep = solid phase particle Reynolds number Rem = mean particle Reynolds number RH = relative humidity, % S = rate of deformation tensor, 1/s Sgij = gas phase shear rate, s−1 Ssij = solid phase shear rate, s−1 Shm = Sherwood number, βdp/Dvap Ṡ depos = rate of droplet deposition on particle, kg/s/m3 Ṡ dry = rate of particle drying, kg/s/m3 Ṡ evap = rate of droplet evaporation, kg/s/m3 Ṡ dep = rate of droplet deposition, kg/s/m3

St = droplet Stokes number St*eff = effective droplet Stokes number t = time, s ttrav = droplet travel time, s Tg = thermodynamic temperature of fluid phase, K Ts = thermodynamic temperature of solid phase, K ud = droplet velocity, m/s umf = minimum fluidization velocity, m/s Vdtot = droplet total volume, m3 Vg = gas volume in drying zone, m3 Vptot = particle total volume, m3 Ugi = fluid phase velocity vector, m/s ug = superficial gas velocity,m/s Usi = solid phase velocity vector, m/s W(x) = W-Lambert function Wbed = bed width, m xn = mass fraction of chemical species n Greek Symbols

α = constant with a value of 1.6, dimensionless γgm = fluid−solids heat transfer coefficient corrected for interphase mass transfer, J/m3·K·s β = mass transfer coefficient, m/s βgs = coefficient for interphase force between fluid phase and solid phase,kg/m3·s δij = Kronecker delta function Δcell = grid size, m εg = volume fraction of fluid phase (void fraction) ε*g = volume fraction of fluid phase in minimum fluidization condition εs,max = packed-bed (maximum) solid volume fraction εs = volume fraction of solid phase η = function of restitution coefficient, (1 + epp)/2 Θ = granular temperature of solid phase, m2/s2 κs = granular energy diffusion coefficient, kg/m·s λ = filtration coefficient, 1/m λrm = solid conductivity function μg = molecular viscosity of fluid phase, kg/m·s μlp = mass loading of particle in solid phase, kg/kg μs = solid phase shear viscosity, kg/m·s μslid = sliding friction coefficient μs,f = frictional shear viscosity of solid phase, kg/m·s μwl = mass loading of liquid in gas phase, kg/kg μwv = mass loading of vapor in gas phase, kg/kg ξ = bulk viscosity, kg/m·s Πs = exchange force in granular energy equation, kg/m·s3 ρg = microscopic (material) density of fluid phase, kg/m3 ρs = microscopic (material) density of mth solid phase, kg/m3 ρw,sat = saturation density of water vapor, kg/m3 σLoD = standard deviation of LoD σT = standard deviation of temperature τy = yield stress, Pa τ̿g = fluid phase stress tensor, Pa τ̿s = solids phase m stress tensor, Pa τ̿s,f = solid phase frictional stress tensor, Pa Φp = particle coating number φ = angle of internal friction, also used as general scalar φs = specularity coefficient ψliq = particle surface coverage Subscripts

air = air d = droplet depos = deposition 12344

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research

(7) Börner, M.; Hagemeier, T.; Ganzer, G.; Peglow, M.; Tsotsas, E. Experimental spray zone characterization in top-spray fluidized bed granulation. Chem. Eng. Sci. 2014, 116, 317. (8) Maronga, S.; Wnukowski, P. Modelling of the three-domain fluidized-bed particulate coating process. Chem. Eng. Sci. 1997, 52, 2915. (9) Maronga, S. J.; Wnukowski, P. The use of humidity and temperature profiles in optimizing the size of fluidized bed in a coating process. Chem. Eng. Process. Process Intensif. 1998, 37, 423. (10) Vanderroost, M.; Ronsse, F.; Dewettinck, K.; Pieters, J. G. Modelling overall particle motion in fluidised beds for top-spray coating processes. Particuology 2013, 11, 490. (11) Toschkoff, G.; Khinast, J. G. Mathematical modeling of the coating process. Int. J. Pharm. 2013, 457, 407. (12) Freireich, B.; Li, J.; Litster, J.; Wassgren, C. Incorporating particle flow information from discrete element simulations in population balance models of mixer-coaters. Chem. Eng. Sci. 2011, 66, 3592. (13) Crites, T.; Turton, R. Mathematical model for the prediction of cycle-time distributions for the Wurster column-coating process. Ind. Eng. Chem. Res. 2005, 44, 5397. (14) Börner, M.; Bück, A.; Tsotsas, E. DEM-CFD investigation of particle residence time distribution in top-spray fluidised bed granulation. Chem. Eng. Sci. 2017, 161, 187. (15) Börner, M.; Peglow, M.; Tsotsas, E. Particle residence times in fluidized bed granulation equipments. Chem. Eng. Technol. 2011, 34, 1116. (16) Fullmer, W. D.; Hrenya, C. M. Quantitative assessment of finegrid kinetic-theory-based predictions of mean-slip in unbounded fluidization. AIChE J. 2016, 62, 11. (17) Ostermeier, P.; DeYoung, S.; Vandersickel, A.; Gleis, S.; Spliethoff, H. Comprehensive investigation and comparison of TFM, DenseDPM and CFD-DEM for dense fluidized beds. Chem. Eng. Sci. 2019, 196, 291. (18) Banaei, M.; Jegers, J.; van Sint Annaland, M.; Kuipers, J. A.; Deen, N. G. Effect of superficial gas velocity on the solid temperature distribution in gas fluidized beds with heat production. Ind. Eng. Chem. Res. 2017, 56, 8729. (19) Heine, M.; Antonyuk, S.; Fries, L.; Niederreiter, G.; Heinrich, S.; Palzer, S. Modeling of the spray zone for particle wetting in a fluidized bed. Chem. Ing. Tech. 2013, 85, 280. (20) Sutkar, V. S.; Deen, N. G.; Patil, A. V.; Salikov, V.; Antonyuk, S.; Heinrich, S.; Kuipers, J. CFD−DEM model for coupled heat and mass transfer in a spout fluidized bed with liquid injection. Chem. Eng. J. 2016, 288, 185. (21) Link, J.; Godlieb, W.; Tripp, P.; Deen, N.; Heinrich, S.; Kuipers, J.; Schönherr, M.; Peglow, M. Comparison of fibre optical measurements and discrete element simulations for the study of granulation in a spout fluidized bed. Powder Technol. 2009, 189, 202. (22) Ozel, A.; Kolehmainen, J.; Radl, S.; Sundaresan, S. Fluid and particle coarsening of drag force for discrete-parcel approach. Chem. Eng. Sci. 2016, 155, 258. (23) Radl, S.; Sundaresan, S. A drag model for filtered Euler−Lagrange simulations of clustered gas−particle suspensions. Chem. Eng. Sci. 2014, 117, 416. (24) Askarishahi, M.; Salehi, M. S.; Radl, S. Full-physics simulations of spray-particle interaction in a bubbling fluidized bed. AIChE J. 2017, 63, 2569. (25) Salehi, M.-S.; Askarishahi, M.; Radl, S. Analytical solution for thermal transport in packed beds with volumetric heat source. Chem. Eng. J. 2017, 316, 131. (26) Duangkhamchan, W.; Ronsse, F.; Depypere, F.; Dewettinck, K.; Pieters, J. CFD study of droplet atomisation using a binary nozzle in fluidised bed coating. Chem. Eng. Sci. 2012, 68, 555. (27) Duangkhamchan, W.; Ronsse, F.; Dewettinck, K.; Pieters, J. CFD study of solids concentration in a fluidised-bed coater with variation of atomisation air pressure. Powder Technol. 2011, 212, 103. (28) Ronsse, F.; Duangkhamchan, W.; Dewettinck, K.; Pieters, J. G. Computational Fluid Dynamics (CFD) modelling of the fluidised bed coating process. In 7th International conference on Simulation and

dry = drying evap = evaporation exp = expanded bed g = gas phase in = inlet init = at time of zero liq = liquid lw = liquid water max = maximum mix = mixture out = outlet p = particle ref = reference temperature s = solid phase sample = sample tot = total vap = vapor wl = water liquid wv = water vapor Superscripts

coh = cohesion force dry = dry basis max = maximum norm = normalized sat = saturation condition Abbreviations

0D = zero-dimensional CFD = computational fluid dynamics DEM = discrete element method DNS = direct numerical simulation EE = Euler−Euler (method) EL = Euler−Lagrange (method) FB = fluidized bed FBC = fluid bed coater FBG = fluid bed granulator KTGF = kinetic theory of granular flow LoD = loss on drying MFiX = Multiphase Flow with Interphase eXchange MFM = multifluid model PBE = population balance equation TFM = two-fluid model WFB = wet fluidized bed



REFERENCES

(1) Suzzi, D.; Toschkoff, G.; Radl, S.; Machold, D.; Fraser, S. D.; Glasser, B. J.; Khinast, J. G. DEM simulation of continuous tablet coating: Effects of tablet shape and fill level on inter-tablet coating variability. Chem. Eng. Sci. 2012, 69, 107. (2) Fries, L.; Antonyuk, S.; Heinrich, S.; Palzer, S. DEM−CFD modeling of a fluidized bed spray granulator. Chem. Eng. Sci. 2011, 66, 2340. (3) Jiménez, T.; Turchiuli, C.; Dumoulin, E. Particles agglomeration in a conical fluidized bed in relation with air temperature profiles. Chem. Eng. Sci. 2006, 61, 5954. (4) Duangkhamchan, W.; Ronsse, F.; Siriamornpun, S.; Pieters, J. Numerical study of air humidity and temperature distribution in a topspray fluidised bed coating process. J. Food Eng. 2015, 146, 81. (5) Turchiuli, C.; Jimenèz, T.; Dumoulin, E. Identification of thermal zones and population balance modelling of fluidized bed spray granulation. Powder Technol. 2011, 208, 542. (6) Smith, P.; Nienow, A. Particle growth mechanisms in fluidised bed granulationI: the effect of process variables. Chem. Eng. Sci. 1983, 38, 1223. 12345

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346

Article

Industrial & Engineering Chemistry Research Modeling in the Food and Bio-industry (Foodsim’2012); Eurosis-ETI: 2012; p 34. (29) Szafran, R. G.; Kmiec, A. CFD modeling of heat and mass transfer in a spouted bed dryer. Ind. Eng. Chem. Res. 2004, 43, 1113. (30) Reuge, N.; Cadoret, L.; Coufort-Saudejaud, C.; Pannala, S.; Syamlal, M.; Caussat, B. Multifluid Eulerian modeling of dense gas− solids fluidized bed hydrodynamics: influence of the dissipation parameters. Chem. Eng. Sci. 2008, 63, 5540. (31) Roy, S.; Singh, A.; Luding, S.; Weinhart, T. Micro−macro transition and simplified contact models for wet granular materials. Comput. Part. Mech 2016, 3, 449. (32) Roy, S.; Luding, S.; Weinhart, T. Towards hydrodynamic simulations of wet particle systems. Procedia Eng. 2015, 102, 1531. (33) Liu, P.; Kellogg, K. M.; LaMarche, C. Q.; Hrenya, C. M. Dynamics of singlet-doublet collisions of cohesive particles. Chem. Eng. J. 2017, 324, 380. (34) Boyce, C.; Ozel, A.; Kolehmainen, J.; Sundaresan, S. Analysis of the effect of small amounts of liquid on gas−solid fluidization using CFD-DEM simulations. AIChE J. 2017, 63, 5290. (35) Halsey, T. C.; Levine, A. J. How sandcastles fall. Phys. Rev. Lett. 1998, 80, 3141. (36) Gu, Y.; Ozel, A.; Kolehmainen, J.; Sundaresan, S. Computationally generated constitutive models for particle phase rheology in gasfluidized suspensions. J. Fluid Mech. 2019, 860, 318. (37) Sherony, D. F. A model of surface renewal with application to fluid bed coating of particles. Chem. Eng. Sci. 1981, 36, 845. (38) Hulburt, H. M.; Katz, S. Some problems in particle technology: A statistical mechanical formulation. Chem. Eng. Sci. 1964, 19, 555. (39) Hussain, M.; Kumar, J.; Peglow, M.; Tsotsas, E. On twocompartment population balance modeling of spray fluidized bed agglomeration. Comput. Chem. Eng. 2014, 61, 185. (40) Souza, A. S.; Freire, J. T.; Béttega, R. Computational Fluid Dynamics Evaluation of the Influence of Cone Geometry on Solids Circulation in Spouted Beds. Ind. Eng. Chem. Res. 2018, 57, 13876. (41) Syamlal, M.; Rogers, W.; O’Brien, T. J. MFIX Documentation Theory Guide; DOE/METC-94/1004; U.S. Department of Energy: 1993. DOI: 10.2172/10145548. (42) Heinrich, S.; Mörl, L. Fluidized bed spray granulationa new model for the description of particle wetting and of temperature and concentration distribution. Chem. Eng. Process. 1999, 38, 635. (43) Schaeffer, D. G. Instability in the evolution equations describing incompressible granular flow. J. Differ. Equations 1987, 66, 19. (44) Salehi, M.-S.; Askarishahi, M.; Godini, H. R.; Schomäcker, R.; Wozny, G. n. CFD Simulation of Oxidative Coupling of Methane in Fluidized-Bed Reactors: A Detailed Analysis of Flow-Reaction Characteristics and Operating Conditions. Ind. Eng. Chem. Res. 2016, 55, 1149. (45) Deen, N. G.; Kriebitzsch, S. H.; van der Hoef, M. A.; Kuipers, J. Direct numerical simulation of flow and heat transfer in dense fluid− particle systems. Chem. Eng. Sci. 2012, 81, 329. (46) Kariuki, W. I.; Freireich, B.; Smith, R. M.; Rhodes, M.; Hapgood, K. P. Distribution nucleation: quantifying liquid distribution on the particle surface using the dimensionless particle coating number. Chem. Eng. Sci. 2013, 92, 134. (47) Kolakaluri, R. Direct Numerical Simulations and Analytical Modeling of Granular Filtration. Ph.D. Dissertation, Iowa State University, 2013. (48) Radl, S.; Kalvoda, E.; Glasser, B. J.; Khinast, J. G. Mixing characteristics of wet granular matter in a bladed mixer. Powder Technol. 2010, 200, 171. (49) Luding, S. Meso-scale transport in sticky granular fluids. J. Fluid Mech. 2019, 864, 1. (50) Askarishahi, M.; Salehi, M.-S.; Dehkordi, A. M. Numerical investigation on the solid flow pattern in bubbling gas−solid fluidized beds: Effects of particle size and time averaging. Powder Technol. 2014, 264, 466. (51) Askarishahi, M.; Salehi, M.-S.; Godini, H. R.; Wozny, G. CFD study on solids flow pattern and solids mixing characteristics in

bubbling fluidized bed: Effect of fluidization velocity and bed aspect ratio. Powder Technol. 2015, 274, 379. (52) Toschkoff, G.; Just, S.; Funke, A.; Djuric, D.; Knop, K.; Kleinebudde, P.; Scharrer, G.; Khinast, J. G. Spray models for discrete element simulations of particle coating processes. Chem. Eng. Sci. 2013, 101, 603.

12346

DOI: 10.1021/acs.iecr.9b01344 Ind. Eng. Chem. Res. 2019, 58, 12323−12346