Abstract

The ingestion of airborne particulate into aircraft engines is an undesirable consequence of their operation, particularly in and out of arid locations that leads to reduced time between overhaul. Predicting the maintenance burden in environments rich in airborne particulate is made difficult by the large number of parameters that influence the likelihood of retention of the particles on nozzle guide vanes. In this contribution, we propose a new, reduced-order model that can predict the probability of particle retention as a function of a reduced set of independent variables relating to both the carrier gas flow and particle. Two-dimensional CFD simulations of particle deposition are performed on the General Electric E3 nozzle guide vane using the existing, energy-based fouling of gas turbines (EBFOG) particle deposition model. Results from the model are compared with experimental observations of particle deposition and show good agreement with the mass fraction retained by a vane. We introduce a function that allows the probability of retention to be calculated for a range of engine operating states and architectures by defining a new dimensionless parameter, the generalized thermal Stokes number. This parameter normalizes the thermal response of a particle for all gas and particle softening temperatures allowing the retention probability function to be applied universally. Finally, we demonstrate a practical use of this model by showing its use in calculating the accumulation factor for a particle size distribution.

1 Introduction

There are several well-documented occurrences of aircraft engine flame-out due to atmospheric dust ingestion. The most well known of these is the 1982 encounter of a British Airways Boeing 747-236B with a volcanic ash cloud emanating from Mount Galunggung in Indonesia [1]. Fortunately, safety-critical events arising from high volcanic ash concentrations are rare, thanks in part to improved forecasts of volcanic ash cloud dispersion. However, the growth of the Middle East as a hub for trans-continental commercial operations has led to an increase in the number of aircraft exposed to concentrations of other mineral dusts ranging from 0.1 to 10 mg m−3 in ordinary conditions. A large commercial aircraft climbing out in such an atmosphere may ingest in the order of tens of grams of particulate dust into the core of each engine. Predicting the engine performance degradation which results from this would enable optimization of engine maintenance schedules within these operating environments.

A proportion of this dust will adhere to internal surfaces along the gas path. The journey to the high pressure turbine is a complex one, with conditioning of the dust taking place along the way, the extent of which is likely to depend on the architecture of the engine. On entry to most turboshaft engines, the dust first meets a particle separator or barrier filter, which removes almost all particles greater than 20 μm [2]. A similar, albeit less effective process is observed on turbofan engines by the centrifuging action of the fan stage [3] which diverts particles into the bypass. Particles unaffected by these separation processes which become entrained with the core gas flow may then be pulverized by high-speed rotating compressor blades, modifying the particle size and shape distributions further [4]. Finally, as the temperature increases through the high pressure compressor and combustor, particles may soften and change phase, producing a dust of markedly different mineralogy and properties to that on entry [5].

When the heated particles enter the high pressure turbine, they must negotiate the nozzle guide vanes (NGVs). If the particle comes into contact with the vane, it may be retained, resulting in a growing accumulation of particles, typically on the pressure surface (Fig. 1). This creates a flow constriction that causes a reduction in the surge margin. When the dust cloud concentration around the airframe is of the order of hundreds of milligrams per cubic meter, this process can cause an engine surge within minutes [1].

Fig. 1
Dust deposition on the leading edge of the nozzle guide vane row from a General Electric F101 turbofan engine [5]
Fig. 1
Dust deposition on the leading edge of the nozzle guide vane row from a General Electric F101 turbofan engine [5]

It was incidents such as that near Mount Galunggung which precipitated the need to understand and develop engineering models to predict the extent of engine damage due to particle ingestion. Recent advances in high performance computing have facilitated accurate computational fluid dynamics (CFD) simulations of particle deposition (see Refs. [68]) but the computation time-scales involved are not conducive to the rapid decision-making required by military operators [9]. Furthermore, owing to the large number of poorly constrained independent variables upon which deposition depends, the general applicability of these deterministic methods is limited until the ranges of these variables are better understood. In the present study, we address this by developing a reduced-order model to predict the probability that particles will be retained by an NGV. This will ultimately enable a generalized prediction of engine performance degradation in harsh environments.

2 Background

2.1 Accumulation Factor.

The extent of damage endured during a dust encounter event is a function of the dose of dust entering the engine core. Mathematically, the dose is defined as the product of the exposure duration and mean dust concentration, giving units of gram seconds per cubic meter. This approach allows calculations of the particle mass flux on the vane to be carried out. This can be expressed mathematically as Eq. (1) [10]
dmdt=C¯WcoreζNGVρf
(1)
where the rate of mass accumulation is a function of the mean dust concentration at the engine inlet C¯, the mass flowrate through the engine core Wcore, and the air density at the engine inlet ρf. These combined terms are equivalent to the rate at which dust is ingested by the engine and are dependent on the aircraft trajectory and therefore time varying. The final term in this equation is the accumulation factor ζNGV, which is defined as the proportion of ingested dust mass that deposits on a vane.

The accumulation factor is a function of the temperature-sensitive mechanical properties of the particle, the vane geometry, and the amount of deposit that is shed by aerodynamic forcing or surface erosion. A model proposed by Clarkson and Simpson [10] used two parameters to describe the overall flux of mass on the NGV: the accumulation factor ζNGV and shedding rate parameter λNGV. These parameters are applied in a simple mass balance, with the mass flux accumulated and mass flux shed compared, giving the change in deposited mass with time. Solution of the equation gives rise to an asymptotic trend of deposited mass with time, indicating that eventually, the deposited mass on the vane will reach a constant value where the accumulating and shed mass fluxes are equal and opposite. It should be noted that the accumulation factor in Eq. (1) is analogous to that used by Clarkson and Simpson [10]; therefore, the mass flux in Eq. (1) does not incorporate any deposit shedding. In an analysis of historic volcanic ash encounter events, the authors treated the accumulation factor and shedding rate parameter as constant for a particular ash composition, engine architecture and operating condition. The model was used to develop a new qualitative safety metric, in which dose limits of volcanic ash relating to operation with no safety implications, possible safety concerns and probable safety implications were defined.

Determining appropriate accumulation factors for volcanic ash has been the subject of considerable work, with composition identified as the key variable [11] driving changes in this parameter. There are a wide range of accumulation factors reported in the literature for different compositions covering volcanic ash, fly-ash, and sands, as summarized in Table 1. Basic rig testing has been used to refine the accumulation factor for volcanic ashes only. An example is the experiment of Dean et al. [12] in which ash particles were flash-heated in a plasma jet and accelerated toward a metal coupon. Four ash types of differing chemical composition exhibited accumulation factors that varied by up to 40%. A similar experiment was used by Giehl et al. [11] who showed that the rates of adhesion and removal of ash particles were broadly dependent upon the proportion of glass phases in the material. The authors concluded that compositional variables that affect the degree of melting (e.g., SiO2 content) have a significant effect on the accumulation factor.

Table 1

Summary of recent experimental studies that report accumulation factors

ReferenceParticle typeζNGV
Kim et al. [13]Ash, Sand0.000–0.048
Bons et al. [14]SynFuel Ash0.050–0.170
Crosby et al. [15]Coal Ash0.065–0.048
Dean et al. [12]Volc. Ash0.050–1.000
Taltavull et al. [16]Volc. Ash0.000–0.500
ReferenceParticle typeζNGV
Kim et al. [13]Ash, Sand0.000–0.048
Bons et al. [14]SynFuel Ash0.050–0.170
Crosby et al. [15]Coal Ash0.065–0.048
Dean et al. [12]Volc. Ash0.050–1.000
Taltavull et al. [16]Volc. Ash0.000–0.500

2.2 Processes Influencing Particle Retention.

When a particle interacts with a vane, it may form a bond and be retained on the surface. The condition for adhesion is that the attractive forces of van der Waals, surface tension or electrostatic charge overcome the elastic rebound force from the energy stored in the particle during impact. There is the additional consideration that as a particle impacts the vane, there is a mechanical response that may not simply be elastic. If the impact energy is sufficiently large, the particle may yield plastically and/or fragment upon impact. These processes consume some of the stored energy reducing the likelihood of the particle rebounding from the surface. The particle-vane impact is therefore a complex process involving exchange of momentum and energy, which depend on the mechanical response of particle and substrate, how this varies with temperature, and the energy required to form a bond between particle and substrate [8].

The deformation process a particle undergoes during impact is a strong function of its temperature, which depends on its thermal inertia and the time it spends in the surrounding gas. The particle temperature directly determines its physical state and indirectly dictates its composition by promoting chemical reactions such as decomposition. All materials consisting of a mixture of glass and crystalline phases exhibit a continuous reduction in viscosity with temperature unless this is interrupted by a phase change which can cause a discontinuity. However, the rate at which this reduces is much greater when the material exceeds its softening temperature [8]. This softening temperature is best defined as that at which the strength of the material starts to decrease sharply with further increases in temperature. Physical models predicting the reduction in viscosity as a function of temperature and chemical composition of the material [17] have been shown to provide a good approximation of this relationship for volcanic ash and coal fly-ash [7].

Some authors have attempted to fit the pattern of mineral particles reaching a softening temperature by analyzing melt behavior of mineral powders [18]. This is the current best estimate, but it is only applicable to a blend of minerals, rather than an individual particle. In reality, the response of crystalline materials during heating to their melting temperature involves a rich array of processes. For example, crystalline phases may undergo changes in crystal structure, decomposition reactions, reactions with other phases present, or none of these. These processes may have a profound effect on both the softening temperature and the mechanical response of a particle during and following impact with a substrate. It is therefore the complexity of these responses to heating which complicates the task of modeling particle retention from first principles. However, just because the underlying physical processes are complex, doesnot mean that rates of particle retention cannot be successfully captured by models with a few well-chosen variables. It is the identification of these dominant parameters which we intend to investigate in this work.

2.3 Particle Deposition Modeling.

Several numerical models have been proposed which attempt to simplify the physical processes involved in the particle-substrate collision. These models can be summarized by the following three categories:

  1. Semi-empirical, elastically based, critical velocity/moment models

  2. Models based on the probability of sticking as a function of viscosity

  3. Models built from a consideration of the mechanical interaction during impact

The first attempts at this modeling used the semi-empirical critical velocity model of Brach and Dunn [19] where the energy stored in the particle during an elastic collision is used to define a critical normal velocity for sticking. The critical normal velocity is a function of material constants relating to the particle and substrate such as their Poisson ratio and elastic modulus along with the particle diameter and work of adhesion. The energy stored in the particle during impact is compared to the work of adhesion required for the particle to be retained. As both of these depend on the normal velocity, a critical velocity can be defined where the energy stored and the work of adhesion are equal and thus the particle is retained. The model is empirically based as it requires the relationship between the contact surface area and particle normal velocity to be experimentally determined in order for the work of adhesion to be calculated. El-Batsh [20] extended this method to include deposit detachment caused by the shear flow in the boundary layer using the critical moment theory [21]. This determines the work done by the shearing force about the adhesion point which is compared with a critical moment, defined as the adhesive energy in the particle-substrate bond. If the moment acting on the particle is found to exceed the critical moment then the particle detaches and is returned to the flow. The main drawback of using the critical normal velocity as a criterion for particle retention is the lack of a dependence on temperature, something that extensive experimental testing has shown the accumulation factor to be a function of Ref. [22]. Ai and Fletcher [6] have since demonstrated that incorporating temperature dependence of the particle elastic modulus in the critical normal velocity calculation produces a reasonable agreement with experimental data for coal ash deposition.

An alternative method commonly used is the critical viscosity model. In this formulation, the temperature variation of the particle viscosity is calculated as a function of its oxide composition using the model defined by Senior and Srinivasachar [17]. This reflects the gradual softening of a particle resident in the high temperature flows typical of gas turbine engines and is based upon the observation that as a particle softens, its probability of retention increases. Sreedharan and Tafti applied this principle to model fly-ash depositing on metal coupons over a range of gas temperatures [7]. The probability of deposition was described as a function of the particle viscosity and hence its temperature. This function is defined such that it gives a sticking probability of zero below the softening temperature, increasing to unity at particle temperatures exceeding the bulk melting temperature. The probability between these bounds is defined as the ratio between the critical viscosity—the viscosity at the particle softening temperature—and the viscosity determined by the current temperature of the particle.

While the results of Sreedharan and Tafti [7] showed good agreement with experimental observations, the test was limited to relatively low gas temperatures (up to 550 K), neglected impact dynamics and was only validated for the deposition of PVC spheres at these temperatures. Singh and Tafti [23] addressed the neglecting of impact dynamics by employing a hybrid critical viscosity—elastic impact energy loss model to predict particle accumulation from a jet onto a flat plate. The probability of a particle sticking was defined as a composite function of the coefficient of restitution and viscosity of the particle and gave good agreement for the sticking probability of a range of ash particles and a particulate described only as sand.

Currently, the model most grounded in physics is that of Bons et al. [8] (referred to as the OSU model here), which incorporates both temperature-dependent elastic properties and plastic deformation. The model assumes a linear stress-strain relation for a cylinder impacting a substrate end-on and elastic compression of the particle until the yield stress is reached. At this point, the model assumes that elastic energy storage ceases and all remaining normal kinetic energy is expended in plastic deformation. This allows a critical energy storage to be defined as a function of particle diameter and the temperature-dependent elastic properties. Particles with impact energies greater than this are treated as plastically deforming. The condition for sticking is then applied as an energy balance, with the sum of the elastic or plastic stored energy and the critical moment removal energy compared to the work of adhesion. When the work of adhesion exceeds these two energy sources then then particle is assumed to be retained.

Good agreement has been shown for the case of two size distributions of a sub-bituminous ash deposited on a flat plate, and a four-vane experimental turbine nozzle cascade. Given that the premise of the model is a mechanical process based prediction of adhesion, incorporating thermo-mechanical properties is simply unavoidable. For example, Bons et al. [8] prescribed a linear reduction in yield strength with temperature in order to capture the temperature dependence of deposition. However, there are concerns over how precisely these properties are known and whether they are appropriate for the complexity of the processes involved [24].

To address some of these concerns, the energy-based fouling of gas turbines (EBFOG) model was proposed in which the probability of particle deposition can be described as a single function and a series of material dependent constants [24,25]. The constants can be found from a multi-parameter curve fitting process to experimental measurements of accumulated mass obtained at different temperatures. Experiments used to derive these accumulation factors have ranged from simple coupon deposition tests to individual NGVs and multiple NGV cascades. The implementation of this model in CFD simulations has also shown good agreement for both the bulk deposit mass and chord-wise distribution on an NGV [25]. The main difference between EBFOG and the OSU model of Bons et al. [8] is that EBFOG does not contain the physics of deposit removal by the boundary layer flow. As such this makes this a pure retention probability compared to the combined retention-detachment probability that the OSU model provides.

The prediction of particle retention via implementation of any of the models described requires the use of CFD, which is a highly deterministic, expensive process. A more useful application of these models would be in a generalized form which eliminates the need to carry out unnecessary CFD simulations.

3 Theory

3.1 Thermal Response Time.

As a particle passes through a gas stream the two phases exchange both momentum and energy with each other. The rate at which the particle responds to changes in the flow temperature is characterized by the thermal response time. This is defined as the time required for a particle to achieve 63% (or 1 − e−1) of a step change in temperature of the carrier phase [26]. It is given mathematically by Crowe et al. [26] as
τth=cp,pρpdp212kf
(2)
where cp,p is the particle specific heat capacity at constant pressure, ρp is the particle density, dp is the particle diameter, and kf is the fluid thermal conductivity. Crowe et al. [26] showed how the change in the particle temperature, Tp, as a function of time, t, can be expressed in terms of the thermal response time:
dTpdt=Nup21τth(TfTp)
(3)
where Nup is the particle Nusselt number and Tf is the temperature of the fluid surrounding the particle. Equation (3) effectively describes the thermal inertia of the particle. Figure 2 shows the particle-fluid temperature and velocity ratios, hence response times, for a 10 μm particle subjected to a step change in temperature and fluid velocity, respectively. The 63% (or 1 − 1/e) particle-fluid ratio that defines the response time is also included as the horizontal dotted line. In the case illustrated in Fig. 2, the thermal response time is an order of magnitude less than the equivalently defined velocity response time
Fig. 2
Inertial and thermal response of a 10 μm particle subjected to equivalent step changes in fluid velocity and temperature. For this case, the particle properties were assumed to be cp,p = 984, ρp=2320 and the fluid properties: Tf = 1500, kf = 8.83 E − 2, cp,f = 1211, μf = 5.26 E − 5, ρf=0.61.
Fig. 2
Inertial and thermal response of a 10 μm particle subjected to equivalent step changes in fluid velocity and temperature. For this case, the particle properties were assumed to be cp,p = 984, ρp=2320 and the fluid properties: Tf = 1500, kf = 8.83 E − 2, cp,f = 1211, μf = 5.26 E − 5, ρf=0.61.

3.2 Thermal Stokes Number.

The inertial Stokes number is useful for predicting the behavior of a particle undergoing a velocity change in the vicinity of an obstacle, in this case an NGV. The Stokes number is defined as the ratio of the respective response times of the particle and fluid to an object of characteristic length L in flow. Similarly, a thermal Stokes number can be helpful in estimating whether a particle will reach the temperature of its carrier gas field in the time taken to translate across the vane. Broadly speaking, if the inertial Stokes number of the particle is greater than its thermal Stokes number then there is a high likelihood that the particle will be in equilibrium with the flow temperature as it traverses the NGV.

Singh and Tafti [27] proposed that the thermal Stokes number can be expressed in terms of the freestream velocity, Uin, and the particle convective heat transfer coefficient, h, defined as the rate of heat transfer between a particle surface and carrier fluid per unit surface area. Similar to this was the thermal Stokes number for NGVs suggested by Bojdo and Filippone [28], in which the heat transfer coefficient was replaced by the non-dimensional particle Nusselt number Nup, which is a function of the particle diameter, dp, and fluid thermal conductivity, kf:
Stkth=ρpcp,pdpUin6hL=ρpcp,pdp26kf1NupUinL
(4)
It is difficult to use these definitions of thermal Stokes number in a way that accommodates changes in state of the particle prior to impact with the vane, such as allowing the particle properties to vary with time or temperature, as they would during passage through the NGV stage of an engine. Incorporation of this behavior would allow the probability of retention to be defined as a function of the thermal Stokes number reflecting the properties of the particle upon impact.

3.3 Reynolds Number Effects.

Particles that exit the highly turbulent, high velocity flow from the combustor outlet experience large velocity gradients as they are accelerated toward the nozzle guide vanes. This can be expressed in terms of the non-dimensional particle Reynolds number:
Rep=ρf|uv|dpμf
(5)
where ρf, u, and μf are the density, velocity, and viscosity of the fluid, respectively. v and dp are the velocity and density of the particle, respectively.
The inertial drag force acting on the particle also increases with the particle Reynolds number, resulting in an increase of the particle drag coefficient as the particle exits the linear, Stokesian drag regime. As the classical Stokes number is derived for a particle whose Reynolds number lies in this linear range, a corrected version, the generalized Stokes number for an NGV has been proposed by Bojdo et al. [29]. The generalized Stokes number can be expressed in terms of engine design and operating state parameters including the particle Reynolds number, Rep, engine core mass flowrate, W1, fluid density, ρf, and the product of the number of vanes, N, and the vane height, h:
Stkgen=τvτf=ρpdp18(Repψ)(Nh)UW1
(6)
The generalized Stokes number uses the Israel and Rosner [30] non-Stokes drag correction factor, ψ, a function of the particle relative Reynolds number that corrects for non-linearity of the drag force for particles outside the Stokesian drag regime [29]. Given that the particle Reynolds number is also related to the particle Nusselt number and therefore the heat transfer rate between particle and fluid, it follows that particles that enter the inertial drag regime will be subject to an enhanced heat transfer rate compared to that in the purely Stokesian regime. This has implications for the impact mechanics as these particles will see reductions in their mechanical properties commensurate with this heating rate. It is this effect which we incorporate into a generalized model in this work.

4 Methodology

The particle fate models previously discussed have routinely been applied to individual particles as user defined sub-routines within high fidelity CFD models to predict whether or not the particle sticks. This condition has to be applied for each particle, with tens of thousands often required per simulation in order to obtain statistically significant results. There is the added drawback that the solutions of these high fidelity models are highly deterministic. For example, to investigate multiple combinations of fluid and particle properties, the CFD model would need to be run repeatedly for all cases. This leads to a significant computational burden and requires high performance computing in order to provide timely solutions. Even with these resources, a solution for a single combination of particulate and engine operating state can take anywhere from hours to days to solve. In order to reduce this computational burden, we therefore seek a generalized model for particle retention that uses inertial and thermal Stokes numbers to predict the likelihood of particle retention using two curve fit equations describing the probabilities of interaction and retention. To do this, we set the following objectives:

  1. Develop a generalized form of the thermal Stokes number to normalize the retention probability of particles on a vane.

  2. Propose a new retention probability to describe the likelihood of particle sticking for a range of engine operating states.

  3. Propose a reduced-order function to predict the retention probability of a particles on a nozzle guide vane.

We perform a numerical simulation on the two-dimensional mid-span section of the nozzle guide vane from the General Electric Energy Efficient Engine (GE-E3) using the CFD code ansys fluent. This code has been chosen due to its ability to allow user customization for the dispersed phase calculations, particularly the particle-wall interaction. More detail regarding the GE-E3 high pressure turbine and its design can be found in Ref. [29]. For each case tested, the continuous phase is solved first in order to provide a converged solution for the fluid flow into which the particles can be introduced. The particle tracking algorithm is then exercised multiple times using the same converged continuous phase solution. A particle fate model is built into these simulations through the use of a user-defined function, which applies a condition for retention during the particle-wall interaction.

4.1 Computational Domain and Mesh.

We have adopted the domain and boundary conditions from the two-dimensional mid-span section of the GE-E3 NGV used in Ref. [29]. The inlet and outlet boundaries are spaced one vane chord length from the leading and trailing edges of the vane, respectively. The vane surface is specified as uniformly smooth with a no-slip condition, and the fluid-wall interface is treated as adiabatic. While this is not representative of the cooled wall surfaces used on turbine vanes, in this study, we have restricted our attention to the effects of the free-stream temperature on the probability of retention only. The influence of vane temperature is the subject of ongoing work.

The computational mesh was generated by the ansys icem cfd meshing software. The resulting grid uses tetrahedral elements and is unstructured in the far field, with a structured 50 layer prism mesh used to resolve the near wall boundary layer flow as shown in Fig. 3. The prism mesh was grown at a specified rate in order to ensure a wall y+ less than unity around the entire surface of the vane for all converged solutions. This eliminated the need to apply wall functions in the turbulent boundary layer as this was always fully resolved. A grid independence study was carried out by repeatedly refining the mesh density throughout the domain and monitoring both the average skin friction coefficient and the skin friction profile around the vane surface. Once there was less than 0.01% change between successive converged solutions, the flow computation was considered to be grid independent. The total cell count for the final computational grid was 298,772 with a solution time for the continuous phase of 1.5 CPU hours.

Fig. 3
The vane leading edge, showing the unstructured far-field mesh and the structured near wall prism mesh
Fig. 3
The vane leading edge, showing the unstructured far-field mesh and the structured near wall prism mesh

4.2 Continuous Phase Modeling.

Pressure and velocity fields around the vane were solved using the Navier–Stokes equations for mass and momentum for each continuous phase set of boundary conditions. These solutions were carried out in a coupled manner with the conservation equations solved simultaneously using a second-order least squares cell-based discretization. Variations in density were accounted for by considering the flow to be an ideal gas. Due to this assumption and the limited change in the fluid temperature through the domain, the fluid viscosity, thermal conductivity, and heat capacity were assumed constant for the continuous phase temperature. The temperature distribution around the vane was determined by solving the energy equation and was then used to calculate the local enthalpy using the equation of state. Turbulent variations in the local velocity field were evaluated by a second order discretization of Menter’s shear stress transport (SST) k–ω form of the Reynolds-averaged Navier-Stokes (RANS) equations. This method was selected due to its good suitability for flows with large positive pressure gradients such as the vane stage of a high pressure turbine. Turbulence intensity and integral length scale were assumed constant throughout, with values of 5% and 2.216 cm chosen to match those determined experimentally by Whitaker et al. [31].

Continuous phase solutions were considered to be converged once the residuals for each conservation equation had reduced by at least five orders of magnitude. The average skin friction, lift, and drag coefficients of the vane were also monitored, and it was ensured that once converged, these variables showed less than 0.001% change over 10 solution iterations. The continuous phase has been validated against the reference data set of Timko [32] for a continuous phase Reynolds number of 176,480 [29].

4.3 Dispersed Phase Modeling.

The dispersed phase particle tracking computations were carried out in a Lagrangian framework using the discrete phase model (DPM) provided with ansys fluent. Having independently solved the flow-field, spherical particles with constant density were injected from 500 points spaced equally across the domain inlet. An injector independence study was also carried out. All particles were introduced from rest and in thermal equilibrium with the continuous phase. The exchange of momentum and energy between the dispersed and continuous phases was assumed to be one-way coupled, with only the effect of the continuous phase on the dispersed phase accounted for. This is an appropriate assumption to make as the maximum particle volume fraction tested of 4.903 × 10−8 was sufficiently low [20].

The trajectory followed by a particle is computed by stepwise integration at a series of discrete time-steps. At each time-step, the motion variables were used to compute new forces. The resulting acceleration from these forces along with the current position and velocity of the particle was then used to solve a set of coupled ordinary differential equations for the new particle velocity and position. The solution of this coupled system was carried out using an implicit trapezoidal discretization of the differential equations for the velocity and position changes of the particle. The time-step used for the trajectory integration was calculated based upon a specified step length factor, which was defined as equal to the maximum number of time-steps allowed divided by the number of mesh cells in the main flow direction, in this case giving a length factor of 50.

The presiding force driving the particle motion was the drag force, this being dominant for particles greater than 1 μm in diameter. The pressure gradient force term was also included due to the non-zero pressure gradient across the domain. The thermophoretic force that arises due to the production of a temperature gradient across the particle was also included. This is particularly important for particles of the order 1 μm in diameter which traverse the vane boundary layer where thermal gradients are the steepest. The gravity force was neglected, as was the virtual mass force due to the large ratio between the dispersed and continuous phase densities. Turbulent dispersion was modeled using the discrete random walk model with 25 particles injected from each injection point to ensure statistically significant results. The model subjects each particle to a Gaussian distributed random velocity fluctuation with a specified time scale constant, in this case defined to be 0.15. The components of this random velocity fluctuation are added to the local Reynolds-averaged velocities giving the instantaneous particle velocities used in the trajectory calculation. As a result, particles do not necessarily enter the domain normal to the inlet boundary.

Heat transfer to and from the particle is governed by Eq. (3). The effect of radiative heat transfer has been neglected in this analysis. The particle temperature was calculated in the same manner as the velocity and position, using implicit trapezoidal integration of Eq. (3) and the temperature of the particle in the previous time-step. To carry out this integration, the particle Nusselt number and therefore the convective heat transfer coefficient must be specified. This was implemented through the use of the Ranz–Marshall correlation.

4.4 EBFOG Fate Model.

The particle fate model applied in this work is the EBFOG model proposed by Casari et al. [24]. The model assumes that the probability of a particle sticking Sp can be described using an Arrhenius-type equation. This type of equation is used to describe the temperature dependence of a wide range of processes in which an energy threshold must be overcome to initiate the process. Casari et al. applied it to the particle deposition process under the assumption that only particles with an energy greater than some threshold value will stick on impact [24].

In this specific application, the probability of sticking Sp depends on a rate constant A and the ratio between an activation energy Eact and a reference energy Eref of the particle. The activation energy defines the energy threshold of the particle-wall adhesion process based upon the physical state of the particle and the reference energy is specified as the kinetic energy of the particle normal to the vane. The sticking probability can therefore be expressed mathematically in the form of Eq. (7).
Sp=Ae(Eact/Eref)
(7)
The ratio between the reference energy Eref and activation energy Eact for a given particle is defined by two empirical constants C1 and C2, the wall normal kinetic energy of the particle, including the particle mass mp and the wall normal particle velocity vp,n2 and a Taylor series expansion of the ratio between the fluid temperature Tf and characteristic transition temperature T* truncated at the first order:
EactEref=C10.5mpvp,n2(1+C2(TfT*))
(8)
Casari et al. [24] showed that C2 is a universal non-dimensional constant equal to 3027 that applies for all particle types. The constants A (Eq. (7)) and C1 (Eq. (8)) both depend upon the composition of the dust and can only be determined by curve fitting experimental results for different particle types. Whether these constants also depend on the particle size is currently unknown. This fitting process was carried out by Casari et al. [24] for the volcanic ash particles, giving C1 = 2.51 × 10−5 J and A = 0.897.

It is through the ratio Tf/T* in Eq. (8) that the temperature dependence of deposition is modeled. The physical state of the particle is associated with the term T*, defined as a characteristic transition temperature of the particle. It is suggested that for amorphous ash particles, this is the glass transition temperature, and the melting temperature for crystalline materials [24]. As the transition temperature of the material reduces, the ratio of the activation energy to the reference energy also reduces, making sticking more probable. The diameter dependence which is key to the OSU model is captured by including the wall normal kinetic energy in the energy ratio calculation (Eq. (8)).

EBFOG has been implemented into our CFD simulations by combining it with a continuous random walk process to emulate the stochastic nature of particle adhesion. To do this, a Metropolis-Hasting algorithm is used [24]. First, a random number in the range 0–1 is generated, against which the calculated sticking probability is compared. If the sticking probability exceeds the randomly generated number then the particle is assumed to stick. When the converse is true, the particle rebounds. In the current work, we do not model deposit growth and its influence on the flow field. This would require a dynamic mesh morphing approach that has been applied previously with some success for computationally simple internal cooling geometries [33]. However, use of this approach for external deposition on an NGV would make production of the data necessary for this study prohibitively expensive.

4.5 Test Matrix.

Continuous phase solutions corresponding to four different bulk temperatures over the range 1400–1700 K have been conducted. We limit our analysis to this range due to the availability of experimentally determined accumulation factors for these temperatures. A summary of the boundary conditions for these cases is provided in Table 2. The continuous phase static pressure gradient was assumed to be constant for all, with the outlet pressure and temperature prescribed to maintain a constant stage reaction for the given turbine inlet temperature.

Table 2

Continuous phase boundary conditions and prescribed fluid properties for the four temperatures tested

Vane inlet temperature, TT,in K]
1400150016001700
PS,in [Pa]179263
TT,out [K]1345144115371630
PS,out [Pa]98595
μf [kg/ms]5.06E-55.26E-55.46E-55.65E-5
cp,f [J/kgK]1200.501211.201220.701229.30
kf [W/mK]8.45E-28.83E-29.20E-29.55E-2
Vane inlet temperature, TT,in K]
1400150016001700
PS,in [Pa]179263
TT,out [K]1345144115371630
PS,out [Pa]98595
μf [kg/ms]5.06E-55.26E-55.46E-55.65E-5
cp,f [J/kgK]1200.501211.201220.701229.30
kf [W/mK]8.45E-28.83E-29.20E-29.55E-2

The dispersed phase was based on the widely studied JBPS sub-bituminous ash, for which the elastic and thermal properties are well documented [6]. All particles have a bulk density of 2320 kg/m3, thermal conductivity of 0.5 W/m K, and specific heat capacity 984 J/kg K. We have also limited our analysis to spherical particles, but expect that particle shape will influence the probability of retention due to the increasing surface area to volume ratio as particles become less spherical. This is expected to enhance the heat transfer rate to the particle and thus the rate at which it reaches its softening temperature and future work should aim to include this in the approach.

For the dispersed phase calculations, a logarithmically spaced set of 30 target particle Reynolds numbers were defined over five orders of magnitude (0.1–1000). For each of these, a corresponding particle diameter was calculated using the bulk inlet velocity from the converged continuous phase calculations, average density, and the viscosity of the gas. As particles were injected from rest, the particle velocity at the inlet is equal to zero, and therefore, a corresponding diameter can be determined. A sub-set of these calculations is shown in Table 3 for continuous phase temperatures of 1400 K and 1600 K. The diameters obtained for each temperature were allowed to vary in order to match the target particle Reynolds number while the flow density and viscosity vary with fluid temperature.

Table 3

Test matrix showing a sub-set of the discrete phase diameters tested for each continuous phase temperature

Temperature, TT,in [K]
14001600
Density, ρf[kg/m3]0.640.56
Viscosity, μf[Pa.s]4.98 × 10−55.38 × 10−5
Rep,in [ − ]dp [μm]
0.110.080.94
0.530.390.45
1.080.800.93
2.211.641.90
4.523.353.88
8.266.137.08
11.728.6910.1
23.9517.820.5
48.9436.342.0
80.0559.468.7
400.00297343
600.00445515
Temperature, TT,in [K]
14001600
Density, ρf[kg/m3]0.640.56
Viscosity, μf[Pa.s]4.98 × 10−55.38 × 10−5
Rep,in [ − ]dp [μm]
0.110.080.94
0.530.390.45
1.080.800.93
2.211.641.90
4.523.353.88
8.266.137.08
11.728.6910.1
23.9517.820.5
48.9436.342.0
80.0559.468.7
400.00297343
600.00445515

4.6 Interaction and Retention Probability.

To determine the retention probability for a given continuous phase temperature and particle Reynolds number, 25 identical particles were injected from the 500 injection locations and their fates upon impact determined by application of the particle fate model described in Sec. 4.4. Retention depends first upon whether the particle hits the vane, and second upon whether it meets the criteria for sticking. The former is dealt with by the interaction probability defined in Ref. [29] and given by Eq. (9). The probability of interaction is a function of the inertial response of the particle to changes in the core gas flow and depends on multiple factors including particle shape and density. Bojdo et al. [29] showed that the number of particles interacting with the vane compared to the number injected into the domain (ninteracted and ninjected) can be expressed in terms of particle-fluid system properties through the particle generalized Stokes number, Stkgen.
ηint=ninteractedninjected=a¯tanh(b¯Stkgenc¯)+d¯
(9)
The use of this generalized parameter produces an interaction probability function (Eq. (9)), which depends on particle shape only. The particle shape has the effect of modifying the drag coefficient for a given particle Reynolds number and thus changes the trajectory followed by the particle around the vane, modifying its probability of interaction. Therefore, the coefficients a¯, b¯, c¯, and d¯ in Eq. (9) depend on the sphericity of the particle concerned. This implies that only a fraction of the particles introduced into the domain will hit the vane, and only a proportion of that fraction will go on to be retained. Thus, we define the retention probability in Eq. (10) as the ratio of the number of particles retained, nretained to the number that interacted with the vane ninteracted. In this work, we only consider particles which are retained upon their first impact with the NGV. For simplicity, we neglect particles that collide for a second time after initially failing to meet the criteria for sticking. Having lost some energy during their initial contact with the vane, such particles may not possess sufficient kinetic energy to deform significantly enough to meet the sticking criteria during a second interaction.
ηret(Stkth)=nretainedninteracted
(10)

5 Discussion

5.1 Generalized Thermal Stokes Number.

As an alternative to the forms given in Eq. (6), the particle thermal Stokes number can be expressed as a ratio of response times, neglecting the Nusselt number or heat transfer coefficient. Instead, the particle thermal response time, τth, is normalized by the velocity response time of the fluid to the presence of the NGV, τf, where the characteristic dimension is chosen to be the throat spacing, Lth, resulting in Eq. (11).
Stth=τthτf=ρpcp,pdp2Uin12kfLth
(11)
The two forms of thermal Stokes number presented in Sec. 3.2 contain the convective heat transfer coefficient and its non-dimensional form the Nusselt number. These are therefore functions of the particle relative Reynolds number and thus depend on the particle trajectory through the NGV passage. From the viewpoint of low-order modeling, requiring each trajectory to be solved becomes a time-consuming process and negates the benefits of fast running generalized models. While the Eq. (11) form of thermal Stokes number eliminates the trajectory-dependent Nusselt number and heat transfer coefficient, it does not provide any indication of the particle’s state on impact. For instance, there is no indication that the particle has exceeded the softening temperature and therefore undergone a considerable reduction in its elastic properties. Given that a softened particle is more likely to deform plastically and adhere, it is logical that the parameter used to characterize the probability of retention contains some reference to this phenomenon.

In seeking a generalized model, we carry out a parameter analysis of the variables contained within the EBFOG model. We choose EBFOG over other particle fate models due to the fact it is already a generalized model for the sticking probability and requires only a few experimentally determined constants in order to be applied. By contrast, higher fidelity models such as that of Bons et al. require an understanding of how the mechanical properties of the particle vary with local fluid temperature [8]. This makes the sticking probability predicted by this model a function of the trajectory followed by the particle through the temperature gradients around the vane. As a result, any generalized thermal Stokes number derived from a model containing these relationships would have to account for this variation of mechanical properties with temperature, something which is currently not well known for any particles undergoing plastic deformation. By contrast the EBFOG model incorporates the effect of mechanical property changes with temperature into the model through the constants A and C1, which Casari et al. [24] have shown can be easily obtained by curve fitting to experimental data.

In this parameter sensitivity analysis, the retention probability corresponding to a range of thermal Stokes numbers has been calculated for all parameters in Eq. (11) along with the EBFOG model constants, A and C1 and the material softening temperature T*. The retention probabilities obtained by varying each of these parameters are shown in Fig. 4. As we expect, the probability of retention increases with increasing freestream temperature, reducing particle transition temperature and heat capacity. The particle density however appears to have a relatively limited effect on the retention probability. The effect of the EBFOG constants are also evident, with an increasing activation energy C1 or reduced maximum retention probability A resulting in a lower probability of retention. It is clear from Fig. 4 that there is a large variation in the number of particles retained, even when different flow conditions and particle properties produce the same thermal Stokes number. This limits the general applicability of the thermal Stokes number in its current form for predicting the retention probability. We therefore seek to obtain a new, reduced-order parameter that can predict the retention probability for different flow conditions and particle properties.

Fig. 4
Parameter study of particle, fluid, and EBFOG model variables showing the range of retention probabilities that can be obtained from varying different parameters. The particle reference values used were based upon JBPS fly-ash: ρp=2320kg/m3, T* = 1197 K, cp,p = 984 J/kgK, and Laki volcanic ash: A = 0.897, C1 = 2.51 × 10−6 J. The reference fluid temperature was Tf = 1600 K.
Fig. 4
Parameter study of particle, fluid, and EBFOG model variables showing the range of retention probabilities that can be obtained from varying different parameters. The particle reference values used were based upon JBPS fly-ash: ρp=2320kg/m3, T* = 1197 K, cp,p = 984 J/kgK, and Laki volcanic ash: A = 0.897, C1 = 2.51 × 10−6 J. The reference fluid temperature was Tf = 1600 K.
We begin by considering the particle heat capacity cp,p and adopt a reference value for a given particle type, this being the heat capacity at standard conditions. Something widely available in the literature for many particle types. Figure 4 shows how the retention probability changes for cases in which the fluid temperature, hence thermal conductivity and the particle density are varied independently. These are all Reynolds number effects due to the variability of the heat transfer coefficient with particle Reynolds number. Inertial Reynolds number effects have previously been accommodated using the generalized Stokes number to normalize the interaction probability [29]. We therefore include the non-Stokes drag correction factor ψ as a means of normalizing the heat transfer Reynolds number effects in the generalized thermal Stokes number. The final parameter that results in variation of the retention probability is the transition temperature of the particle T*. We accommodate this in the generalized thermal Stokes number by including its ratio with the gas temperature as in Eq. (7). The resulting generalized thermal Stokes number can therefore be expressed as the product of the ratio between the particle thermal response time and the fluid velocity response time (Eq. (11)) expressed in terms of engine design and operating state parameters, the non-Stokes drag correction factor (ψ) and the temperature ratio (Tf/T*).
Stkth,gen=τthτf=cp,pdp2ρpW112kfNhρf(TfT*)ψ
(12)

5.2 Retention Probability.

The next step in the production of a generalized retention model is to account for the effect that the constant A has in the EBFOG model. As the exponent in Eq. (7) is negative, this term will always be less than unity and therefore A represents the maximum achievable sticking probability. This constant depends on the mechanical properties of the particles in the dust which in turn depend upon the particle composition; thus, there is an indirect dependence of A on the composition of the dust. As with the activation energy constant C1, A can only be determined by curve fitting experimental data. By definition, the maximum achievable probability of retention that can be obtained is unity. To account for this in the model and normalize all curves to have retention probabilities spanning the range 0–1, we multiply the retention probability by the reciprocal of A. As such, the equation to calculate the retention probability from a sigmoidal curve fit function becomes
ηret=1A[a+ba[1+(Stkth,genc)d]e]
(13)
When we apply the correction factor to the results shown in Fig. 4, neglecting the results generated by varying C1 and cp,p, we obtain the result in Fig. 5 in which the retention probabilities for each thermal Stokes number collapse onto a single line. This can be described by a single sigmoidal curve fit function for which the coefficients are given in Table 4. This curve is valid for all freestream temperatures, particle densities, and particle transition temperatures. The curve fit is geometry-specific to the GE-E3 vane and is dependent on the model constant C1 and the particle heat capacity, cp,p. Thus, we can expect different fitting coefficients for different particle types. However, with suitable experimental data to provide the constants A and C1, it is possible to use this technique to develop a generalized curve for each material (combination of cp,p, A, C1) which is valid for all continuous phase conditions and therefore a range of engine operating states. The constant a in Eq. (13) defines the maximum sticking probability for a given dust and therefore is equivalent to the EBFOG model constant A.
Fig. 5
Results from Fig. 4 normalized using the generalized thermal Stokes number and generalized retention probability
Fig. 5
Results from Fig. 4 normalized using the generalized thermal Stokes number and generalized retention probability
Table 4

Sigmoidal curve fit parameters for the generalized retention probability function (Eq. (13))

Particulateabcde
JBPSA0.000f(C1)f(C1)f(C1)
Particulateabcde
JBPSA0.000f(C1)f(C1)f(C1)

Note: The functions that can be used to calculate coefficients c, d, and e are given by Eqs. (14)(16), respectively.

Finally, we seek to generalize the effect that the EBFOG activation energy constant, C1 has on the probability of retention. As previously discussed, this is a material-dependent constant and will therefore allow the sticking probability predicted by the EBFOG model to vary for different particulate types. Currently, its applicability to mineral dusts is not yet proven due to a lack of experimental data. In order to account for this in our reduced-order model, we define the coefficients c, d, and e in the sigmoidal function describing the relationship between generalized thermal Stokes number and retention probability as functions of the constant C1. These functions are given by Eqs. (14)(16), respectively, and can be found for any value of C1 which is determined experimentally (see Ref. [24]).
c=7.85×107C12+5.86×104C1+7.89×101
(14)
d=0.984C10.114
(15)
e=3.400×103ln(C1)+4.101×101
(16)

When appropriate experimental accumulation factor data which are suitable for deriving the constants A and C1 become available, these coefficients can be simply predicted and the reduced-order retention probability model applied for the given dust.

5.3 Validation.

Data that can be used to validate the retention probability model are scarce for all particulate types and geometries. There is not a well-defined validation case for which all properties of the particle-vane-fluid system are known. However, if a particle size distribution and sufficient detail regarding the continuous and dispersed phase properties is available, then results from the retention probability model can be reasonably validated.

Data from the experimental deposition studies by Whitaker et al. [31] and Bonilla et al. [34] of JBPS coal fly-ash on a set of CFM-56 NGVs using the Ohio State University TuRFR rig, shown in Table 5, allow a comparison to be made with accumulation factor results from the reduced-order model for these studies. Each study used a range of size distributions with different mean diameters, as noted in Table 5. We take the quoted particle size distributions used in each experiment and determine generalized inertial and thermal Stokes numbers for each discrete particle diameter in the distribution. This allows interaction and retention probabilities corresponding to each inertial and thermal Stokes number to be determined. The interaction probability is calculated using Eq. (9) with the curve fit parameters from Ref. [29] and the Sacco et al. [35] form of the momentum Stokes number. The retention probability is calculated based on the curve fit parameters in Table 4. These probabilities are multiplied to give the accumulation factor for each particle diameter (ζmodel).

Table 5

Predicted values of the interaction probability, retention probability, and accumulation factor from the generalized particle retention model compared with experimentally determined accumulation factors from Whitaker et al. [31] and Bonilla et al. [34]

TfUindp¯Rep,in¯Stkgen¯Stkth,gen¯ηintηretζmodelζexp.
Ref.[K][m/s][μm][-][-][-][-][-][%][%]
Whitaker et al. [31]1400704.631.620.270.280.350.00033.26%1.90–3.20%
1400706.482.260.520.540.750.00264.54%2.50–4.00%
Bonilla et al. [34]1400703.801.330.190.190.210.00012.55%6.80%
1400706.402.230.510.520.730.00243.77%7.50%
14007011.84.121.611.671.000.099128.40%10.70%
14007015.85.512.782.871.000.333930.50%25.00%
TfUindp¯Rep,in¯Stkgen¯Stkth,gen¯ηintηretζmodelζexp.
Ref.[K][m/s][μm][-][-][-][-][-][%][%]
Whitaker et al. [31]1400704.631.620.270.280.350.00033.26%1.90–3.20%
1400706.482.260.520.540.750.00264.54%2.50–4.00%
Bonilla et al. [34]1400703.801.330.190.190.210.00012.55%6.80%
1400706.402.230.510.520.730.00243.77%7.50%
14007011.84.121.611.671.000.099128.40%10.70%
14007015.85.512.782.871.000.333930.50%25.00%

Comparing experimental results of Whitaker et al. [31] with the retention probability model, it is clear that the model tends to predict the accumulation factor well for size distributions with mean diameters less than 10 μm. For the experimental deposition results of Bonilla at al. [34], we see that the retention probability model under-predicts the accumulation factor for relatively small mean diameters and then over-predicts as the diameter increases. An over-prediction of the accumulation factor is achieved for size distributions with mean particle diameters greater than 10 μm. Given that the EBFOG particle fate model tends to predict low retention probabilities for low thermal Stokes numbers (Fig. 5) and larger retention probabilities as the thermal Stokes number increases, this is an understandable result. The over-prediction of the accumulation factor for relatively large average diameters is in part due to the fact that EBFOG does not account for particle removal, a mechanism whose significance increases with particle diameter.

The under-prediction of accumulation factor for size distributions with sub-10 μm average diameters in Table 5 may be explained by some of the second-order effects which are not captured in the reduced-order approach. These include the effects of vane surface condition and particle non-sphericity on particle-fluid heat transfer, melting rate and thus retention probability. We would expect the ashes used in the experiments of Bonilla et al. [34] and Whitaker et al. [31] to exhibit a range of particle shapes and for constituent particles to melt more rapidly than their spherical counterparts due to a reduced thermal response time. Thus, for a fixed equivalent particle diameter, we would expect an increased generalized Stokes number as the transition temperature T* reduces (Eq. (12)). With reference to Fig. 5 this would result in a larger probability of retention for the non-spherical particle and an increase in the accumulation factor of the bulk ash. Incorporating these second-order effects into the reduced-order approach should therefore be a priority of future work.

In the context of predicting engine performance loss, assessing the required accuracy of the reduced-order model relies on understanding the sensitivity of engine performance degradation to the accumulation factor. An estimate of the sensitivity to this parameter can be made using the Clarkson and Simpson [10] approach for the deposition rate of ash on an NGV (Eq. (1)). It is thought that a critical deposit mass exists, which will ultimately result in a surge event. Given that the rate at which deposit mass builds up is directly proportional to the accumulation factor, we expect a 10% over-prediction of accumulation factor would result in a 10% under-prediction of the ash dose (see Sec. 2.1) required to surge the engine. Determining an acceptable level of accuracy in the accumulation factor prediction is therefore contingent on an improved understanding of the sensitivity between engine performance degradation and mass deposition rate.

5.4 Application.

The dusts ingested by gas turbine engines contain a range of particle sizes that can be represented in a particle size distribution by mass, volume or number. For each of these, the proportions of particles within a given discrete size bin are represented as fractions of the total mass, volume or number of particles. We now demonstrate how the reduced-order model can be applied to approximate the accumulation factor of a hypothetical size distribution by mass of a bulk dust. To do this, we retain the generic log-normal, bi-modal, two-dust distribution used in Ref. [29] for which determination of the interaction factor was demonstrated. The geometric mean and standard deviation (μgeo and σgeo) of these distributions are summarized in Table 6, with 200 logarithmically spaced points in the diameter interval 0.1–100 μm used to define the distribution.

Table 6

Size distribution function geometric mean, standard deviation and mass weightings for the generic bi-modal dust blend

Mineral:ABA + B
μgeoμm1.02.0
σgeoμm0.40.3
ρp,bulkkg/m325003000
Weighting(–)40.0%60.0%100.0%
χNGV(–)4.45%40.6%45.0%
ζNGV(–)1.75%29.6%31.6%
ξNGV,max(–)10.0%42.3%52.3%
Mineral:ABA + B
μgeoμm1.02.0
σgeoμm0.40.3
ρp,bulkkg/m325003000
Weighting(–)40.0%60.0%100.0%
χNGV(–)4.45%40.6%45.0%
ζNGV(–)1.75%29.6%31.6%
ξNGV,max(–)10.0%42.3%52.3%

Note: The corresponding interaction, retention, and accumulation factors for the individual dusts and the combined blend are shown.

We can extend the analysis carried out in Ref. [29] to calculate the accumulation factor of the same bulk dust. To do this, we first re-cast the interacted distribution in terms of the generalized thermal Stokes number as defined by Eq. (12) for a given set of engine operating conditions. We assume that the EBFOG constants A and C1 are 0.897 and 2.51 × 10−5, respectively, while T* is defined to be 1132 K as defined by Suman et al. [36] for Laki vent volcanic ash. For each discrete particle size in the interacted distribution we then calculate the corresponding probability of retention using Eq. (13). This is analogous to applying the retention probability function (Eq. (13)) across the whole particle size distribution and multiplying the interaction and retention probabilities for each discrete diameter band. In this case, the retention probability obtained represents the maximum amount of particulate that would be retained, if all particles were assumed to interact with the vane.

The final step is to determine the bulk proportion of the size distribution which is ultimately retained by the vane. To do this, we simply multiply the corresponding retention probability for a given particle size by its corresponding mass fraction in the size distribution (nd(dp)). Integrating under the resulting retained curve and dividing by the corresponding area under the particle size distribution represents the accumulation factor, (ζNGV). This is represented by the retained shaded area in Fig. 6 for the hypothetical dust distribution. In this figure, the retained area represents the accumulation factor, ζNGV, of the dust under the given set of conditions. Thus, we can mathematically describe the accumulation factor of a size distribution under a particular set of engine operating conditions to be
ζNGV=dlowerduppernintηretd(dp)
(17)
where nint represents the mass fraction of a particular diameter band which interacts with the vane.
Fig. 6
Normalized fraction by mass of a hypothetical bi-modal dust distribution that accumulates on the GE-E3 NGV stage at design point conditions specified in Ref. [29]
Fig. 6
Normalized fraction by mass of a hypothetical bi-modal dust distribution that accumulates on the GE-E3 NGV stage at design point conditions specified in Ref. [29]

For the hypothetical distribution, we can also define the corresponding interaction (χNGV), accumulation (ζNGV), and maximum retention factors (ξNGV,max) for the two constituents of the bi-modal size distribution. The maximum retention factor represents the worst case scenario of accumulation factor for the given dust, that is, what the accumulation factor would be if all dust were to interact with the vane. This is summarized in Table 6.

6 Conclusions

This work bridges the gap between the two extremes of previous studies—empiricism and high fidelity deterministic CFD – by proposing a reduced-order, generalized particle retention model which can be applied for a range of particle types, engine architectures, and operating states. The EBFOG particle retention model has been applied in 2D CFD simulations of particle deposition on nozzle guide vanes and used to derive a universal sigmoidal function. This can be used to describe the retention probability over a range of engine operating conditions and for different particles as a function of a new non-dimensional property of the particle, the generalized thermal Stokes number and a single experimentally determined constant, A from the EBFOG particle fate model. The generalized thermal Stokes number corrects for both the non-linearity of the particle drag force due to changing thermal properties of the fluid, different particle transition temperatures and has been expressed in terms of parameters related to the operating state of the engine and its geometry.

The coefficients of the sigmoidal curve fit function have subsequently been shown to be functions of the empirical constant C1 only. Therefore, if suitable experimental constants A and C1 are determined for different particle types, the coefficients in the sigmoidal retention probability curve can be obtained for each. These can then be used to calculate the retention probability for particles with different generalized thermal Stokes numbers arising from a range of engine operating states, dust composition and engine architecture, giving the model wide applicability. It has been demonstrated how the retention probability as a function of generalized thermal Stokes number can be used with the probability of interaction, providing a two equation model to predict the capture probability of a particle. When applied over a range of particles within a size distribution this allows the accumulation factor of the distribution to be calculated.

Acknowledgment

Part of this work was joint funded by the UK Engineering and Physical Sciences Research Council (EPSRC) and the Defence Science and Technology Laboratory (DSTL), who are funding the doctoral studies of lead author, Matthew Ellis, as part of EPSRC grant EP/P510579/1. The dataset generated during the current study is available in the Mendeley repository, 10.17632/hptxjdjss6.1.

Conflict of Interest

There are no conflicts of interest.

Data Availability Statement

Data provided by a third party listed in Acknowledgment.

Nomenclature

     
  • d =

    diameter

  •  
  • h =

    heat transfer coefficient

  •  
  • k =

    thermal conductivity

  •  
  • m =

    mass

  •  
  • n =

    number of particles

  •  
  • r =

    vane height

  •  
  • t =

    time

  •  
  • u =

    fluid velocity

  •  
  • v =

    particle velocity

  •  
  • A =

    EBFOG maximum sticking probability

  •  
  • E =

    energy

  •  
  • L =

    characteristic length scale

  •  
  • N =

    number of vanes in NGV stage

  •  
  • P =

    pressure

  •  
  • S =

    sticking probability

  •  
  • T =

    temperature

  •  
  • U =

    characteristic fluid velocity

  •  
  • W =

    mass flow rate

  •  
  • c¯ =

    vane chord

  •  
  • C¯ =

    concentration

  •  
  • cp =

    specific heat capacity

  •  
  • C1 =

    EBFOG activation energy constant

  •  
  • C2 =

    EBFOG universal constant

  •  
  • Wa =

    surface energy

  •  
  • T* =

    transition temperature

  •  
  • Nu =

    Nusselt number

  •  
  • Re =

    Reynolds number

  •  
  • Stk =

    Stokes number

  •  
  • ζ =

    accumulation factor

  •  
  • ηcap =

    capture probability

  •  
  • ηint =

    interaction probability

  •  
  • ηret =

    retention probability

  •  
  • μ =

    viscosity

  •  
  • μgeo =

    geometric mean

  •  
  • ξ =

    retention factor

  •  
  • ρ =

    density

  •  
  • σgeo =

    geometric standard deviation

  •  
  • τ =

    response time

  •  
  • χ =

    interaction factor

  •  
  • ψ =

    Non-Stokes drag correction factor

  •  
  • [.]core =

    core flow

  •  
  • [.]crit =

    particle yield point

  •  
  • [.]drag =

    pertaining to the drag force

  •  
  • [.]f =

    fluid

  •  
  • [.]NGV =

    nozzle guide vane

  •  
  • [.]p =

    particle

  •  
  • [.]p,in =

    particle at inlet

  •  
  • [.]p,th =

    particle temperature

  •  
  • [.]p,soft =

    particle softening point

  •  
  • [.]pv =

    particle velocity

  •  
  • [.]0 =

    standard state conditions

  •  
  • [.]S =

    static conditions

  •  
  • [.]T =

    total conditions

  •  
  • [.]throat =

    vane throat

References

1.
Clarkson
,
R. J.
,
Majewicz
,
E. J.
, and
Mack
,
P.
,
2016
, “
A Re-Evaluation of the 2010 Quantitative Understanding of the Effects Volcanic Ash Has on Gas Turbine Engines
,”
Proc. Inst. Mech. Engrs, Part G: J. Aerosp. Eng.
,
230
(
12
), pp.
2274
2291
.
2.
Bojdo
,
N
,
2012
, “
Rotorcraft Engine Air Particle Separation
,”
Ph.D. thesis
,
University of Manchester
,
Manchester, UK
.
3.
Vogel
,
A
,
2018
, “
Volcanic Ash Properties, Atmospheric Properties & Effects on Aero Engines
,”
Ph.D. thesis
,
University of Oslo
,
Oslo
.
4.
Goodwin
,
J. E.
,
Sage
,
W.
, and
Tilly
,
G. P.
,
1969
, “
Study of Erosion by Solid Particles
,”
Proc. Inst. Mech. Engrs
,
184
(
1
), pp.
279
292
.
5.
Dunn
,
M.
,
2012
, “
Operation of Gas Turbine Engines in An Environment Contaminated with Volcanic Ash
,”
ASME J. Turbomach.
,
134
(
5
), p.
051001
.
6.
Ai
,
W.
, and
Fletcher
,
T.
,
2012
, “
Computational Analysis of Conjugate Heat Transfer and Particulate Deposition on a High Pressure Turbine Vane
,”
ASME J. Turbomach.
,
134
(
4
), p.
041020
.
7.
Sreedharan
,
S.
, and
Tafti
,
D.
,
2011
, “
Composition Dependent Model for the Prediction of Syngas Ash Deposition in Turbine Gas Hotpath
,”
Int. J. Heat Fluid Flow
,
32
(
1
), pp.
201
211
.
8.
Bons
,
J. P.
,
Prenter
,
R.
, and
Whitaker
,
S.
,
2017
, “
A Simple Physics-Based Model for Particle Rebound and Deposition in Turbomachinery
,”
ASME J. Turbomach.
,
139
(
8
), p.
081009
.
9.
2019
,
Gas Turbine Engine Environmental Particulate Foreign Object Damage
.
Technical Report, TR-AVT-250, NATO STO
.
10.
Clarkson
,
R.
, and
Simpson
,
H.
,
2017
,
Maximising Airspace Use During Volcanic Eruptions: Matching Engine Durability against Ash Cloud Occurrence
.
Technical Report, STO-MP-AVT-272, NATO STO
.
11.
Giehl
,
C.
,
Brooker
,
R.
,
Marxer
,
H.
, and
Nowak
,
M.
,
2017
, “
An Experimental Simulation of Volcanic Ash Deposition in Gas Turbines and Implications for Jet Engine Safety
,”
Chem. Geol.
,
461
, pp.
160
170
.
12.
Dean
,
J.
,
Taltavull
,
C.
, and
Clyne
,
T.
,
2016
, “
Influence of the Composition and Viscosity of Volcanic Ashes on Their Adhesion Within Gas Turbine Aeroengines
,”
Acta. Mater.
,
109
, pp.
8
16
.
13.
Kim
,
J.
,
Dunn
,
M. G.
,
Baran
,
A. J.
,
Wade
,
D. P.
, and
Tremba
,
E. L.
,
1993
, “
Deposition of Volcanic Materials in the Hot Sections of Two Gas Turbine Engines
,”
ASME J. Eng. Gas Turbines Power
,
115
(
3
), pp.
641
651
.
14.
Bons
,
J.
,
Crosby
,
J.
,
Wammack
,
J.
,
Bentley
,
B.
, and
Fletcher
,
T.
,
2007
, “
High-Pressure Turbine Deposition in Land-Based Gas Turbines From Various Synfuels
,”
ASME J. Eng. Gas Turbines Power
,
129
(
1
), pp.
135
143
.
15.
Crosby
,
J.
,
Lewis
,
S.
,
Bons
,
J.
,
Ai
,
W.
, and
Fletcher
,
T.
,
2008
, “
Effects of Temperature and Particle Size on Deposition in Land Based Turbines
,”
ASME J. Eng. Gas Turbines Power
,
130
(
5
), p.
051503
.
16.
Taltavull
,
C.
,
Dean
,
J.
, and
Clyne
,
T.
,
2016
, “
Adhesion of Volcanic Ash Particles Under Controlled Conditions and Implications for Their Deposition in Gas Turbines
,”
Adv. Eng. Mater.
,
18
(
5
), pp.
803
813
.
17.
Senior
,
C. L.
, and
Srinivasachar
,
S.
,
1995
, “
Viscosity of Ash Particles in Combustion Systems for Prediction of Particle Sticking
,”
Energy Fuels
,
9
(
2
), pp.
277
283
.
18.
Krisak
,
M
,
2015
, “
Environmental Degradation of Nickel-Based Superalloys Due to Gypsiferous Desert Dusts
,”
Ph.D. thesis
,
Air Force Institute of Technology
,
Wright-Patterson AFB, OH
.
19.
Brach
,
R.
, and
Dunn
,
P.
,
1992
, “
A Mathematical Model of the Impact and Adhesion of Microspheres
,”
Aerosol. Sci. Technol.
,
16
(
1
), pp.
51
64
.
20.
El-Batsh
,
M. H
,
2001
, “
Modeling Particle Deposition on Compressor and Turbine Blade Surfaces
,”
Ph.D. thesis
,
Vienna University of Technology
,
Vienna
.
21.
Soltani
,
M.
, and
Ahmadi
,
G.
,
1994
, “
On Particle Adhesion and Removal Mechanisms in Turbulent Flows
,”
J. Adhes. Sci. Technol.
,
8
(
7
), pp.
763
785
.
22.
Dunn
,
M.
,
Baran
,
A.
, and
Miatech
,
J.
,
1996
, “
Operation of Gas Turbine Engines in Volcanic Ash Clouds
,”
ASME J. Eng. Gas Turbines Power
,
118
(
4
), pp.
724
731
.
23.
Singh
,
S.
, and
Tafti
,
D.
,
2015
, “
Particle Deposition Model for Particulate Flows At High Temperatures in Gas Turbine Components
,”
Int. J. Heat Fluid Flow
,
52
, pp.
72
83
.
24.
Casari
,
N.
,
Pinelli
,
M.
,
Suman
,
A.
,
Di Mare
,
L.
, and
Montomoli
,
F.
,
2017
, “
An Energy-Based Fouling Model for Gas Turbines: EBFOG
,”
ASME J. Turbomach.
,
139
(
2
), p.
021002
.
25.
Casari
,
N.
,
Pinelli
,
M.
,
Suman
,
A.
,
Montomoli
,
F.
, and
di Mare
,
L.
,
2017
, “
EBFOG: Deposition, Erosion and Detachment on High Pressure Turbine Vanes
,”
ASME J. Turbomach.
,
140
(
6
), p.
061001
.
26.
Crowe
,
C.
,
Sommerfeld
,
M.
, and
Tsuji
,
Y.
,
1998
,
Multiphase Flows with Droplets and Particles
, 1st ed.,
CRC Press
,
Boca Raton, FL
.
27.
Singh
,
S.
, and
Tafti
,
D.
,
2016
, “
Prediction of Sand Transport and Deposition in a Two-Pass Internal Cooling Duct
,”
ASME J. Eng. Gas Turbines Power
,
138
(
7
), p.
072606
.
28.
Bojdo
,
N.
, and
Filippone
,
A.
,
2019
, “
A Simple Model to Assess the Role of Dust Composition and Size on Deposition in Rotorcraft Engines
,”
Aerospace
,
6
(
4
), p.
44
.
29.
Bojdo
,
N.
,
Ellis
,
M.
,
Filippone
,
A.
,
Jones
,
M.
, and
Pawley
,
A.
,
2019
, “
Particle-Vane Interaction Probability in Gas Turbine Engines
,”
ASME J. Turbomach.
,
141
(
9
), p.
091010
.
30.
Israel
,
R.
, and
Rosner
,
D. E.
,
1982
, “
Use of a Generalized Stokes Number to Determine the Aerodynamic Capture Efficiency of Non-Stokesian Particles From a Compressible Gas Flow
,”
Aerosol. Sci. Technol.
,
2
(
1
), pp.
45
51
.
31.
Whitaker
,
S. M.
,
Prenter
,
R.
, and
Bons
,
J. P.
,
2015
, “
The Effect of Freestream Turbulence on Deposition for Nozzle Guide Vanes
,”
ASME J. Turbomach.
,
137
(
12
), p.
121001
.
32.
Timko
,
L.
,
1984
,
Energy Efficient Engine High Pressure Turbine Component Test Performance report
.
NASA Technical Report NASA-CR-168289, NASA
.
33.
Forsyth
,
P.
,
Gillespie
,
D.
, and
McGilvray
,
M.
,
2018
, “
Development and Applications of a Coupled Particle Deposition–Dynamic Mesh Morphing Approach for the Numerical Simulation of Gas Turbine Flows
,”
ASME J. Eng. Gas Turbines Power
,
140
(
2
), p.
022603
.
34.
Bonilla
,
C.
,
Webb
,
J.
,
Clum
,
C.
,
Casaday
,
B.
,
Brewer
,
E.
, and
Bons
,
J. P.
,
2012
, “
The Effect of Particle Size and Film Cooling on Nozzle Guide Vane Deposition
,”
ASME J. Eng. Gas Turbines Power
,
134
(
10
), p.
101901
.
35.
Sacco
,
C.
,
Bowen
,
C.
,
Lundgreen
,
R.
,
Bons
,
J. P.
,
Ruggiero
,
E.
,
Allen
,
J.
, and
Bailey
,
J.
,
2018
, “
Dynamic Similarity in Turbine Deposition Testing and the Role of Pressure
,”
ASME J. Eng. Gas Turbines Power
,
140
(
10
), p.
102605
.
36.
Suman
,
A.
,
Casari
,
N.
,
Fabbri
,
E.
, and
Pinelli
,
M.
,
2018
, “
Gas Turbine Fouling Tests: Review, Critical Analysis and Particle Impact Behavior Map
,”
ASME J. Eng. Gas Turbines Power
,
141
(
3
), p.
032601
.