*2.2. Numerical Method*

In this work, commercial software (Ansys Fluent version 16.0) is adopted to perform this simulation. The continuous phase of drying gas is treated by an Eulerian approach, and a standard *k-ε* model is utilized for the turbulence description.

Continuity and momentum equations [23]:

$$\frac{\partial \rho}{\partial t} + \nabla \left(\rho \stackrel{\rightarrow}{\vec{v}}\right) = S\_m \tag{1}$$

$$\frac{\partial}{\partial t} \left( \rho \stackrel{\rightarrow}{\boldsymbol{\upsilon}} \right) + \nabla \cdot \left( \rho \stackrel{\rightarrow}{\boldsymbol{\upsilon}} \stackrel{\rightarrow}{\boldsymbol{\upsilon}} \right) = -\nabla p + \nabla \cdot \left[ \mu \left( \boldsymbol{\nabla} \stackrel{\rightarrow}{\boldsymbol{\upsilon}} + \boldsymbol{\nabla} \stackrel{\rightarrow}{\boldsymbol{\upsilon}} \right) - \frac{2}{3} \mu \boldsymbol{\nabla} \cdot \stackrel{\rightarrow}{\boldsymbol{\upsilon}} \stackrel{\rightarrow}{\boldsymbol{I}} \right] + \rho \stackrel{\rightarrow}{\boldsymbol{\xi}} + \stackrel{\rightarrow}{F} \tag{2}$$

where *ρ* and *v* are drying gas density and velocity, and *Sm* is the mass source term. *ρ* → *g* is the gravitational force, <sup>→</sup> *F* is the sum forces exerted by particles on the gas phase, *μ* is the drying gas effective viscosity, *I* is the unit tensor.

The standard *k-ε* is adopted to describe the flow. The turbulence kinetic energy and its rate of dissipation are obtained from the following transport equations [24]:

$$\frac{\partial}{\partial t}(\rho k) + \frac{\partial}{\partial \mathbf{x}\_{\dot{l}}}(\rho k u\_{\dot{l}}) = \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left[ \left( \mu + \frac{\mu\_{\rm f}}{\sigma\_{\rm k}} \right) \frac{\partial k}{\partial \mathbf{x}\_{\dot{j}}} \right] + \mathbf{G}\_{\dot{k}} + \mathbf{G}\_{\dot{b}} - \rho \varepsilon - \mathbf{Y}\_{\rm M} + \mathbf{S}\_{\rm k} \tag{3}$$

$$\frac{\partial}{\partial t}(\rho \varepsilon) + \frac{\partial}{\partial \mathbf{x}\_i}(\rho \varepsilon \boldsymbol{u}\_i) = \frac{\partial}{\partial \mathbf{x}\_j} \left[ \left( \mu + \frac{\mu\_t}{\sigma\_\varepsilon} \right) \frac{\partial \varepsilon}{\partial \mathbf{x}\_j} \right] + \mathbf{C}\_{1t} \frac{\varepsilon}{k} (\mathbf{G}\_k + \mathbf{C}\_{3t} \mathbf{G}\_b) - \mathbf{C}\_{2t} \rho \frac{\varepsilon^2}{k} + \mathbf{S}\_\varepsilon \tag{4}$$

where *Gk* represents the generation of turbulence kinetic energy due to the mean velocity gradients, *Gb* is the generation of turbulence kinetic energy due to buoyancy, *YM* is the generation of turbulence kinetic energy due to buoyancy.*C*1*ε*, *C*2*ε*, and *C*3*<sup>ε</sup>* are constants. *σ<sup>k</sup>* and *σε* are the turbulent Prandtl numbers for *k* and *ε*. *Sk* and S*<sup>ε</sup>* are user-defined source terms.

The trajectory of a discrete phase droplet integrates the force balance on the particle, which is written in a Lagrangian reference frame [25,26]. This force balance equates the droplet inertia with the forces acting on the droplet and can be written as:

$$\frac{d\stackrel{\rightarrow}{\boldsymbol{\mu}}\_{p}}{dt} = \boldsymbol{F}\_{\boldsymbol{D}}(\stackrel{\rightarrow}{\boldsymbol{\mu}} - \stackrel{\rightarrow}{\boldsymbol{\mu}}\_{p}) + \frac{\stackrel{\rightarrow}{\boldsymbol{\mathcal{G}}}(\boldsymbol{\rho}\_{p} - \boldsymbol{\rho})}{\rho\_{p}} + \stackrel{\rightarrow}{F} \tag{5}$$

$$F\_D = \frac{18\mu}{\rho\_p d\_p^2} \frac{C\_d \text{Re}\_p}{24} \tag{6}$$

$$\mathrm{Re}\_p = \frac{\rho d\_p \left| \stackrel{\rightarrow}{\tilde{\mathcal{U}}}\_p - \stackrel{\rightarrow}{\tilde{\mathcal{U}}} \right|}{\mu} \tag{7}$$

where <sup>→</sup> *u <sup>p</sup>* is the droplet velocity, *FD*( → *<sup>u</sup>* <sup>−</sup> <sup>→</sup> *u <sup>p</sup>*) is the drag force per unit droplet mass, <sup>→</sup> *u* is the fluid phase velocity, *<sup>ρ</sup><sup>p</sup>* is the density of the droplet, *<sup>ρ</sup>* is the fluid density, <sup>→</sup> *F* is an additional acceleration (force/unit droplet mass) term, *μ* is the molecular viscosity of the fluid, *Re* is the relative Reynolds number of gas and liquid. *Cd* is the drag coefficient.

The droplet temperature is updated according to a heat balance that relates the sensible heat change in the droplet to the convective and latent heat transfer between the droplet and the continuous phase [27]:

$$m\_p c\_p \frac{dT\_p}{dt} = hA\_p \left(T\_{\infty} - T\_p\right) \qquad T\_p < T\_{vap} \tag{8}$$

$$m\_p c\_p \frac{dT\_p}{dt} = hA\_p \left(T\_\infty - T\_p\right) + \frac{dm\_p}{dt} h\_{l\lg} \quad T\_{vap} \le T\_p < T\_{bp} \tag{9}$$

$$\frac{d(d\_p)}{dt} = \frac{4k\_\mathcal{g}}{\rho\_p c\_{p,\mathcal{g}} d\_p} (1 + 0.23\sqrt{\text{Re}\_p}) \ln\left[1 + \frac{c\_{p,\mathcal{g}}(T\_\infty - T\_p)}{h\_{l\_\mathcal{g}}}\right] \quad T\_{bp} \le T\_p \tag{10}$$

where *mp*, *cp*, *Tp*, and *Ap* are the mass, specific heat capacity at constant pressure, temperature, and surface area of droplet particles, *T*<sup>∞</sup> is the temperature in the flue gas, *Tvap* is the droplet vaporization temperature, *Tbp* is the boiling temperature of the droplet, *hlg* is the latent heat of the droplet vaporization, *cp,g* is the heat capacity of the gas, *ρ<sup>p</sup>* is the droplet density, *kg* is the thermal conductivity of the gas. *h* is the convective heat transfer coefficient, which is calculated with a modified *Nu* number as follows [28]:

$$Nu = \frac{hd\_p}{k\_{\infty}} = 2.0 + 0.6 \text{Re}\_d^{1/2} Pr^{1/3} \tag{11}$$

where *k*<sup>∞</sup> is the thermal conductivity of the continuous phase, *Red* is the Reynolds number, and *Pr* is the Prandtl number of the continuous phase.

The gas phase is described using the species transport model. The atomization model of the droplet adopts the hollow cone model. The atomization angle and inner diameter are 89◦ and 0.23 m, respectively. The droplet after atomizing follows the Rosin–Rammler distribution with the constant distribution coefficient (1.2). Different parts are meshed separately to attain a high-quality grid, as shown in Figure 3. After performing gridindependence tests, the mesh with 650,000 cells was finally used in the simulation study.

The simulated cases are listed in Table 2. Cases 1, 2, and 3 are used to study the influence of flue gas flow rate on evaporation performance. Cases 1, 4, and 5 are used to discuss the effect of flue gas temperature on evaporation performance. Cases 1, 6, and 7 are used to discuss the effect of wastewater flow rate on evaporation performance. Cases 1, 8,

and 9 are used to discuss the effect of the initial temperature of wastewater on evaporation performance. Cases 10, 11, and 12 are used to discuss the effect of the droplet size of wastewater on evaporation performance. Cases 13–30 are the orthogonal test cases of evaporation performance.

**Figure 3.** Grid division of this spray drying tower without deflectors.


