Generalized Predictions of Particle-Vane Retention Probability in Gas Turbine Engines

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⁻³ 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].

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. [6–8]) 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.

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., SiO₂ content) have a significant effect on the accumulation factor.

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.

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.

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.

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-E³) 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-E³ 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.

We have adopted the domain and boundary conditions from the two-dimensional mid-span section of the GE-E³ 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.

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].

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⁻⁸ 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.

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 S_p 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].

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.

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.

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/m³, 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.

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 C₁, 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 C₁ 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 C₁ 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.

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

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$⁠).

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.

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.

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 C₁ are 0.897 and 2.51 × 10⁻⁵, 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.

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.

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 C₁ only. Therefore, if suitable experimental constants A and C₁ 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.

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.

There are no conflicts of interest.

Data provided by a third party listed in Acknowledgment.

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

c_p =

specific heat capacity

C₁ =

EBFOG activation energy constant

C₂ =

EBFOG universal constant

W_a =

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

[.]₀ =

standard state conditions

[.]_S =

static conditions

[.]_T =

total conditions

[.]_throat =

vane throat