Nomenclature

AERD

average erosion rate density (mg/s.mm²)

FPR

fan pressure ratio

HBTFE

high bypass turbofan engine

HEM

hourly eroded mass (mg/hr)

IGV

inlet guide vane

OGV

outlet guide vane

LE

leading edge

PS

pressure side

SS

suction side

TE

trailing edge

Symbols

$A$

area (m²)

c

chord (m)

${C_D}$

drag coefficient

${C_{SL}}$

Saffman lift coefficient

c _p

specific heat at a constant pressure (J/kg.K)

d

diameter of particle (m)

E

erosion rate density (mg/s.mm²)

f

force by mass (N/kg)

g

gravity (m/s²)

h

specific enthalpy (J/kg)

i

turbulence intensity level (%)

k

turbulent kinetic energy (m²/s²)

m

mass (kg)

$\dot{m}$

mass flow rate (kg/s)

$\vec n$ , $\vec t$

normal and tangential unit vectors

$N,\;{N_d}$

rotational speed, design speed (rpm)

p

pressure (Pa)

r

radius (m)

$R{e_p}$

particle Reynolds number

$R{e_s}$

shear flow Reynolds number

T

temperature (K)

$t$

time (s)

$\vec V$

velocity vector (m/s)

Greek symbols

$\beta$

impact angle (deg)

$\epsilon$

local erosion rate (mg/g)

$\varepsilon$

turbulence rate of dissipation (m²/s³)

$\epsilon$

local erosion rate (mg/g)

$\phi$

particle shape factor

ρ

density (kg/m³)

μ

dynamic viscosity (kg/m.s)

$\overrightarrow \Omega $

blade speed of rotation (rad/s)

$\overrightarrow {{\omega _f}} $

fluid rotation (s⁻¹)

$\theta $

tangential component (rad)

Subscripts

abs

absolute

co

corrected

f

fluid

i

inlet

o

outlet

p

particle

r

radial

rel

relative

st

standard atmosphere

$\theta$

tangential

t

total

z

axial

1, 2

impact/rebound

2.0 Fan stage configuration

Current commercial aircraft are powered by HBTFEs, which rely on a large fan as an essential component to achieve the overall engine performance. This fan is characterised by a low aspect ratio and highly loaded thin blades of transonic velocities. The studied components are from an old CF6-80A engine, with typical cruise flight operating parameters listed in Table 1, adapted from [37]. The rotor consists of 38 highly twisted transonic blades (Fig. 1(a)) made from titanium. A splitter located downstream of the fan rotor directs the majority of the air flow to discharge through the bypass duct via 80 OGVs, while the lower sections of the fan rotor are connected to the core of the engine through 78 IGVs, as illustrated in Fig. 1(b). The determination of the geometry of the fan blade and vanes used an optical scanning machine. The CAD models are shown in Fig. 2(a) for the rotor and in Fig. 2(b) for the complete fan stage. Detailed geometric data are listed in Table 2.

3.0 Flow solution

In dilute particulate flows, the gas phase is assumed to be a continuum phase while the particles are dispersed throughout. The flow field in this large transonic fan of low aspect ratio and highly twisted blades is 3-D, compressible and very complex. It involves flow features such as shock formation, boundary layer transition, separation, shock-boundary layer interactions, secondary flows, tip leakage, and vortices formation, as well as significant flow fluctuations. Moreover, the upper sections are prone to flow separation and often operate close to stall flutter in off-design conditions.

Table 1. Operating parameters

Figure 1. Fan stage components: (a) Fan blade (b) OGVs and IGVs.

Figure 2. CAD models: (a) Fan rotor (b) Fan stage components.

Table 2. Geometry parameters

Figure 3. Computational domain.

Although the flow complexity can be high, CFD simulations often rely on Reynolds-Averaged Navier-Stokes (RANS) approximations, since the Direct Numerical Simulation (DNS), capable of accurately representing all scales of turbulence, requires significant computating power. On the other hand, the Large Eddy Simulation (LES) remains prohibitively expensive and not commonly used in major engineering design problems due to computational resource limitations. Despite their inability to capture large separations, RANS models are widely used in industry due their simplicity and ability to achieve reliable results. The k- $\omega$ based SST (shear stress transport) turbulence model [Reference Menter38] combines the robustness of the standard k- $\omega$ model near walls with the capabilities of the k- $\varepsilon$ model, offering a good balance between grid resolution and flow complexity. Moreover, it can predict boundary layer thickness and separations based on wall function, reducing the need for large mesh size and making it a viable option for many industrial scenarios [Reference Menter, Kuntz and Langtry39]. This model has been successfully verified in various instances [Reference Pecnik, Witteveen and Iaccarino40], demonstrating good agreement between simulations and experimental results. Validation tests conducted on the NASA transonic Rotor 35 [Reference Balasubramanian, Barrows and Chen41] demonstrated a high level of agreement between the predictions of the k- $\omega$ based SST turbulence model and LDV data collected at the design and 80% speeds. Another assessment focusing on ultra-high-lift turbine blades [Reference Zhou, Hodson and Himmel42] concluded that the k- $\omega$ SST model accurately predicted aerodynamic performance. In the present research, the flow solution is based on the RANS equations with k- $\omega$ based SST as the turbulence model implemented in the CFX solver.

3.2 Mesh generation

Hexahedral mesh blocks were distributed to secure fine meshing around the fan blade, spinner, vanes, hub and shroud, in addition to the Pitot intake. Near the walls, 12 meshlines guaranteed a minimum volume distortion and to solve for the boundary layers. In the tip clearance, 15 meshlines allowed to capture the leakage flow. The mesh obtained for the fan blade is shown in Fig. 4. The first layer of nodes near walls corresponds to a distance $\;\Delta y$ targeting ${y^ + } = 2$ . The dimensionless parameter $\;{y^ + } = \Delta y\;\frac{\rho }{\mu }{U_t}\;$ considers the friction law of a flat plate $\;{U_t} = {V_\infty }\sqrt {{c_f}/2}$ , where ${c_f} = 0.027Re_c^{ - 1/7}$ . ${V_\infty }$ is the free velocity outside the boundary layer. The Reynolds number $R{e_c} = \rho {V_\infty }c/\mu$ is based on the mean chord of blade or vane. At the cruise operating conditions the values of $R{e_{c\;}}$ of the rotor blade, IGVs and OGVs are equal to 2.5×10⁶, 5.3×10⁵ and 3.3×10⁵, respectively. The turbulence model k- $\omega$ based SST implements an automatic near wall treatment to switch between the low-Reynolds model for finer mesh and the wall-function for coarser mesh. Figure 5 demonstrates that the values of ${y^ + }$ on the surfaces of the fan blade, hub and spinner are below 4.5, with the highest values observed near the LE and blade tip. In the IGVs the values of ${y^ + }$ do not exceed 3.2, and that of OGVs is below 2.9, with a maximum occurring on the SS from the LE.

Figure 4. Rotor mesh: (a) Hub (b) Mid-span (c) Tip (d) Blade.

Figure 5. Contours of y⁺ on the surfaces of: (a) Fan blade (b) IGVs (c) OGVs.

The grid size independence verification carried out in cruise point shows (Fig. 6) that the fan pressure ratio (FPR), bypass ratio (BPR) and total-to-total isentropic efficiency tend to stabilise beyond the fourth mesh. Additionally, the evaluated hourly eroded mass per fan blade for the particle concentration of 2mg/m³ indicates (Fig. 6) slight variation beyond the fourth mesh size. As a result, separate grids of a total number of 1,611,925 nodes and 1,525,900 elements were generated and adopted in this study.

Figure 6. Grid size independence verification.

3.3 Flow results

The flow solutions based on the RANS equations and the stage interface allowed to compute the transonic fan stage performance maps and characterize the complex flow structures through the front components of this HBTFE. Moreover, the flow field provided the necessary data for calculating the trajectories of ash particles and erosion.

3.3.1 Performance maps

The mass averaging of flow properties at the inlet and outlet of the fan stage allowed calculating $FPR = \frac{{{P_{{t_o}}}}}{{{P_{{t_i}}}}}$ and the total-to-total isentropic efficiency ${\eta_{is\;t - t}} = {c_p}_i\;{T_{ti}}\left( {FP{R^{\frac{{{\rm{\gamma }} - 1}}{{\rm{\gamma }}}}} - 1} \right)/\left( {{h_{to}} - {h_{ti}}} \right)$ . These performance parameters are shown in Figs. 7(a) and 7(b), plotted against the corrected mass flow rate ${\dot m_{co}} = \dot m\sqrt {{T_{{t_i}}}/{T_{st}}} /\left( {{P_{{t_i}}}/{P_{st}}} \right)$ and speed ${\dot N_{co}} = N/\sqrt {{T_{{t_i}}}/{T_{st}}}$ , with ${P_{st\;}}$ and ${T_{st}}$ are the static values at standard sea level. The fan rotational speeds are 70, 80, 90, $\;100$ , $102.7$ and 110 ${\rm{\% }}{N_d}$ . From the computed performance (Fig. 7), the design point is determined for ${N_d} = 3,432{\rm{rpm}}$ at the corrected mass flow rate of 637.5kg/s corresponding to the true mass flow rate of 250kg/s. At that point the computed performance gives an FPR = 1.711 and a total-to-total isentropic efficiency $\eta_{is\;t - t} = 82.4$ 1%. Moreover, the cruise operating point corresponds to the rotational speed $\;N = 3,525{\rm{rpm}}$ and a corrected mass flow rate of 654.7kg/s equivalent to the true mass flow rate of 254.15kg/s. At the cruise point, the predicted values are FPR = 1.733, BPR = 4.287, and $\eta_{is\;t - t} = 79.96\% $ . These obtained values align with the data of the operating manual [37]. This latter states that in the flight cruise conditions, the true air mass flow rate is 257.6kg/s equivalent to a corrected mass flow rate of 663.5kg/s, while the FPR = 1.745 and BPR = 4.24.

Figure 7. Fan stage performance: (a) FPR (b) Total-to-total isentropic efficiency.

Figure 8. Vectors of flow velocity at span: (a) 10% (b) 50% (c) 75% (d) 98%.

Figure 9. Streamlines and blade tip flow structure.

Figure 10. Mach number at different blade spans: (a) 10% (b) 50% (c) 75% (d) 98%.

3.3.2 Flow field characterisation

In order to replicate the main flow patterns, the flow field was solved in the cruise flight conditions, Mach 0.8, altitude of 10,650m, and fan speed of $\;102.7{\rm{\% }}{N_d} = 3,525{\rm{rpm}}$ . This large fan is characterised by transonic flow velocities for more than 60% of the blade span. Figure 8 displays the flow velocity vectors at different locations, such as 10% span, mid-span, 75% span and near the tip of blade. Upstream and around the fan blade LE, the relative flow velocities are extremely high and increasing with the blade height. High flow velocities, reaching up to 465.9m/s, are observed on the SS of the fan blade, leading to the formation of shock structure. Towards the blade TE, the lower velocities are due to diffusion and wake formation. At the interface plane between the rotor and IGVs/OGVs, the flow exhibits drastic changes in flow direction, and the blade’s wake disturbs the angle-of-attack around the LE of the stationary vanes. This is clearer at the LE of IGVs since are closer to the fan blade. Across IGVs the flow accelerates, reaching high velocities on the front part of SS. In contrast, the OGVs, located farther, dissipate more the wake issued from the fan blade, hence the flow is well-aligned with the vanes’s LE.

The most intricate flow patterns are observed above the blade tip and through the tip gap (Fig. 9). In this typical region, the leakage flow is clearly depicted by the vectors of flow velocity crossing from the PS to SS, with their directions varying from the LE to TE of the blade. Figure 9 depicts the formation of the leading edge vortex (LEV) from the LE, and extending along the blade tip. The LEV rotates opposite to the tip leakage flow vortex (TLFV) and mixes out with the annulus boundary layer. At the exit from the rotor blade, a swirl velocity is established to ensure a radial equilibrium balanced by the axial velocity profile.

In this transonic fan the pressure loading is most pronounced at the front of the blade, and the diffusion process is greatly influenced by the radial position. The peak of static pressure occurs at the blade LE, while a noticeable depression is visible on the SS starting from mid-span onward. The formation of shock waves is responsible for the main increase in static pressure across the fan blade passage, with the highest aerodynamic loading at the upper sections from the mid-span.

The distributions of the relative Mach number are displayed in Fig. 10 at 10% span, mid-span, 75% span and near blade tip. As the flow moves from the inlet of the Pitot intake towards the fan rotor, the convergence imposed by the conical shape of the spinner caused an increase in Mach number before reaching the fan blade. Near the hub, the relative Mach number is practically subsonic and the flow velocity decreases along the blade passage. A supersonic bubble appears (Fig. 10(a)) at the SS close to the blade LE and extends towards 25% of the span. It is perceived that the transonic flow velocities dominate over more than one-third of the SS, while the subsonic velocities prevail on the remaining of the fan blade. As depicted in Fig. 10(b), a passage shock region develops from the LE over approximately the third of the SS and may reach the opposite fan blade. Downstream of the shock, the Mach number decreases until reaching a value around 0.6 at the exit of the fan blade. The distinctive characteristics of this transonic fan blade include an oblique shock sitting near the LE (Fig. 10(c)) and crossing the blade passage at 50% of the span, which interacts with the neighbouring blade. In the upper sections where the flow is entirely supersonic over the SS, a passage normal shock sits near the TE of the blade (Fig. 10(d)), with a sonic layer extending up to a third or two-thirds of the blade chord on the opposite PS. The strength of this shock is indicated by a jump in the Mach number from a value of 1.43 to 0.83. Near the TE of the fan blade, the boundary layer of the SS becomes considerably thicker compared to that from the PS, and eventually merges to form a thick wake diffusing downstream.

4.0 Particle trajectory

There are two common approaches in particle-laden turbulent flows simulations, which differ in the treatment of the dispersed phase. The Eulerian approach considers both the continuous and dispersed phases as interpenetrating and solves the time-averaged transport equations to predict the field distributions of momentum and volume fraction for each phase [Reference Zhou45]. In the present computations, as the volume fraction at the highest concentration is about 1.5×10⁻⁹ below the limit of 4×10⁻⁷ [Reference Sommerfeld46], the dispersed phase is assumed to be dilute, and the Lagrangian particle tracking model has to be considered. This model is used to track large numbers of particles using a stochastic approach to achieve a significant solution [Reference Chen and Pereira47]. Compared to the Eulerian approach, the Lagrangian approach is more cost-effective [Reference Abbas, Koussa and Lockwood48] and better at handling particle-wall interactions [Reference Shirolkar, Coimbra and Queiroz49]. After obtaining the flow field, the flow data are exported to the trajectory code to compute the particle trajectory and impacts in step by step.

4.1 Particle motion modelling

Particle motion is governed by the balance between the rate of momentum change and the superposition of external forces. The left side of the vector equation of particle motion (1) takes into account inertia, centrifugal force, and Coriolis force, while the right side sums the external forces reduced by the mass of particle, including the drag force $\overrightarrow {\;{f_D}} $ , combined gravity and buoyancy force $\overrightarrow {\;{f_{GB}}} \;$ , pressure gradient force $\overrightarrow {{f_P}} $ , and Saffman force $\overrightarrow {{f_S}} $ .

(1) \begin{align}\left[ {\frac{{{d^2}{r_p}}}{{d{t^2}}} - {r_p}{{\left( {\frac{{d{\theta _p}}}{{dt}} + \omega } \right)}^2}} \right]\overrightarrow {{e_r}} + \left[ {{r_p}\frac{{{d^2}{\theta _p}}}{{d{t^2}}} + 2\frac{{d{r_p}}}{{dt}}\left( {\frac{{d{\theta _p}}}{{dt}} + \omega } \right)} \right]\overrightarrow {{e_\theta }} + \frac{{{d^2}{z_p}}}{{d{t^2}}}\overrightarrow {{e_z}} = \overrightarrow {{f_D}} + \overrightarrow {{f_{GB}}} + \overrightarrow {{f_P}} + \overrightarrow {{f_S}} \end{align}

The dominant external force is the drag force, expressed in terms of the drag coefficient and Reynolds number $R{e_p} = \frac{{{\rho _f}}}{{{\mu _f}}}{d_p}\left| {\overrightarrow {{V_f}} - \overrightarrow {{V_p}} } \right|$ , which when reduced by particle mass becomes:

(2a) \begin{align}\overrightarrow {{f_D}} = \frac{3}{4}\frac{{{\mu _f}R{e_p}}}{{{\rho _{pd_p^2}}}}{C_D}\left( {\overrightarrow {{V_f}} - \overrightarrow {{V_p}} } \right)\end{align}

At very low $R{e_p} \lt 0.1,$ the drag coefficient corresponds to the Stokes regime $\;{C_D} = 24/R{e_p}$ . For the wider range 0.01 $ \le {\rm{R}}{{\rm{e}}_{\rm{p}}} \le $ 2.6×10⁵, the four-parameter drag coefficient [Reference Haider and Levenspiel50] is considered as follows:

(2b) \begin{align}{C_D} = \;\frac{{24}}{{R{e_p}}}\;\left( {1 + aRe_p^b} \right) + \frac{c}{{1 + \frac{d}{{R{e_p}}}}}\end{align}

The constants $\;a$ , $b$ , $c$ and $d$ are revaluated by Haider and Levenspiel [Reference Haider and Levenspiel50] based on the shape factor.

(2c) \begin{align}a &= \exp\left( {2.3288 - 6.4581\phi + {{2.4486}}\phi}^2 \right)\nonumber\\[5pt]b &= 0.0964 + 0.5565\phi\nonumber\\[5pt] c &= \exp\left( {4.905 - 13.8944\phi + {{18.4222}}{\phi}^2 - {{10.2599}}\phi}^3 \right)\nonumber \\[5pt] d &= \exp\left( {1.4681 + 12.2584\phi - {{20.7322}\phi}^2 + {{15.8855}}\phi}^3 \right)\end{align}

Since ash particles are of irregular shapes, the constants $a$ , $b$ , $c$ and $d$ are based on the shape factor $\phi = {A_S}/{A_P}$ ( ${A_S}$ area of a sphere of equivalent volume, and $\;{A_P}\;$ actual area of the particle). De-Giorgi et al. [Reference De Giorgi, Campilongo, Ficarella, Coltelli, Pfister and Sepe51] by analysing the Etna volcano ash particles with a scanning electron microscopy found that the shape factor is equal to 0.68, in the same range of 0.6–0.8 recommended by Riley et al. [Reference Riley, Rose and Bluth52] based on three samples of volcanic ash. In this work the value of 0.7 is opted for.

The combined gravitational force and buoyancy force per unit of particle mass is as follows:

(3) \begin{align}\overrightarrow {{f_{GB}}} = \left( {1 - \frac{{{\rho _f}}}{{{\rho _p}}}} \right)\vec g\end{align}

The pressure gradient force relates the effect of the local fluid pressure gradient around a particle, which when reduced by particle mass becomes:

(4) \begin{align}\overrightarrow {{f_P}} = - \frac{1}{{{\rho _p}}}\vec \nabla p\end{align}

The slip shear layers induced force is due to non-uniform relative flow velocities. Saffman [Reference Saffman53] established the first expression, which when developed in 3-D and reduced by particle mass becomes:

(5a) \begin{align}\overrightarrow {{f_{SL}}} = \frac{3}{4}\frac{{{\rho _f}}}{{{\rho _p}}}\;{C_{LS}}\left( {\overrightarrow {{V_f}} - \overrightarrow {{V_p}} } \right) \times \overrightarrow {{\omega _f}} \end{align}

With the lift coefficient [Reference Sommerfeld46] given by:

(5b) \begin{align}\;{C_{LS}} = \frac{{4.1126}}{{\sqrt {R{e_s}} }}f\left( {\;R{e_p},\;R{e_s}} \right)\end{align}

The correction function $f\left( {\;R{e_p},\;R{e_s}} \right)$ is proposed by Mei [Reference Mei54] for particle Reynolds number in the range of $0.1 \le R{e_p} \le 100$ and $\beta = 0.5\frac{{\;R{e_s}}}{{\;R{e_p}}}$ .

(5c) \begin{align}\begin{array}{l@{\quad}l}f\left( {\;R{e_p},\;R{e_s}} \right) = \left( {1 - 0.3314\sqrt \beta \;\;} \right){e^{ - R{e_p}/10}} + 0.3314\sqrt \beta & R{e_p} \le 40\\[5pt]f\left( {\;R{e_p},\;R{e_s}} \right) = \;0.0524\sqrt {\beta R{e_p}} & R{e_p} \gt 40\end{array}\end{align}

The shear flow Reynolds number ( $R{e_s}$ ) is given by:

(5d) \begin{align}R{e_s} = \frac{{{\rho _{fd_p^2}}}}{{{\mu _f}}}\|\overrightarrow {{\omega _f}}\| \;\;\;\;\;\;\;\;\;{\rm{where}}\quad \overrightarrow {{\omega _f}} = \nabla \times \overrightarrow {{V_f}} \end{align}

4.2 Turbulence effect

The stochastic random walk implemented within the Lagrangian particle tracking model is based on the eddy interaction model of Gosman and Ioannides [Reference Gosman and Ionnides55]. This latter assumes that a particle and an eddy interact as long as the interaction time $\Delta {t_{in}} = {\rm{min}}\left( {\Delta {t_e},\,\,\Delta {t_r}} \right)$ [Reference Brown and Hutchinson56] is shorter than the eddy lifetime $\Delta {t_e} = 0.37\frac{k}{\epsilon }\;$ and the particle displacement is smaller than the eddy’s length ${l_e} = 0.3{{{k^{3/2}}} \over \epsilon }$ (Shirolkar et al. [Reference Shirolkar, Coimbra and Queiroz49]). The transit time $\Delta {t_r}$ is the duration taken by a particle to cross through an eddy, which has a linearised form [Reference Gosman and Ionnides55] function of ${l_e}$ and relaxation time ${\tau _p}$ as follows:

(6) \begin{align}\Delta {t_r} = - {\tau _p}Ln\left[ {1 - {{{l_e}} \over {{\tau _p}\|\overrightarrow {{V_f}} - \overrightarrow {{V_p}}\| }}} \right]\;\,{\rm{and}}\,\,{\tau _p} = {4 \over 3}{{{\rho _p}} \over {{\rho _f}}}{{{d_p}} \over {{C_{D \| \overrightarrow {{V_f}} - \overrightarrow {{V_p}} \| }}}}\end{align}

When a particle is trapped inside an eddy the turbulent components $u',v',w' = \xi\sqrt {\frac{2}{3}k} $ are added to the mean velocity components of the fluid through a Gaussian random number $\xi $ . The turbulence effect is updated when the particle either passes into a new mesh cell or exceeds the eddy lifetime which is a function of the specific dissipation rate.

Figure 11. Impact conditions.

4.3 Boundary conditions

The true impact on a wall is determined by half the diameter of a particle (Fig. 11) to avoid non-physical trajectories and impacts. Interpolation is used to determine a more precise time step, and the calculation is repeated. At an impact point, the local surface unit vectors $\vec n\;$ and $\vec t\;$ (Fig. 11) are determined and used to calculate the normal and tangential components of impact velocity. Tabakoff et al. [Reference Tabakoff, Grant and Ball57] first studied the post-impact properties of solid particles using high-speed photography to analyse the rebound correlations for annealed 2024 aluminium alloy. Then followed the work of Grant et al. [Reference Grant, Ball and Tabakoff58] carried out in an erosion wind tunnel. By using the transit laser anemometer, Tabakoff et al. [Reference Tabakoff, Hamed and Murugan59] determined the restitution factors for typical turbomachinery materials such as 2024 Aluminium, TiA1-4, AM355 and Rene 41. The particle impact is influenced by its hardness, irregular shape and surface roughness, leading to the rebound characteristics which are dispersed around the averaged values, considering the values of standard deviation [Reference Tabakoff, Hamed and Murugan59]:

(7a) \begin{align}{e_v} = {V_{p2}}/{V_{p1}} = \mathop \sum \limits_{i = 0}^4 {a_i}\beta _i^i \quad {\rm{and }}\quad {\sigma _{\left( {{e_v}} \right)}} = \mathop \sum \limits_{i = 0}^4 {a_{\sigma i}}\beta _i^i\end{align}

(7b) \begin{align}{e_\beta } = {\beta _{p2}}/{\beta _{p1}} = \mathop \sum \limits_{i = 0}^4 {b_i}\beta _i^i\; \quad {\rm{and}}\quad {\sigma _{\left( {{e_\beta }} \right)}} = \mathop \sum \limits_{i = 0}^4 {b_{\sigma i}}\beta _i^i\end{align}

At the symmetry boundary condition, such as the axis of rotation, the velocity components are updated by pure reflection. When it comes to the periodic boundary condition, the particle coordinates and velocity components are transferred to the opposite side using the pitch angle. A FEM-based search algorithm is used to determine the natural coordinates and update for the cell number. At the interface between the stationary domain and rotating domain, the coordinates of a particle are transferred locally to the next frame and the FEM algorithm searches for the next cell number. The particle velocity is updated using the velocity composition.

(8) \begin{align}{\overrightarrow{V_{p_{{}_{abs}}}}}= {\overrightarrow {V_{p_{{}_{rel}}}}} + \overrightarrow {\Omega} \times \overrightarrow {r\;} \end{align}

4.5 Particle seeding

Particle trajectory and erosion calculations were performed considering four samples of volcanic ash collected from various volcanic regions around the globe, and close to the eruption source. The first sample was collected from the Eyjafjallajokull volcano in Iceland (North Europe, 2010), and particle sizes determined through sieving (Fig. 12(a)) to mimic steady-state flight conditions [Reference Wylie, Bucknell, Forsyth, Mcgilvray and Gillespie60]. The second sample (Fig. 12(b) was collected from the Chaiten volcano in Chile (South-America, 2008) [Reference Watt, Pyle, Mather, Martin and Matthews61]. The third sample (Fig. 12(c)) was collected from the Etna volcano in Italy (South Mediterranean, 2002) [Reference De Giorgi, Campilongo, Ficarella, Coltelli, Pfister and Sepe51]. Lastly, the fourth sample (Fig. 12(d)) corresponds to the Kelud volcano in Indonesia (South Asia, 2014) [Reference Sahay, Pati, Singh, Niyogi, Chakarvorty, Prakash and Dwivedi62]. The cumulative frequency curves of particle sizes of the four samples are shown in Fig. 13. The mean size and standard deviation (Table 3) of each ash sample are calculated from the size density derived from the cumulative frequency curve. The density of volcanic ash (Table 3) is influenced by the quartz concentration, as reported by Lau et al. [Reference Lau, Grainger and Taylor63]. The main component of volcanic ash is quartz, ranging from 47 to 73.2% [Reference Lau, Grainger and Taylor63]. The Chaiten sample, from a rhyolitic eruption, has the highest SiO₂ percentage of 73.9%, while the samples from Eyjafjallajokull, Kelud, and Etna have percentages of 57.8, 54.07, and 47.1% respectively, which are mainly consisting of basaltic material.

Figure 12. SEM images of different volcanic ash samples: (a) first sample from Eyjafjallajokull [Reference Wylie, Bucknell, Forsyth, Mcgilvray and Gillespie60] (b) second sample from Chaiten [Reference Watt, Pyle, Mather, Martin and Matthews61] (c) third sample from Etna [Reference De Giorgi, Campilongo, Ficarella, Coltelli, Pfister and Sepe51] (d) fourth sample from Kelud [Reference Sahay, Pati, Singh, Niyogi, Chakarvorty, Prakash and Dwivedi62].

Figure 13. Cumulative frequency of ash particle sizes.

The maximum concentration in the vicinity of an aircraft is at least equal to 4mg/m³ [Reference Davison and Rutke9, Reference Clarkson, Majewicz and Mack11], although the level of confidence in the estimated concentrations is low and uncertainties of an order of magnitude are possible [Reference Witham, Webster, Hort and Thomson64]. In this study, the number of ash particles and seeding positions, radially and tangentially, in the Pitot intake are randomly set according to the concentration of particles, size distribution (Fig. 13), intake geometry and flow conditions. The mass rate of volcanic ash particles is obtained by multiplying the air volume flow rate and particle concentration of successively 1, 2 and 4mg/m³. A subprogram continuously iterates on the particle sizes and release positions until the convergence of the total particle mass rate.

5.0 Erosion model

Early studies on erosion have concluded that this phenomenon is affected by various factors, such as the characteristics of erodent particles, size, concentration, impingement velocity and angle, and the properties of the target material. Finnie [Reference Finnie65] established the basis for the sand erosion, which was later improved by Bitter [Reference Bitter66] to account for the effect of plastic deformation. Grant and Tabakoff [Reference Grant and Tabakoff67] developed successful erosion models for aerospace materials, beginning with aluminium alloys and quartz particles. Their correlation of erosion rate incorporates two mechanisms predominant at low impact angle and high impact angle, and expressed as material removed (milligram) per mass of impacting particles (in grams). Tabakoff et al. [Reference Tabakoff, Kotwal and Hamed68] developed further erosion models for aluminium, stainless steel, and titanium impacted by ash particles, which are expressed as follows:

(9) \begin{align}\epsilon = {k_1}{\left[ {1 + ck.{k_{12}}{\rm{sin}}\left( {90/{\beta _0}} \right){\beta _1}} \right]^2}V_{p1}^2{\cos ^2}{\beta _1}\left[ {1 - R_\theta ^2} \right] + {k_3}{\left( {{V_{p1}}\sin {\beta _1}} \right)^4}\end{align}

The tangential restitution ratio: ${R_\theta } = 1 - 0.0016{V_{p1}}\sin {\beta _1}$ , and the constant $ck = \left\{ {\begin{array}{*{20}{c}}{1\;\;if\;\;{\beta _1} \le 3{\beta _o}}\\{0\;\;\;if\;\;\;{\beta _1} \gt 3{\beta _o}}\end{array}} \right.$

Table 3. Volcanic ash particle size and density

The Pitot intake and spinner are made from aluminium alloy, whereas the fan blades, IGVs and OGVs are made from titanium. For the 2024 aluminium alloy and ash as erodent particles, the angle of maximum erosion $\;\;{\beta _0} = 20\;{\rm{deg}}\;$ and the material constants are: $\;{k_1} = 0.156988 \times {10^{ - 5}}$ , ${k_{12}} = $ 0.3193, ${k_3} = 2 \times {10^{ - 12}}$ . For the titanium alloy 6Al-4V and ash as erodent particles, the angle of maximum erosion $\;\;{\beta _0} = 20\;{\rm{deg}}\;$ and the material constants are: $\;{k_1} = 0.1564951 \times {10^{ - 5}}$ , ${k_{12}} = $ 0.173636, and ${k_3} = 3 \times {10^{ - 12}}$ .

At each impact point, the eroded mass is calculated using the erosion rate evaluated from Equation (9). An area around the concerned node, denoted centred area $\;{A_{node}}$ , is formed by taking a quarter of area of the four faces that share this node. The values of eroded mass per unit of time (second) are summed and divided by the area $\;{A_{node}}$ to calculate the equivalent erosion rate $E$ (mg/s.mm²) according to Equation (10) and distributed among the nodes of impacted faces.

(10) \begin{align}E = \frac{{\mathop \sum \nolimits_1^{ni} {\epsilon _i}{m_{pi}}}}{{{A_{node}}}}\end{align}

6.0 Code structure and validation

The organisation and functionalities of the FORTRAN code PARTRAJ, depicted in Fig. 14, are summarised as follows: COMMON BLOC contains all the arrays and matrices allocations, and the declarations and common variables. INPART reads the control and initial parameters. The geometry, grids, and flow data are imported from the flow solver. PGEODATA reads the geometry and topology, while FLODATA reads the flow data. ARDIST distributes particles at the release plane according to the selected injection mode, size distribution, and concentration. BOUNDARY sets all boundary conditions. INTVEL interpolates for the flow data. EXFORCE calculates all activated external forces. RKUTTA integrates the particle motion equations. JACOBI computes the shape function and the Jacobian matrix. TRANSF converts the global coordinates to the local coordinates. LOCATE searches for particle position and updates for the cell number. REBOUND calculates the rebound properties. STATISTIC generates random numbers and calculates the mean and variance for the given distributions. SPLINT interpolates functions. EROSION calculates erosion rate and cumulates eroded mass using CUMEROS , and calls DISTEROS to calculate erosion rate density and its distribution. DEGGEOM calculates depths of penetration and deteriorated geometry. OUTPUT generates files for the plots of trajectories, erosion and deteriorated geometry.

Figure 14. PARTRAJ flowchart.

This code has been validated based on a high-speed axial fan stage of 11,300rpm, made of aluminium and used for aircraft equipment ventilation, as detailed in the author’s thesis [Reference Ghenaiet69]. Silica sand MIL-E 5007E of particle sizes ranging from 0 to 1,000 microns was used to test and simulate particle trajectories and erosion. Figure 15(a) shows sand particle paths deviating from the streamlines due to inertia and radial spread of particles. Large particles follow ballistic trajectories after hitting walls, while the small ones closely follow the streamlines. In IGVs, the impacts mainly occur at the aft upper corner of the PS, while the remaining particles traverse completely without collisions. Figure 15(b) shows frequent impacts around the blade’s LE caused by sand particles at high velocity and incidence. High frequency of impacts result in intense erosion spreading in a triangular region on the PS, and over the tip of the blade. Particles migrating from the blade’s PS to the SS, through the tip clearance, cause clear erosion of the tip. Figures 15(c) and 15(d) illustrate the computed erosion rate densities on the PS and SS of the blade, which are consistent with the paint stripping test results. Figures 16(a) and 16(b) reveal the measured performance degradation based on the tested eroded blade profiles compared to the predicted performance degradation based on the computed eroded blade profiles. Initially, at the design point defined for the mass flow rate of 0.8kg/s, the fan stage had an efficiency of 0.802 and a pressure rise coefficient of 0.395. After nine hours of sand ingestion, the measured performance droped to the successive values of 0.748 and 0.381, while the predicted performance declined to the values of 0.740 and 0.360, respectively.

Figure 15. Sample of trajectories (a) sample of impacts (b) computed erosion rate density on the blade PS (c) and SS (d) tested erosion patterns showing the PS (e) and SS (f).

Figure 16. Compared performance degradation: (a) Pressure rise coefficient (b) Efficiency.

Figure 17. Samples of trajectories of particle size 1μm coloured by velocity: (a) Meridional view (b) Top view of fan blade, OGVs and IGVs.

7.0 Particle trajectories results

The results of the particle trajectory calculations include information about the particles coordinates, velocities, collision conditions, and impact frequency. The behaviour of ash particles trajectories can be seen in different components, including the Pitot intake, spinner, fan blade, splitter, IGVs, OGVs and outlet ducts. It is evident that most of released particles do not reach beyond 80% of the fan blade span due to the shape of the Pitot intake lip and inner contour. Many ash particles are observed to deflect downward as they collide with the intake lip and strike the spinner, rotor’s hub, and lower sections of the fan blade. The paths of small ash particles closely align with the streamlines, as revealed by the trajectories of particles circulating around the intake lip. Small particles released centrally accelerate along the intake duct and impact the spinner while others move away alongside the spinner, and by entering the rotor they become strongly affected by the centrifugal force. On the other hand, large ash particles of higher inertia than the drag, acquire an outward radial velocity after striking the spinner, while others after collisions bounce off and curve upwards to collide with the upper sections of the fan blade and shroud. Generally, the deviation of particle paths from the streamlines increases with the particle size. Heavy ash particles are more likely to strike the blade surfaces and outer walls at high velocities, bouncing, and bending trajectories.

Figure 17(a) displays the paths of 1-micron ash particles, with their velocities indicated by colours. These particles generally follow the flow streamlines and are influenced by the turbulence and secondary flows. These tiny particles circulate and collide repeatedly around the Pitot intake lip and travel alongside the inner contour. Evidence from Fig. 17(a) reveals that these particles adhere to the streamlines and closely follow the flow details, to reach velocities up to 466m/s. As observed from the top view (Fig. 17(b)), small ash particles travelling around the blade are entrained by the leakage flow over the blade tip, thereby causing multiple impacts and erosion. Figure 18(a) depicts the trajectories of 10 microns ash particles. Many of these particles are seen to circulate around the intake lip a short distance before inducing multiple collisions and being diverted slightly. After recovering kinetic energy from the flow, they move directly towards the fan blade. The same particles when released near the center of intake, by reaching the front of the fan they circulate away from the spinner and only induce a few collisions. Owing to the high peripheral blade speed and flow velocity combined with the blade twist, many high-velocity particles impact the blade’s LE, spanning from the root upward, as depicted in Fig. 18(b). Added to that, there are multiple impacts spreading over the PS of the fan blade. Figure 18(b) shows particles reaching the top of the blade at high velocities up to 457m/s, which when they hit the blade PS at low-impact angles become practically aligned with the streamlines. Other particles at the inner sections of the fan blade, are seen to impact and bounce off at high rebound angles but at relatively lower velocities. Due to the shape of the intake lip and inner contour pushing ash particles away, the majority of them reach the fan blade along approximately 80% of the span, while the rest exhibits scattered impacts. Upon leaving the rotor blade, the particles change their directions and pass to the IGVs and OGVs. The tiny ash particles such as 1 micron are observed to travel through the stationary vanes and follow closely the streamlines while exchanging kinetic energy from the fluid. The large particles of 10 microns seem to deviate notably (Fig. 18(b)) from the streamlines and impact more frequently the PS of the stationary vanes. As also observed, the small ash particles issued from the lower sections of the fan blade are redirected by the secondary flow to reach the LE of IGVs at high incidence and travel farther while inducing several impacts.

Figure 18. Samples of trajectories of particle size 10μm coloured by velocity: (a) Meridional view (b) Top view of fan blade, OGVs and IGVs.

Figure 19. Samples of trajectories of particle size 100μm coloured by velocity: (a) Meridional view (b) Top view of the fan blade, OGVs and IGVs.

Figure 19 illustrates the trajectories of large ash particles (i.e. 100 microns), which are shown to collide repeatedly around the intake lip and diverte downward in sort of ballistic trajectories and impact the lower sections of the fan blade. As shown in Fig. 19(a), the particles released centrally travel along the intake, strike the spinner, bounce back, and then turn to follow paths in a parabolic manner. Moreover, owing to the influence of flow structure and imparted centrifugation, large particles are forced to move upwards and reach the upper parts of the fan blade. The high flow velocities in this transonic fan, leading to high drag force, combined with the significant centrifugal force, cause the large particles to impact the LE of the fan blade from both sides and rebound to reach the opposite PS at the aft of the surface. As observed in Fig. 19(b), the particle velocity increases significantly while moving across the fan blade, reaching up to 434m/s. These large particles are seen to impact the PS of the blade near the optimum impact angle and rebound off from the surface, and after multiple collisions they exit the fan blade as coalesced particles. By crossing the rotor interface, their directions change significantly. Through IGVs or OGVs, they seem to travel almost straight, as depicted in Fig. 19(b), and impact the stationary components at high velocities and angles, primarily at the LE from the SS and more frequently on the PS of the vanes. It should be noted that the majority of ash particles, after undergoing significant centrifugation by the fan blade, experience higher frequency of impacts on the OGVs' surfaces and the shroud compared to the small particles. Furthermore, owing to intense centrifugation imparted by the fan blade, a large amount of ash particles bypass outward, preventing them from traversing the IGVs and entering the engine core.

To demonstrate the behaviours of ash particles of various sizes randomly released, trajectory simulations were performed for the ash samples Eyjafjallajokull (0.4–80 microns) and Kelud (0.5–300 microns), which are characterized by relatively small and large particle sizes, respectively. These particles, randomly injected upstream of the Pitot intake combine all the discussed aspects associated with small, medium and large ash particles. Their paths along the Pitot intake and through the fan blade, OGVs and IGVs, and outlet ducts are depicted in Fig. 20, and coloured by the particle diameter. As observed in Fig. 20(a), the small particles belonging to the first sample tend to bypass the spinner without collisions and move away from the rotor’s hub, and then deviate upward to reach the upper sections of the fan blade. On the other hand, the large particles released centrally are shown to strike the spinner, rebound, and deflect upward to reach the upper sections of the fan blade. Even particles moving near the lower sections of fan blade are redirected upwards and projected against the shroud, due to the high centrifugation. The trajectories of particles belonging to the fourth ash sample, of medium and large sizes (Fig. 20(b)), appear more diverted while crossing the intake. After multiple collisions with the intake lip, these particles follow ballistic paths towards the lower sections of the fan blade. When ash particles collide with the spinner, they rebound backward and then curve upward in parabolic trajectories, striking the upper regions of the fan blade and shroud. For the same reasons, larger particles (in red colour) repeatedly collide with the shroud and are more likely to infiltrate through the tip gap. Upon leaving the fan blade, ash particles enter the OGVs/IGVs at high incidence angles, while others deflect upward and move alongside the shroud. Heavy ash particles are predominantly diverted upward by the fan rotor, pass through the OGVs and are finally expelled through the bypass nozzle. On the other hand, the light particles are less diverted and may reach the engine core. Consequently, this large fan operates as a separator for the large-size particles.

Figure 20. Trajectories coloured by particle diameter: (a) Sample 1 Eyjafjallajokull (b) Sample 4 Kelud.

8.0 Erosion results and discussions

The local erosion rates and patterns are obtained in cruise flight conditions, considering ash particle concentrations 1–4 ${\rm{mg}}/{{\rm{m}}^3}$ and the four samples collected from the major volcanic regions worldwide. Due to high number of impacts on each face of a mesh element, the local erosion rates $\left( {{\rm{mg}}/{\rm{g}}} \right)$ determined from Equation (9) were used to evaluate the erosion rate density ( ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2})$ based on Equation (10). To facilitate the discussion, erosion patterns of the four ash samples are presented only for the concentration of 2mg/m³. Contours of erosion rate density are displayed with the same scales to compare between patterns of different samples. As noticed, higher erosion rates are obtained with the Kelud’s ash sample, followed by that of the Etna volcano.

Figure 21. Erosion rate density on the fan blade caused by volcanic ash, at concentration of 2mg/m³: (a) Sample 1 Eyjafjallajokull (b) Sample 2 Chaiten (c) Sample 3 Etna (d) Sample 4 Kelud.

Figure 22. Erosion rate density on the fan hub and shroud caused by volcanic ash, at concentration of 2mg/m³: (a) Sample 1 Eyjafjallajokull (b) Sample 2 Chaiten (c) Sample 3 Etna (d) Sample 4 Kelud.

Figure 23. Erosion rate density in OGVs and IGVs caused by volcanic ash, at concentration of 2mg/m³: (a) Sample 1 Eyjafjallajokull (b) Sample 2 Chaiten (c) Sample 3 Etna (d) Sample 4 Kelud.

Figure 24. Erosion rate density in Pitot intake, spinner and fan blades caused by volcanic ash, at concentration of 2mg/m³: (a) Sample 1 Eyjafjallajokull (b) Sample 2 Chaiten (c) Sample 3 Etna (d) Sample 4 Kelud.

The patterns of erosion, depicted in Figs. 21–24, are globally similar, except the slight variations in erosion rates, spread, and some specific details between the four ash samples. The first and second samples from the Eyjafjallajokull and Chaiten volcanoes caused dense erosion spreading more on the fan blade’s PS, but with lesser erosion intensity on the SS. On the other hand, the third and fourth samples from the Etna and Kelud volcanoes resulted in scattered erosion on both sides of the fan blade due to dispersed high-velocity impacts, leading to higher local erosion rates compared to the first and second ash samples. As depicted in Fig. 21(a)–(d), there are noticeable areas of erosion spreading on the PS of the fan blade, typically from the hub to 80% of the span, where the frequency of impacts and velocities are very high. The remaining of the top sections of the fan blade from the PS are characterized by fewer impacts and subsequently lower erosion rates. The maximum of local erosion rates on the PS of the fan blade are equal to 1.352×10⁻⁵, 1.974×10⁻⁵, 1.106×10⁻⁴ and 1.791×10⁻⁴ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ for the first, second, third and fourth ash samples, respectively. It is evident, that the fourth and third ash samples caused the highest erosion densities, due to higher proportion of large size particles. The LE and front of the fan blade were directly impacted by high-velocity particles, resulting in a band of extreme erosion at the LE that extends to both sides. Figures. 21(a)–(d) exhibit pronounced erosion extending upward from the hub and covering a large portion of the PS until the TE, where the third and fourth samples led to extensive erosion rates. The scattered areas of dense erosion observed at the top sections of the fan blade are caused by the centrifuged heavy particles added to the high-velocity particles reflected upward by the spinner. Moreover, the SS of the fan blade reveals impacts spreading from the LE up to two-thirds of the surface and resulting in scattered erosion. The front of the blade’s SS presents concentrated erosion, but of lower rates compared to the PS. As observed, the spread of erosion and rates are the highest with the third and fourth ash samples, characterised by a broader range of particle sizes, higher density and quartz content. According to Figs. 22(a)–(b), the rotor’s hub is severely eroded from its entrance alongside the PS whereas the opposite side has scattered erosion due to irregular impacts and particles deflected away around the blade root. Furthermore, the rotor’s shroud presents scattered erosion due to irregular impacts caused by the centrifuged and reflected heavy particles from the spinner and blade root sections, and subsequently, the cumulated erosion appears more scattered compared to the hub.

Erosion of the IGVs, as illustrated in Figs. 23(a)–(d), extends across both the PS and SS of the vanes, and are characterized by erosion rate densities reaching up to 5.108×10⁻⁵, 9.601×10⁻⁵, 1.426×10⁻⁴ and 1.351×10⁻⁴ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , between the first, second, third and fourth ash samples, respectively. However, the third and fourth ash samples caused higher erosion rates at the aft of the vane’s PS until the TE, while the SS is eroded across the entire surface until the TE, but at lower rates. The IGVs’ hub is eroded slightly since ash particles are diverted upwards due to preceding centrifugation while the shroud is more eroded near the splitter. The front junction between IGVs and shroud appears the region more prone to erosion wear. On the other side, the OGVs display patterns of erosion spreading across the entirety of the PS, with the highest rates at the rear of the PS, while the lowest rates are observed near the root. Moreover, the main eroded areas of the SS are concentrated on the LE and the front part of OGVs, while the upper sections seem unaffected by erosion and the same for the aft of this surface. The OGVs’ shroud displays erosion patterns, particularly, at the front part, generated by previously centrifuged particles and those following ballistic trajectories. The hub, on the other hand, is eroded predominantly from the splitter entrance, and onward the erosion becomes weak.

Figures 24(a–d) display the front and back views of erosion patterns formed in the Pitot intake, spinner, fan blades, hub and shroud. The highest erosion rates in the Pitot intake is observed on the intake’s lip, unlike the inner contour of much lower erosion. The Pitot intake experiences erosion rate density of up to 2.811×10⁻⁸, 2.805×10⁻⁸, 7.469×10⁻⁸, and 7.865×10⁻⁸ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , for the first, second, third and fourth ash samples, respectively, but of lower rates compared to the spinner and the fan blades. The spinner displays erosion spreading more near the fan hub due to the high frequency of impacts, and subsequently, the erosion rate densities are high and equal to 1.248×10⁻⁶, 1.191×10⁻⁶, 3.291×10⁻⁶ and 9.005×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ between the first, second, third, and fourth ash samples, respectively. As noticed, the fourth ash sample caused more erosion of the spinner. The back view shows that the fan blades suffer more erosion mainly from the blades’ PS, with extreme erosion seen at the rear part and the TE. Accordingly, the fourth ash sample induced the highest erosion rate density of up to 1.791×10⁻⁴ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^{2{\rm{\;}}}}$ . On the other side, the front view of the components displays erosion extending across the intake lip. Also, reveals scattered erosion on the blades’ SS, with the highest rates observed near the LE. One may conclude that the most of erosion scatters, of higher rates, are induced by the fourth ash sample.

The average erosion rate density (AERD), in milligrams per second and square millimetre (mg/s.mm²), is another metric for assessing the overall erosion and characterising how erosion develops in each component. This parameter is determined using the local erosion rates, mass of impacting particles, and component’s area. To evaluate the AERD, the erosion rate density $E$ based on Equation (10) is integrated over the whole surface of the concerned component and then divided by its area. In addition, the hourly eroded mass (HEM), expressed in milligrams per hour (mg/hr), provides a comprehensive assessment of the material removal from a component. The HEM is calculated by multiplying the value of AERD by the component’s area, for duration of 3,600 seconds.

Table 4. Erosion parameters caused by four volcanic ash samples

Table 4 summarises the obtained results in terms of ash particle concentration, particle mass rate, AERD, and HEM. At the highest ash particle concentration, the fan blades’ AERD, reached the maximum values of 0.161×10⁻⁶, 0.162×10⁻⁶, 0.266×10⁻⁶ and 0.276×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ between the first, second, third and fourth ash samples, respectively. It is evident that the fourth and third samples led to the highest values of AERD. Moreover, the HEM (mg/hr), for the totality of 38 blades, is equal to 5.871×10³, 5.910×10³, 9.705×10³, and 10.067×10³mg/hr, where the highest values occur for the third and fourth ash samples. The fan’s hub exhibits AERD values of 0.175×10⁻⁶, 0.178×10⁻⁶, 0.335×10⁻⁶, and 0.352×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , for the first, second, third and fourth ash samples, respectively. Accordingly the HEM values are 0.514×10³, 0.523×10³, 0.983×10³ and 1.033×10³mg/hr, where the highest values are caused by the third and fourth samples. The shroud reveals AERD values of 0.064×10⁻⁶, 0.583×10⁻⁶, 0.525×10⁻⁶, and 1.430×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , between the first, second, third and fourth ash samples, respectively, leading to HEM of 0.4952×10³, 4.487×10³, 4.041×10³ and 11.013×10³mg/hr, with the highest values corresponding to the fourth ash sample. Overall, the fourth ash sample induced the most severe erosion of the fan blades, followed by the third sample. The rotor’s shroud, despite its lower erosion rates, has a higher HEM than the hub because of larger surface.

Figure 25. HEM (g/hr) of the complete fan rotor: (a) blades (b) hub (c) shroud.

Figure 26. HEM (g/hr) of complete: (a) IGVs (b) OGVs.

Figure 27. HEM (g/hr) of complete: (a) spinner (b) Pitot intake.

The entire IGVs are characterized by AERD values equal to 0.404×10⁻⁶, 0.386×10⁻⁶, 0.506×10⁻⁶, and 0.514×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , respectively, and the subsequent HEM are equal to 2.337×10³, 2.233×10³, 2.923×10³, and 2.968×10³mg/hr, respectively. Moreover, the IGVs’ shroud, presents HEM values of 0.522×10³, 0.551×10³, 0.983×10³, and 1.008×10³mg/hr, while the erosion of the hub is negligible. The fourth and third ash samples exhibit higher values compared with the two other samples. As for the OGVs, the values of AERD are respectively equal to 0.095×10⁻⁶, 0.098×10⁻⁶, 0.136×10⁻⁶, and 0.141×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , and their respective HEM values are 2.277×10³, 2.341×10³, 3.245×10³, and 3.375×10³mg/hr. Furthermore, the OGVs’ shroud depicts values of HEM equal to 0.164×10³, 0.165×10³, 0.565×10³, and 0.817×10³mg/hr, respectively, while the hub has insignificant erosion.

The spinner is characterized by values of AERD equal to 0.341×10⁻⁶, 0.375×10⁻⁶, 0.973×10⁻⁶, and 1.099×10⁻⁶ ${\rm{mg}}/{\rm{s}}.{\rm{m}}{{\rm{m}}^2}$ , and the subsequent HEM values are 0.736×10³, 0.811×10³, 2.10×10³, and 2.373×10³mg/hr between first, second, third and fourth ash samples, respectively. Despite the low erosion rate density in the Pitot intake, its large geometry resulted in high values of HEM equal to 0.115×10³, 0.113×10³, 0.161×10³, and 0.163×10³mg/hr between the four ash samples, respectively, where the fourth sample produced more erosion.

Finally, Figs. 25–27 compare between the values of HEM of the Pitot intake, spinner, entire fan blades, and entire IGVs and OGVs, depicting essentially a marked increase with ash concentration. The highest values of HEM correspond to the third and fourth ash samples. As outlined in Table 5, in these front components of an HBTFE, the fourth volcanic ash sample (Kelud), followed by the third sample (Etna), caused the most severe erosion due to their relatively broader size range (0.5–300 microns) and (0.5–140 microns), combined with higher density and quartz content.

Table 5. Compared HEM at high ash concentrations

The continual removal of material from the surfaces of this transonic fan stage, particularly the rotor blades, causes a rapid decline in aerodynamic performance and lifespan. By simulating various operating situations, ash particle types, and concentration, a correlation can be formed between particle type and concentration, operating conditions, blade geometry deterioration, performance degradation, and lifespan.