A comparison of turbulence models and two and three dimensional meshes for unsteady CFD ash deposition tools

This is a version of a publication in Please cite the publication as follows: DOI: Copyright of the original publication: This is a parallel published version of an original publication. This version can differ from the original published article. published by A comparison of turbulence models and two and three dimensional meshes for unsteady CFD ash deposition tools Garcia Perez Manuel, Vakkilainen Esa García Pérez, M., Vakkilainen, E. (2019). A comparison of turbulence models and two and three dimensional meshes for unsteady CFD ash deposition tools. Fuel, Vol. 237, Issue 1 February 2019, p. 806–811. DOI: 10.1016/j.fuel.2018.10.066 Publisher's version Elsevier Fuel 10.1016/j.fuel.2018.10.066 © Elsevier 2019 Contents lists available at ScienceDirect Fuel journal homepage: www.elsevier.com/locate/fuel Full Length Article A comparison of turbulence models and two and three dimensional meshes for unsteady CFD ash deposition tools Manuel García Pérez ⁎ , Esa Vakkilainen Lappeenranta University of Technology, Energy Technology, P.O. Box 20, FIN-53851 Lappeenranta, Finland ARTICLE INFO Keywords: Fouling Combustion Ash deposition Computational fluid dynamics Turbulence model ABSTRACT This work aims to assess the adequacy of the often made two dimensional mesh simplification in ash deposition models. Little information is available regarding its validity due to the heavy computational costs that a proper three-dimensional grid model would entail. We have implemented a case study (a deposition probe in a kraft recovery furnace) with 2D and 3D mesh models in order to compare their results regarding the ash deposition and the fluid flow. An additional simulation has been carried out to compare the results between URANS and DES turbulence models. For the particular case studied in this article, the two-dimensional simplification is justified as the results did not vary notably whereas entailing remarkably smaller computational costs. Nonetheless, the usage of DES turbulence model yielded moderately different results, qualitatively closer to deposit observations, justifying perhaps the three-dimensional approach when accuracy is needed for the deposition of fine particles on the lee edges of the tubes. 1. Introduction Fly ash impaction and deposition on boiler tubes entail erosion and corrosion of the heat exchange surfaces, as well as downtimes and overall performance penalties in boilers of any kind, up to a point that these may have an impact on the overall design [1]. The physical and chemical phenomena involving ash generation, growth, transport, and deposition is thus of a major concern [1,2]. A considerable effort is being made for a better knowledge of these phenomena, CFD modeling being a particularly popular and affordable approach. Unfortunately these tools are still at an early stage [3] and would ideally be enhanced with grid-independency studies (as, e.g., [4]) and/or empirical validations (as, e.g., [5,6]), often omitted due to their challenging nature [1]. In the work of Li et al. [7] the deposition on a rhombic heat transfer tube array was modeled and validated satisfactorily. Mavridou et al. [8] analyzed the effect of the usage of tubes of different diameter within a row. Han et al. [9] compared the deposition on circular versus elliptical tube arrays. All the studies mentioned in the previous paragraph used 2D meshes for their models. Notably fewer three-dimensional simulations are typically found in literature. For instance, Wang et al. [10] studied the flow patterns and deposition over a tandem of two H-type finned tubes. The deposition and erosion on a superheater corner was modeled by Li et al. [11]. Leppänen et al. [12] predicted the fume ash formation and deposition onto the superheater surface with a 3D mesh comprising the entire furnace and the backpass of a kraft recovery boiler. All these works considered study domains with geometries which were not amenable to a two-dimensional approach, resulting into considerably heavy 3D meshes. More often than not, the works which used a 3D approach reviewed by the authors of the present study did not attempt to perform an unsteady flow simulation or a grid convergence analysis. Typically, the two-dimensional approach is executed when the modeled domain shows an appropriate periodicity and/or symmetry [7]. This way, the calculation time may be reduced by a few orders of magnitude. Frequently, for long unsteady calculations within a parametric study, this may be the only reasonable way of simulating. However, the adequacy of this very common simplification has seldom been assessed or even considered. Many deposition surfaces under study are tube arrays which experience three-dimensional turbulence eddies [13,14] which could be affecting particle trajectories. Greifzu et al. [15] observed differences between the tracked particle trajectories for a 2D and a 3D mesh of the same domain; attributing the best fitof the 3D case to the experimental results to a better simulation of the three dimensional particle dispersion caused by the turbulent eddies. Li et al. [13] compared the deposition modeled over a tube with both 2D and 3D meshes. It was concluded that the 3D approach yielded https://doi.org/10.1016/j.fuel.2018.10.066 Received 23 January 2018; Received in revised form 27 September 2018; Accepted 10 October 2018 ⁎ Corresponding author. E-mail address: manuel.garcia.perez@lut.fi(M. García Pérez). Fuel 237 (2019) 806–811 0016-2361/ © 2018 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/BY-NC-ND/4.0/). T moderately better results according to empirical measurements, although a rather coarse mesh was used, possibly resulting into an overestimation of the front deposits for the smallest particles [3]. Unfortunately, both the simulation and the experimental set-up were poorly detailed [13].Wefind fitting to consider how, in general, results on ash deposition could be affected by the 2D mesh simplification. The present study aims to confirm the validity of this 2D assumption by modeling a case study in both twoand three-dimensional meshes. The case study selected for the unsteady particle deposition was a water-cooled deposition probe inserted in the last superheater area of a kraft recovery boiler. The simulations were carried out in Ansys FLUENT 18.0 enhanced with the discrete phase model and user-defined functions. Different mesh resolutions, dimensions (in the sense of spatial coordinates), ash particle diameters and turbulence models are tested and compared. 2. Model description The ash deposition around a 4-cm (outer diameter) water-cooled probe with an outer wall temperature of 39 °C inserted in kraft recovery boiler is simulated with FLUENT. Typical kraft recovery values and gas parameters have been selected [16–18] to match the ones of the superheater region. The upstream gas comes at 686 °C with a velocity of 3.5 m/s and at a pressure of 92 kPa (outlet). The flue gas has a dynamic viscosity of − 3.82·10 5 kg/(m s), a thermal conductivity of 0.0663 W/(m K), a specific heat of 1248 J/(kg K) and a molecular weight of 28.97 kg/ kmol. Different meshes of the same case study are implemented in order to compare their results. The grid resolution is varied to ensure that the numerical convergence is reached. In addition, for one of the 3D meshes, the turbulence models URANS −k ω SST and DES are compared. Conclusions regarding the validity of the 2D simplification and turbulence models shall be drawn. In this study, the direction along the probe is often referred to as the third coordinate, which is the simplified one in the 2D simulations. The locations on the perimeter over the circular probe are determined with the angular coordinate θ, where =θ0 corresponds to the probe lee and =±θ π corresponds to the probe wind. 2.1. Computational domain and meshes The simulations executed in this study use a total of four different meshes for the same case study. The two-dimensional domains consist of a rectangle with a circular hole (4 cm diameter) which represents the probe. The flue gas comes from left to right, thus the left edge is set as a velocity inlet (located 0.2 m, or 5 times the probe diameter upstream the probe center) and the right edge is set as a pressure outlet (located 0.4 m, or 10 times the probe diameter downstream the probe center). The upper and lower edges are set as periodical boundaries so as not to constrain the flow (located each one at a distance from the probe of 0.14 m or 3.5 times the probe diameter). These 2D domains are meshed with triangular-paved schemes following two different resolution requirements to ensure the grid independency of the results. The coarse (fine) 2D meshing of the domain proceeds as follows. Firstly, the probe perimeter is divided in 600 (800) elements, accounting for the numeric accuracy guidelines proposed by Weber et al. [3,19]. A size function is then implemented to control the size of the triangular cells as they are paved further away from the probe perimeter with a growing ratio of 1.40 (1.03). The maximum allowed cell size area 35 mm 2 for both meshes. The resulting 2D grid is composed of 25256 (62246) cells. Fig. 1 shows the resulting fine mesh. The three dimensional meshes are built based on the corresponding two-dimensional ones by ’stacking’or replicating slices of them along the third coordinate. The thickness of the 3D domain is twice the probe diameter, i.e., 8 cm. The calculation time of the 3D mesh increases more than proportionally with the third dimension, and it is thus prohibitive to increase it much further. The coarse (fine) mesh is composed of 50 (70) slices of the two-dimensional mesh, resulting into a slice thickness of 1.60 (1.14) mm and a total of 1.26 (4.24) million cells. These are remarkably heavy meshes for a transient study. Fig. 2 illustrates the three dimensional coarse mesh. 2.2. Ashes and discrete particle tracking The discrete phase model available in the software package is implemented to track and calculate the motion and trajectories of kraft ash particle parcels. If a parcel impacts onto the probe surface, the userdefine routine DEFINE_DPM_EROSION [20] is called to perform the sticking-rebound submodel briefly described below. Three different ash particle size diameters are injected in the domain through the gas inlet and are calculated independently. The particle diameters under consideration are 0.7, 4.0 and 40 μ m. With this aim, multiple independent injections are set up. The dust concentration in the flue gas is 8 g/m 3 for each particle diameter, according to typical values [21]. By injecting multiple particle parcels per each inlet face (pppif) it is possible to have a sufficiently large number of parcels being tracked, and thus, to obtain a more statistically robust and less biased result after simulating the flow during a limited number of von Kármán oscillations [22]. Consequently, 10 pppif were injected at each time-step for all the simulations except for the DES case, which was instead limited to 3 pppif (thus, each individual parcel represents a larger number of particles accordingly for the same gas ash concentration) due to the particularly longer DES simulation time and the heavy calculation cost of the trajectories of over 25 million parcels in the domain (once reduced to 3 pppif). The sticking model approach used in this work is based on the mechanistic model of van Beek [23] for regular particle impactions, which has been implemented by a number of authors with satisfactory results, some of which have already been cited [4,9,10,24]. In addition, the Konstandopoulos criterion for oblique impacts [25] has been considered with the use of the empirical rebound correlations of Brach, Dunn and Li [26,27]. This combination of approaches has been used in previous simulations by the present authors. A full description of the sticking-rebound routines implemented in this study is somewhat long and fully detailed elsewhere [6]. Fig. 1. Two dimensional domain with the fine mesh resolution. Fig. 2. Three dimensional domain with the coarse mesh resolution. M. García Pérez, E. Vakkilainen Fuel 237 (2019) 806–811 807 Kraft fine ash is made mainly of inert alkali sulfates [21,28] and it has been set to the following properties: density 2664 kg/m 3 , specific heat 902 J/kg K, and thermal conductivity 0.0608 W/mK. Regarding the mechanical properties for the mechanistic stick-rebound routines, the properties of K 2 SO 4 deposits have been selected [6,23]: Young modulus 3·1010 Pa, Poisson’s ratio 0.3, Yield stress 4 .10·108Pa, friction coefficient 0.7, and surface energy 0.15 J/m 2 . The particle parcels are tracked with a customized drag law which accounts for the spatial variability of the particle Knudsen number Kn and the Cunningham correction factor C c [29,30] =μ ρd M RT Kn 2.533 p w (1) =+ ⎡ ⎣⎛ ⎝ −⎞ ⎠+⎛ ⎝ −⎞ ⎠ ⎤ ⎦ C1 Kn 1.205exp 0.0026 Kn 0.425exp 0.74 Kn c(2) where μ ρM,, wand Tare respectively the gas viscosity, density, molecular weight and temperature. dp is the particle diameter, and = R 8.31 4 J/(mol · K). The drag coefficient in this study is thus computed as − CC Dc 1 , where C D is the drag coefficient yielded by the law of Morsi and Alexandre [31]. In addition, thermophoresis simulation is enabled [32] as it has been observed to be a major deposition mechanism for the smallest types of particles [5,22,33]. 2.3. Cases under study Five different simulations are carried out in this study. Four of them are identical except for their mesh, in order to perform a proper comparison among their results. Due to the low turbulence ( = R e 122 5 )an URANS −k ω SST model is suggested [34,33]. Additionally, a fifth case is executed using the coarse 3D mesh but with a Detached-Eddy Simulation (DES) turbulence model where the near-wall region turbulence is simulated with the −k ω SST. The inlet turbulence intensity is 7%, and the viscosity ratio is 10. For the DES case, a spectral synthesizer which takes into account the inlet turbulent intensity has been enabled to generate adequately fluctuating velocity components in the inlet for a more realistic simulation. The simulation parameters are summed up in Table 1. The simulated flow time for the DES case needs to be particularly longer in order to capture a less biased particle impaction sample, since the flow was not strictly periodic as it could be appreciated in Fig. 3. The other cases may be performed with a relatively low integer number of flow oscillations due to their periodicity and the large quantity of independent particle injections implemented. 2.4. Solver and outline of the model execution The model is executed in ANSYS Fluent 18.0. The double-precision solver with default discretization schemes are used. The time-step chosen for all simulations is − 1 0 4 s, calculating 35 iterations per time step. The lift force over the probe is monitored as the flow with the particles is simulated to the point where its oscillations become quasistable. Only from that point onwards, data about each particle parcel impaction (such as geometric coordinates of the impaction point, impaction velocities and rebound velocities if applicable) on the probe is collected and registered during a time span of t(see Table 1) for ulterior postprocessing. This strategy has been followed in previous work [22,33]. 3. Results and discussion 3.1. On the computational costs Fig. 4 highlights the wall-clock time needed to compute one timestep of flow, to track the particle parcels, and to calculate the rebound of ones that require it. All the time-steps executions for this figure were carried out in the same machine, an HP Proliant SL230s G8 with two processors Intel Xeon E5-2660 (setting up Parallel Fluent with 16 threads, no hyperthreading was used), running with Linux CentOS 7. It should be noted that the large number of parcels being tracked (over 167 million parcels for the 3DF case, 89 million for the 3DC case and 26 million for the DES) require a large amount of available RAM memory, which is greater than 128 GiB for the 3DF case. The duration of a time step calculation highlights the somewhat high impracticability of using 3D models in this kind of cases due to the usually large amount of time steps to calculate. At higher flow velocities, the time step duration might have to be decreased even further, entailing an even larger number of steps to compute. 3.2. Deposition around the probe perimeter Figs. 5–7show the collected deposition rates over the probe perimeter (averaged along the third coordinates for the cases using a 3D mesh). It can be noted how for the first four cases the deposition is, at least qualitatively, equivalent. Thus, if URANS turbulence models are Table 1 Summary of executed simulations. Tis the flow oscillation period. tstands for the flow time of particle impaction logging after the flow stabilization. Sim. Mesh T[ms] t[s] Turbulence 2DC 2D, coarse 53.87 =T 5 0.2675 s URANS 2DF 2D, fine 53.62 =T 5 0.268 1 s URANS 3DC 3D, coarse 53.87 =T 5 0.2693 s URANS 3DF 3D, fine 54.03 =T 5 0.270 1 s URANS DES 3D, coarse 49.84 =T25 1.24 6 sDES Fig. 3. Temporary evolution of the aerodynamic lift force on the probe surface. Fig. 4. Average calculation time required per time-step (20 time-steps were measured for each case), in seconds. M. García Pérez, E. Vakkilainen Fuel 237 (2019) 806–811 808 used, the use of 3D meshes may not be justified. However, the DES turbulence model seems to show intuitively better results as the deposition in the leeward side of the probe is somewhat more uniform than in the other cases. Although the lee deposits are less known and understood than the wind deposits (typically less effort has been put on measuring those [1]), qualitative observations and models have reported that the rear deposits acquire a relatively uniform or flat shape [20,35]. Therefore when accuracy is needed for small particle deposition in the rear sides of tubes, the usage of 3D meshes enhanced with DES turbulence modeling might be necessary. This difference in the results of DES versus the URANS cases is ultimately caused by the flow field phenomena occurring in the probe lee. All these 0.7 and 4.0 μ m particles showed an impaction efficiency over a 99.9% consistently around the entire perimeter, meaning that the DES case is bringing fewer particles to the lee surface. Thermophoresis should be discarded as an explanation for this phenomenon, since its propensity (as proposed in previous studies [33,35]) was showed to be marginally larger for the DES case than it was for the URANS case used for comparison (3DC). Moreover, if thermophoresis happened to be significantly less intense in the DES case, then the deposition differences between the DES versus the URANS cases would be, against observations (Figs. 5 and 6), smaller for the 4.0 μ m particles than what they were with the 0.7 μ m particles. Therefore the reasons behind smaller arrival rates is explained by a smaller dust concentration in the probe wake of the DES case as a consequence of the different flow fields and vortex patterns observed when comparing the two cases at the peak of a flow oscillation, as observed in Fig. 8: in the DES case the particlecarrying flow is reattaching the wake at further location downstream form the probe and not surrounding lee-nearby vortices, nor turning away towards the probe surface. As it can be deduced also from the flow lifts in Fig. 3, the DES case presents a smaller transversal component of the flow velocity, implying again that fewer particles are dragged to the vortex-region to be then driven thermophoretically to deposition. Observe those peaks of deposition appearing at angles of ≈± απ 0.55 for the biggest particles tested in this study. Those peaks are a consequence of the upper and lower periodic boundary conditions selected: as these large particles experience a rebound after impacting in the lee of the probe, they still possess an important inertia and thus they are not much dragged by the flow as they travel nearly perpendicularly to it. After they reach the upper (lower) domain boundary, due to the periodical configuration of the model, they impact again the probe on its lower (upper) side. Such an effect may not be observed for the smaller particles since most of them do not rebound (over 99% of sticking probability has been computed for the two smaller particle diameters), and still those which do rebound are then effectively dragged and carried by the gas, not reaching those domain boundaries. 3.3. Deposition along the third coordinate It is illustrative to analyze the particle arrival observed along the third coordinate (only applicable to the cases with 3D meshes). As an example, observe Fig. 9 which highlights the number of 4.0 μ m parcel impactions onto the probe surface. It can be noted that the dependency with the third coordinate is negligible. However, it should be stressed that this fact would not constitute a sufficient condition by itself to warranty that ash deposition may be properly modeled with two-dimensional approaches. Nonetheless it is a condition that should be fulfilled. For all the other possible cases and particle diameters, an identical independency on the third coordinate was also observed. 3.4. Final comments on the limitation of the study The Finite Volume Method often implemented in CFD makes use of a discretization scheme to solve the unsteady Navier–Stokes equation system. At each iteration, these equations are applied using the discretization schemes in order to build an algebraic system of equations. The size of this system is directly related to the number of mesh cells. Thus, qualitatively, the cost of solving a mesh may grow approximately with the cube of its size. On the other hand, this size of a 3D mesh may be typically around a couple of orders of magnitude larger than that of a 2D mesh. Considering also that the calculations must be unsteady for a proper determination of the smallest particle trajectories, the oftenmade 2D assumption is frequently regarded as necessary in literature, especially with limited computing resources. However some studies might be taking this decision in cases where it is not valid. The simulations of this study could be executed in a cluster of computers with a reasonably good performance. Still, it took several months of calculation to complete due to the very high computational demand of the 3D cases. We would have preferred to conduct a more complete study with a sensitivity analysis on the number of tubes in a row, the depth of the domain along the third coordinate, and the inlet Fig. 5. Deposition rates over the probe perimeter for 0.7 μm particles. Fig. 6. Deposition rates over the probe perimeter for 4.0 μm particles. Fig. 7. Deposition rates over the probe perimeter for 40 μm particles. M. García Pérez, E. Vakkilainen Fuel 237 (2019) 806–811 809 fluid velocity. Unfortunately it is yet not reasonable to perform a study which would comprise all these analyses with the computational capability available in a typical research organization. The conclusions drawn in this study are valid for the particular implemented scenario. It would be desirable, for instance, to compare the 2D and 3D cases with higher gas velocities; but as the velocity increases, both the required mesh and the time-step resolution must become even finer [19] entailing an increase on the computational costs proportionally, approximately, to the gas velocity raised to the fifth power. Additional work should be carried on to test the validity of this assumption in a variety of typically encountered cases to have a more proper scientific confidence on the validity of the 2D or URANS simplifications. 4. Conclusions The 2D mesh simplification of the ash deposition CFD models is found frequently in literature. This study has made an attempt to assess its validity by performing CFD models over the same case. The contribution of this work is underlined by: •comparing 2D and 3D results over the same case study, •performing unsteady flow simulations, •using fine grids and confirming the grid-independency of the results, and •implementing a mechanistic particle stick–rebound routine; There is a lack of studies in literature with models analyzing the ash deposition with this combination of features. Unfortunately we did not have any empirical data at disposal for this particular probe to contrast the results. Nonetheless this work was aimed more towards highlighting the possible differences among modeling approaches than towards the development of a new one. Fig. 8. Instantaneous gas velocity fields and vectors for a peak of an oscillation with minimum (peak) lift for the 3DC (top) and DES (bottom) simulations, in the domain middle plane of the third coordinate. Note how for the 3DC case, the gas stream experiences a more sudden direction change towards the wake of the probe carrying particles into the vortices; whereas for the DES simulation this wake entrainment occurs further downstream, with a lower velocity and softer vortices. Fig. 9. Number of particle parcel impacts on the three-dimensional probe surface as a function of the probe angle and the third coordinate (along the probe). Case 3DF, 4 μ m particles. The mass of each impaction particle parcel for this case and diameter is − 1 .776·10 12 kg. M. García Pérez, E. Vakkilainen Fuel 237 (2019) 806–811 810 For the particular case study of a deposition probe in a kraft recovery furnace, it has been found that the two-dimensional simplification itself did not yield important differences in the results, but the usage of a more elaborated, three-dimensional turbulence model (DES versus URANS) did have an effect on the leeward side of the deposit. It is thus possible that the 2D simplification is justified for cases with big particles (at least 10 μ m) and/or intense turbulence. However, if precise accuracy is essential, a 3D grid with a DES turbulence model may be required for fume ash particles and mild turbulence. Therefore it is of high interest to carry out additional work corroborating the validity of the 2D simplification within a variety of geometries, inlet flow velocities and turbulence intensities. References [1] Li MJ, Tang SZ, Wang FL, Zhao QX, Tao WQ. Gas-side fouling, erosion and corrosion of heat exchangers for middle/ low temperature waste heat utilization: a review on simulation and experiment. Appl. Ther. Eng. 2017;126:737–61. [2] Baxter LL. Ash deposit formation and deposit properties. A comprehensive summary of research conducted at Sandia’s combustion research facility, tech. rep., Sandia National Labs., Albuquerque, NM (US); Sandia National Labs., Livermore, CA (US), Livermore, California; 2000. [3] Weber R, Mancini M, Schaffel-Mancini N, Kupka T. On predicting the ash behaviour using computational fluid dynamics. Fuel Process. Technol. 2013;105:113–28. [4] Balakrishnan S, Nagarajan R, Karthick K. Mechanistic modeling, numerical simulation and validation of slag-layer growth in a coal-fired boiler. Energy 2015;81:462–70. [5] Beckmann AM, Mancini M, Weber R, Seebold S, Müller M. Measurements and CFD modeling of a pulverized coal flame with emphasis on ash deposition. Fuel 2016;167:168–79. [6] García Pérez M. Modeling the effects of unsteady flow patterns on the fireside ash fouling in tube arrays of kraft and coal-fired boilers. Lappeenranta University of Technology; 2016. [Ph.D. thesis]. [7] Li J, Du W, Cheng L. Numerical simulation and experiment of gas-solid two phase flow and ash deposition on a novel heat transfer surface. Appl. Therm. Eng. 2017;113:1033–46. [8] Mavridou SG, Bouris DG. Numerical evaluation of a heat exchanger with inline tubes of different size for reduced fouling rates. Int. J. Heat Mass Transf. 2012;55(19–20):5185–95. [9] Han H, He YL, Tao WQ, Li YS. A parameter study of tube bundle heat exchangers for fouling rate reduction. Int. J. Heat Mass Transf. 2014;72:210–21. [10] Wang F-L, Tong Z-X, Songzhen T. Real-time fouling characteristics of a typical heat exchanger used in the waste heat recovery systems. Int. J. Heat Mass Transf. 2017;104(January):774–86. [11] Li B, Brink A, Hupa M. CFD investigation of slagging on a super-heater tube in a kraft recovery boiler. Fuel Process. Technol. 2013;105:149–53. [12] Leppänen A, Tran H, Taipale R, Välimäki E, Oksanen A. Numerical modeling of fine particle and deposit formation in a recovery boiler. Fuel 2014;129:45–53. [13] Li X, Zhou H, Cen K. Influences of various vortex structures on the dispersion and deposition of small ash particles. Fuel 2008;87:1379–82. [14] Zdravkovich MM. Flow around circular cylinders: volume 2: Applications, vol. 2. Oxford University Press; 2003. [15] Greifzu F, Kratzsch C, Forgber T, Lindner F, Schwarze R. Assessment of particletracking models for dispersed particle-laden flows implemented in OpenFOAM and ANSYS FLUENT. Eng. Appl. Comput. Fluid Mech. 2016;10(1):30–43. [16] Tuomenoja M, Vakkilainen EK, Saviharju K. Effect of tube arrangement on superheater heat transfer. Proceedings of the 1999 Engineering Conference. Tappi: Ahlstrom Machinery Corporation; 1999. [17] Vakkilainen EK. Offdesign operation of kraft recovery boiler [Ph.D. thesis]. Lappeenranta University of Technology; 1993. [18] Vakkilainen E. Chemical pulping part 2: Recovery of chemicals and energy, vol. 2. Helsinki Paper Engineers’Association/Paperi ja Puu Oy; 2008. p. 28–33. [19] Weber R, Schaffel-Mancini N, Mancini M, Kupka T. Fly ash deposition modelling: requirements for accurate predictions of particle impaction on tubes using RANSbased computational fluid dynamics. Fuel 2013;108:586–96. [20] Tomeczek J, Wacławiak K. Two-dimensional modelling of deposits formation on platen superheaters in pulverized coal boilers. Fuel 2009;88(8):1466–71. [21] Vakkilainen E. Kraft recovery boilers –principles and practice. Helsinki, Suomen Soodakattilayhdistys; 2005. [22] García Pérez M, Vakkilainen E, Hyppänen T. Unsteady CFD analysis of kraft recovery boiler fly-ash trajectories, sticking efficiencies and deposition rates with a mechanistic particle rebound-stick model. Fuel 2016;181:408–20. [23] van Beek MC. Gas-side fouling in heat-recovery boilers. PhD thesis [Ph.D. thesis]. Technische Universiteit Eindhoven; 2001. [24] Zhan F, Zhuang D, Ding G, Tang J. Numerical model of particle deposition on fin surface of heat exchanger. Int. J. Refrig 2016;72:27–40. [25] Konstandopoulos AG. Particle sticking/rebound criteria at oblique impact. J. Aerosol Sci. 2006;37(3):292–305. [26] Brach RM, Dunn PF, Li X. Experiments and engineering models of microparticle impact and deposition. J. Adhesion 2000;74(1–4):227–82. [27] Li X, Dunn PF, Brach RM. Lycopodium spore impacts onto surfaces. Atmos. Environ. 2000;34(10):1575–81. [28] Adams T, Frederick J, Grace TM, Hupa M, Iisa K, Jones AK, Tran H. Kraft recovery boilers. Atlanta: TAPPI Press; 1997. [29] García Pérez M, Vakkilainen E, Hyppänen T. A brief overview on the drag laws used in Lagrangian tracking of fly ash particle trajectories for boiler ash deposition CFD models. In: Proceedings of the 24nd International Conference on Impacts of Fuel Quality on Power Production, Prague, Czech Republic; 2016. [30] Allen MD, Raabe OG. Re-evaluation of millikan’s oil drop data for the motion of small particles in air. J. Aerosol Sci. 1982;13(6):537–47. [31] Morsi S, Alexander A. An investigation of particle trajectories in two-phase flow systems. J. Fluid Mech. 1972;55(2):193–208. [32] Talbot L, Cheng RK, Schefer RW, Willis DR. Thermophoresis of particles in a heated boundary layer. J. Fluid Mech. 1980;101(4):737–58. [33] García Pérez M, Vakkilainen E, Hyppänen T. The contribution of differently-sized ash particles to the fouling trends of a pilot-scale coal-fired combustor with an ash deposition CFD model. Fuel 2017;189. [34] ANSYS Inc. (US), ANSYS Fluent 17.0 Theory Guide; 2016. [35] García Pérez M, Fry A, Vakkilainen E, Whitty K. Ash deposit analysis of the convective section of a pilot-scale combustor firing two different sub-bituminous coals. Energy Fuels 2016;30(10). M. García Pérez, E. Vakkilainen Fuel 237 (2019) 806–811 811