A Case Study: Sediment Erosion in Francis Turbines Operated at the San Francisco Hydropower Plant in Ecuador

Abstract

The operation of various types of turbomachines is importantly affected by sediment erosion. Francis turbines used for power generation typically suffer said effects due to the fact that they are used in sediment-laden rivers and are usually operated disregarding the long-term effect of the erosion on turbine performance. This investigation seeks to study the erosion rate for the main components of the turbines located at San Francisco hydropower plant in Pastaza, Ecuador. A sediment characterization study was performed in order to determine the properties of the particles present in Pastaza River and accurately predict their effect on the turbine flow passages. A numerical approach combining liquid–solid two-phase flow simulation and an erosion model was used to analyze the erosion rates at different operating conditions and determine wear patterns in the components. As expected, the results indicated that an increase in the erosion rate was obtained for higher intake flows. However, a dramatic increase in the erosion rate was observed when the turbine was operated at near-full-load conditions, specifically when guide vane opening exceeded a 90% aperture.

1. Introduction

Hydroelectric power is a renewable energy source and a significant component of worldwide electricity production. Around 17% of the total consumed electricity is produced through hydraulic energy sources [1,2], and almost 65% of the total electricity produced in Latin America is generated by hydroelectric power plants (around 709 TWh/y) [3,4]. However, most of the total technical hydraulic potential (2859 TWh/y) of the region is not harnessed by its installed hydropower capacity. Several large-scale hydropower projects are being currently studied and developed in the Andean region in the hopes of increasing the installed capacity and harnessing a larger portion of the available hydraulic potential. One of the most crucial factors that needs to be taken into consideration during the development of the aforementioned projects is the fact that hard particles are present in almost all rivers of the Andean and Himalayan region, causing considerable the erosion, mechanical wear, and failure of turbine components [3,5].

Sediments flowing through the river deposit in the dam’s reservoir, reduce the reservoir’s capacity, and increase the erosion wear of critical turbine components, such as: the spiral casing, guide vanes, runner, and draft tube. This phenomenon reduces the lifespan of the turbine and decreases its efficiency, which increases the cost of maintenance over time, leading to economic losses [6,7]. Erosion wear depends on several factors, including particle concentration, velocity, composition, size, and shape. Other variables, namely turbine materials and operating conditions, also have an effect on the erosion rate. Therefore, erosion reduction strategies can only work effectively after in-depth analyses making a holistic assessment of all the variables involved [8,9].

Extensive research has been conducted on erosion in Francis turbines. In 2013, Singh and Banerjee performed an analysis on the erosion of the runner blades, guide vanes, and labyrinth seals of the Maneri Bhali Stage-II hydroelectric power plant in India. Data collection of sediments at relevant locations and measurements of turbine efficiency were performed during three years to determine the effect of silt erosion on the efficiency of turbines [10]. In 2016, Koirala used a computational analysis coupled with field observations to determine the erosion patterns on the guide vanes of Kaligandaki hydroelectric power plant in Nepal and proposed erosion protection methods [11]. A year later, Masoodi and Harmain presented a detailed comparison of two sediment-laden rivers and their effect on the runner blades of Himalayan hydroelectric power plants in India. A new erosion model was proposed in this study [12]. Most recently, in 2020, Qian et al. executed a study on the erosion wear of the runner blades of a Francis turbine in Jhimruk Hydroelectric Center in Nepal using numerical simulations and comparing the results with the damage of the runners. He proposed changing the opening of the guide vanes to improve turbine efficiency and reduce the erosion rate [13]. Moreover, Noon and Kim discussed and analyzed the latest experimental and numerical techniques to determine sediment and cavitation erosion on different turbine components using baseline data from the Tarbela Dam hydroelectric project in Pakistan [14]. However, all the aforementioned studies were performed in Asia, and no research on the topic has been performed on South America, where similar erosion issues are found.

This study focuses on the analysis of sediment erosion in the Francis turbines of San Francisco hydroelectric power plant in Ecuador. The turbines of this power station suffer erosion wear damage, and to date, no effective strategies have been proposed to reduce the damage. A sediment characterization of the Pastaza River was conducted for this study in order to perform a numerical analysis of the turbines with sediment properties set as close as possible to the real conditions. Finally, a study on the erosion rate and pattern in different components of the turbine was carried out to better understand this phenomenon.

2. Case of Study: San Francisco Hydropower Plant

San Francisco hydro-power plant (SFH) is one of the largest energy generation centers in Ecuador, producing around 1140 GWh of electricity per annum, which represents 12% of the energy demand in the country. SFH is a 230 MW run-of-river hydropower plant located along Pastaza River, which consists of two vertical Francis turbines, each one running at 327.27 rpm under a net head of 213.4 m and a flow rate of 58 m³/s. Since the plant began operations in 2007, it has suffered erosion problems, especially at the guide vanes and the outlet band of the runner. Figure 1 presents a georeference of the plant.

2.1. Sediment Characterization at Pastaza River

The Pastaza River originates in the Andes mountains, where irregular geography and the presence of soft sediments due to high volcanic activity contribute to the high sediment content of South American rivers [15]. These sediments pass through the turbines of hydraulic power plants, resulting in sediment erosion on exposed turbine components. To tackle the problems derived from erosion wear, an analysis of the particulate matter that flows through each power plant becomes necessary for a proper assessment of the erosion in a particular station. On this basis, the samples for the analysis were collected from the following zones of the hydropower plant:

• The reservoir after the desilting chamber;
• The outlet of the discharge gate;
• The outlet of the draft tube.

2.1.1. Sediment Analysis

In order to characterize the particles of the Pastaza River at San Francisco hydropower plant, the collected sediment samples were analyzed at the Soil Mechanics and Materials Testing Laboratory, Escuela Politecnica Nacional, performing a sieve and composition analysis. The results of this analysis were considered during the study.

2.1.2. Particle Size and Distribution

A similar approach as the one followed by Koirala et al. [11,16] was used for this study. The sieve analysis was carried out under the ASTM D422-63 (2007) standard. Five sieve measurements of 4.75 mm, 2.00 mm, 0.85 mm, 0.425 mm, and 0.075 mm were used on five 120 g sediment samples.

Additionally, since particle roundness (R) and sphericity (S) affect most macroscale mechanical properties of the particle such as strength, compressibility, and shear wave velocity, it is necessary to estimate these parameters to increase the fidelity of the simulation [17]. Roundness is described as the ratio between the average radius of curvature of the particle corners and the radius of the maximum inscribed circle, while sphericity is defined as the ratio of the particle width to particle length [18]. The sphericity (spherical shape factor) and roundness of the particles were estimated using the Krumbein–Sloss chart [19].

2.1.3. Mineral Composition Analysis

The mineral composition analysis for the study was performed through a particle count method using a D8 ADVANCE X-ray diffractometer and Diffrac plus software.

2.1.4. Sediment Concentration

Pastaza River is the third-largest river in Ecuador with an average annual flow 144.4 m³/s and a precipitation of 3255 mm. Sediment concentration was determined based on data from San Francisco hydropower plant. The samples were collected daily 1 km upstream of the reservoir, in the reservoir itself, in the desilting chamber, and at the discharge during a month, as shown in Figure 2. Since heavier particles tend to deposit on the reservoir bed and mostly only suspended particles are drawn by the intake, water samples were taken at a constant depth of about 2 m from the surface. The average concentration was estimated considering daily values from all sampling points throughout the sampling period.

3. Numerical Analysis of the Erosion in Francis Turbines

3.1. Governing Equations

3.1.1. Liquid Phase Mathematical Model

Fluids were calculated using a Eulerian approach. The general form of the equations involved in the calculations is presented. The mass continuity equation has the following form:

$$\frac{\partial \rho }{\partial t}+\frac{\partial \left(\rho u_i\right)}{\partial x_i}=0$$

(1)

where:

• $u=$ fluid phase velocity;
• $\rho =$ fluid phase density.

The momentum conservation equation is shown in Equation (2).

$$\frac{\rho \partial \left(u_i\right)}{\partial t}+\frac{\rho \partial \left(u_i{u}_j\right)}{\partial x_j}=-\frac{\partial p}{\partial x_i}+\frac{\partial }{\partial x_j}\left[\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)\right]+f_i$$

(2)

where:

• $p=$ pressure;
• $\mu =$ dynamic viscosity;
• $f_i=$ external forces.

3.1.2. Turbulence Model

A realizable k-$ϵ$ turbulence model was selected for its ability to correctly capture the turbulent nature of the flow in Francis turbines [20,21,22]. This model was selected for its robustness and its improved boundary-layer-solving capacity under strong pressure gradients and flow separation compared to the standard k-$ϵ$ model [23]. In addition, the k-$ϵ$ turbulence model has a low computational expense when compared to k-$\omega$ SST. The transport equations for the realizable k-$ϵ$ take the following form:

$$\frac{\partial \left(\rho k\right)}{\partial t}+\frac{\partial \left(\rho k{u}_j\right)}{\partial x_j}=\frac{\partial }{\partial x_j}\left[\left(\mu +\frac{{\mu }_t}{{\sigma }_k}\right)\frac{\partial k}{\partial x_j}\right]+G_k+G_b+\rho ϵ-Y_M+S_k$$

(3)

$$\frac{\partial \left(\rho ϵ\right)}{\partial t}+\frac{\partial \left(\rho ϵu_j\right)}{\partial x_j}=\frac{\partial }{\partial x_j}\left[\left(\mu +\frac{{\mu }_t}{{\sigma }_{ϵ}}\right)\frac{\partial ϵ}{\partial x_j}\right]+\rho C_1{S}_{ϵ}-\rho C_2\frac{{ϵ}^2}{k+\sqrt{vϵ}}+C_{1ϵ}\frac{ϵ}{k}{C}_{3ϵ}{G}_b+S_{ϵ}$$

(4)

$$C_1=max\left[0.43,\frac{\sqrt{2S_{i,j}{S}_{i,j}\frac{k}{ϵ}}}{\sqrt{2S_{i,j}{S}_{i,j}\frac{k}{ϵ}}+5}\right]$$

(5)

where:

• $G_k=$ turbulent kinetic energy due to velocity gradients;
• $G_b=$ turbulent kinetic energy due to buoyancy;
• $Y_M=$ contribution of compressible fluctuations to the dissipation rate;
• $C_{1ϵ},C_2=$ constants;
• ${\sigma }_k,{\sigma }_{ϵ}=$ Prandtl numbers;
• ${\mu }_t=$ Eddy dynamic viscosity;
• $S_k,S_{ϵ}=$ user-defined terms.

The turbulence viscosity ${\mu }_t$ is computed by:

$${\mu }_t=\rho C_{\mu }\frac{k^2}{ϵ}$$

(6)

Model variable $C_{\mu }$ is defined by:

$$C_{\mu }=\frac{1}{4.04+\sqrt{6}cos\varphi \frac{k{U}^{*}}{ϵ}}$$

(7)

$$\varphi =\frac{1}{3}co{s}^{-1}\left(\sqrt{6}W\right)$$

(8)

$$W=\frac{S_{ij}{S}_{jk}{S}_{ki}}{\widetilde{S^3}}$$

(9)

$$\tilde{S}=\sqrt{S_{ij}{S}_{ij}}$$

(10)

$$S_{ij}=\frac{1}{2}\left(\frac{\partial u_j}{\partial x_i}+\frac{\partial u_i}{x_j}\right)$$

(11)

$$U^{*}=\sqrt{S_{ij}{S}_{ij}+\widetilde{{\mathsf{\Omega }}_{ij}}\widetilde{{\mathsf{\Omega }}_{ij}}}$$

(12)

$$\widetilde{{\mathsf{\Omega }}_{ij}}={\mathsf{\Omega }}_{ij}-2{ϵ}_{ijk}{\omega }_k$$

(13)

$${\mathsf{\Omega }}_{ij}=\overline{{\mathsf{\Omega }}_{ij}}-{ϵ}_{ij}{\omega }_k$$

(14)

where $\overline{{\mathsf{\Omega }}_{ij}}$ is the tensor for the mean rate of rotation in a reference frame rotating at an angular velocity ${\omega }_k$.

3.1.3. Solid Phase Mathematical Model

Solid particles were simulated using a Lagrangian approach and were treated as if their volume fraction were low compared to that of the continuous phase. Equation (15) was derived from the force balance on the Lagrangian reference frame.

$$\frac{\partial {\left(v_p\right)}_i}{\partial t}=\frac{18\mu }{{\rho }_p{d}_p^2}\frac{C_{D}Re}{24}+\left(\frac{\rho }{{\rho }_p}\right){\left(v_p\right)}_i\frac{\partial u_i}{\partial x_i}+\left(1+\frac{\rho }{{\rho }_p}\right)g-\frac{1}{2}\frac{\rho }{{\rho }_p}\frac{\partial }{\partial t}(u_i-v_{pi})+F_z$$

(15)

where:

• $v_p=$ particle velocity;
• ${\rho }_p=$ particle density;
• $C_D=$ drag coefficient;
• $Re=$ Reynolds number;
• $F_z=$ other interaction forces per unit mass.

The term $F_z$ activates force terms in situations in which multiple reference frames are used and frame or mesh rotation is activated. The drag coefficient was estimated using the equations for the nonspherical drag law shown in Equation (16).

$$C_D=\frac{24}{R{e}_{sph}}(1+b_{1}R{e}_{sph}^{b_2})+\frac{b_{3}R{e}_{sph}}{b_4+R{e}_{sph}}$$

(16)

$$b_1=exp(2.3288+6.4581\varphi +2.4486{\varphi }^2)$$

(17)

$$b_2=0.0964+0.5565\varphi$$

(18)

$$b_3=exp(4.905-13.8944\varphi +18.4222{\varphi }^2-10.2599{\varphi }^3)$$

(19)

$$b_4=exp(1.4681+12.2584\varphi -20.7322{\varphi }^2+15.8855{\varphi }^3)$$

(20)

where:

• $R{e}_{sph}=$ Reynolds number of an equivalent sphere;
• $\varphi =$ spherical shape factor.

Particle dispersion caused by turbulent flows can be estimated using the stochastic tracking model [24,25]. This method, also known as the discrete random walk model, takes into account the effect of turbulent velocity fluctuations on the trajectories of particles. Instantaneous fluid velocity, as shown in Equation (21), was used to integrate the particle trajectory equations along their path to predict the turbulent dispersion of particles. Random velocity fluctuations $u^{\prime }$ were determined through Equation (22), where $\zeta$ is a normally distributed random number. Particle diffusivity was estimated using Equation (23), where the integral time scale as defined in Equation (24) describes the time the particle remains in turbulent motion along a path $ds$.

$$u=\overline{u}+u^{\prime }\left(t\right)$$

(21)

$$u^{\prime }=\zeta \sqrt{2k/3}$$

(22)

$$D\left(t\right)=\overline{u_i^{\prime }{u}_j^{\prime }}T$$

(23)

$$T={\int }_o^{\infty }\frac{v_p^{\prime }\left(t\right)v_p^{\prime }(t+s)}{\overline{{v^{\prime }}_p^2}}ds$$

(24)

3.1.4. Erosion Model

The erosion model developed by Oka, Okamura, and Yoshida [26,27] was used to determine the erosion for the present case. This model is one of the most frequently used to determine the erosion in CFD analyses where solid particles suspended in a liquid medium are present [28,29,30]. The equation developed by Oka is the following:

$$E\left(\alpha \right)=g\left(\alpha \right)E_{90}$$

(25)

where:

• $E_{\left(\alpha \right)}=$ erosion damage in $m{m}^{3}k{g}^{-1}$;
• $g\left(\alpha \right)=$ impact angle dependence of the normalized erosion;
• $E_{90}=$ erosion damage at a normal angle.

and:

$$g\left(\alpha \right)={\left(sin\alpha \right)}^{n1}{(1+Hv(1-sin\alpha ))}^{n2}$$

(26)

$$E_{90}=K{\left(aHv\right)}^{k_{1}b}{\left(\frac{v}{v^{\prime }}\right)}^{k_2}{\left(\frac{D}{D^{\prime }}\right)}^{k_3}$$

(27)

$$k_2=2.3{\left(Hv\right)}^{0.038}$$

(28)

$$n1,n2=2.3{\left(Hv\right)}^{0.038}$$

(29)

where Hv is the material Vickers hardness in GPa, k${}_2$ is a velocity exponent, k${}_3$ is a diameter exponent, and constants n1 and n2 are model exponents used to calculate the impact angle influence on the erosion rate. D${}^{\prime }$ and v${}^{\prime }$ are the reference diameter and velocity, respectively. The calibrated values of these parameters for the present case, in which sand particles and stainless steel were considered, are the ones presented in

Table 1. Oka model parameters.

Parameter	Units	Value
$k_1$	-	−0.12
$k_2$	-	2.36
$k_3$	-	0.19
$n1$	-	0.78
$n2$	-	1.27
a	-	0.0221
b	-	0.45
$E_{90}$	mm${}^3$ kg${}^{-1}$	3.53

.

3.2. Geometry and Conditions

Details of the general specifications of the turbine are presented in

Table 2. Turbine specifications.

Parameter	Value
Runner inlet diameter (mm)	1530.8
Number of runner blades $Z_b$	13
Height of the guide vane (mm)	540.4
Number of guide vanes	20
Number of stay vanes	20

. In order to reduce the computational effort while using a more precise mesh for the numerical analysis, only one period of the turbine was simulated, which is comprises one stay vane, one guide vane, one runner blade, and an outlet domain representing the draft tube.

The computational domain is presented in Figure 3. The geometry was obtained performing a 3D scanning of the turbine. The obtained profiles were reconstructed using ANSYS BladeGen to obtain a smoother geometry optimized for mesh construction using TurboGrid.

3.3. Operating Conditions

Table 3. Operating conditions.

Parameter	Equation	Unit	Case 1	Case 2	Case 3	Case 4
Guide vane opening	-	%	55.17	77.91	89.92	93.46
Volumetric flow rate	-	m${}^3$ s${}^{-1}$	33.2	50.8	59.7	62.4
Specific speed, $n_{QE}$	$n{Q}^{0.5}{H}^{-0.75}$	-	0.56	0.69	0.75	0.77
Discharge coefficient, $Q_{nD}$	$Q{n}^{-1}{D}^{-3}$	-	3.31	5.06	5.95	6.22
Energy coefficient, $E_{nD}$	$H{n}^{-2}{D}^{-2}$	-	4.78	4.78	4.78	4.78
Speed factor, $n_{ED}$	$nD{H}^{-0.5}$	-	0.46	0.46	0.46	0.46

shows the operating conditions that were used to perform the simulations. These conditions were translated to mass flow inlet and pressure outlet boundary conditions. The atmospheric pressure from the region was used to define the outlet pressure.

The characteristics of the sediments used for this study are shown in

Table 4. Sediment characteristics.

Characteristic	Unit	Value
Density	kg m${}^{-3}$	2650
Size	$\mathsf{\mu }$m	62
Average concentration	kg m${}^{-3}$	0.334

. Particle mass flow rate was calculated as a function of the volumetric water flow rate of each case using the average particle concentration.

3.4. Mesh

The mesh was generated using the Turbogrid module, which employs a high-fidelity hexahedric structured mesh with a uniform distribution. $y^{+}$ values were calculated for stay vanes, guide vanes, and runner blades through the following equation:

$$y^{+}=\frac{u_{t}y}{\nu }$$

(30)

where $u_{\tau }$ is the friction velocity, y is the distance to the nearest wall, and $\nu$ is the kinematic viscosity. The first cell height of each domain was calculated to correctly compute the boundary layer with the selected turbulence model, as shown in Figure 3, by applying a locally refined region near the domain walls. The obtained $y^{+}$ values ranged between 45 and 125, which are considered appropriate for the selected turbulence model [31]. The quality of the mesh was also evaluated using the orthogonal quality model. The orthogonal quality of a cell was estimated as follows.

$$min\left(\frac{\overrightarrow{A_i}\overrightarrow{C_i}}{|A_i\left|\right|C_i|},\frac{\overrightarrow{A_i}\overrightarrow{f_i}}{|A_i\left|\right|f_i|}\right)$$

(31)

where:

• $A_i=$face normal vector;
• $f_i=$vector from the centroid of the cell to the centroid of the face;
• $C_i=$vector from the centroid of the cell to the centroid of the adjacent cell.

A minimum orthogonal quality of 0.269 was obtained in 0.0000002% of mesh elements. The average orthogonal quality in all meshes was 0.93. The walls of the fully structured mesh are shown in Figure 4

Additionally, mesh independence studies were performed considering the pressure drop between the domain inlet and outlet. The structured mesh distribution was modified in all domain directions. Three different mesh resolutions were analyzed for each of the four individual operating conditions evaluated in the study. The results shown in Figure 5 indicate that the difference between the calculations of the fine and medium mesh resolutions is negligible.

3.5. Solver

The numerical simulation was performed using the commercial software ANSYS Fluent. The present simulation used a RANS approach for the liquid phase through a realizable $k-ϵ$ model. The dispersed phase was estimated using the discrete phase model in the commercial CFD software ANSYS Fluent. The numerical simulation was performed coupling all the subdomains with the following imposed boundary conditions:

• The total mass flow inlet was designated at the inlet of the stay vane and the nonoverlapping interfaces of the runner;
• A rotating frame was applied for the runner domain, and other regions were considered to be in a stationary frame;
• Periodic repeats interfaces were created between the stay vane outlet and guide vane inlet and between the outlet of the runner and the inlet of the outlet domain;
• A standard interface was created between the guide vane outlet and the inlet of the runner;
• The total pressure was designated at the outlet of the runner;
• Solid walls were set as nonslip boundary conditions.

A setup used previously in [32] was used for the solid phase, where the injection was applied at the domain inlet and fully elastic collision was assumed at the walls. An analysis of the adequate number of particle injections was also carried out to generate a statistically meaningful sampling. One-hundred stochastic tracking tries were determined to be adequate to ensure that erosion on the walls of the turbine flow passage was independent of the number of injected particles.

The steady-state simulations were carried out using spatial derivatives discretized through a second-order upwind scheme. Full pressure–velocity coupling was enabled using the SIMPLE algorithm. Further, double precision was considered to improve the computational accuracy. A quantitative assessment of the discharge difference was made between the inlet and the outlet, which was lower than the order of $9\times {10}^{-3}$.

The postprocessing phase was carried out in ANSYS CFD-Post, obtaining estimations of the erosion rate on the surfaces of the studied components. Additionally, the pressure and velocity of the flow were determined. Turbine efficiency was calculated based on these results.

3.6. Validation

The study was validated by reproducing the numerical experiment of Nguyen [28], where a wet erosion test rig was used to discharge and project sand particles. A stainless steel plate specimen with a 196 Vickers’ hardness was used for the experiment. The parameters of the experimental setup are shown in

Table 5. Details of the experimental setup for the validation.

Parameter	Units	Value
Particle velocity	m s${}^{-1}$	30
Particle diameter	$\mathsf{\mu }$m	150
Nozzle diameter	mm	6.4
Plate dimensions	mm	25 × 25 × 5
Standoff distance	mm	12.7

.

The results shown in Figure 6 were obtained from the numerical assessment. A satisfactory agreement between the erosion pattern and the results of the experimental and computational tests was observed.

The chart shows the material removal in the specimen caused by sediment erosion, where the center of the chart is aligned with the center of the nozzle. An inverted “W” shape was obtained for the erosion pattern caused by an expected stagnation point in the zone directly below the nozzle. The highest erosion rate was observed right outside this stagnation zone.

4. Results and Discussion

This section presents the results from the sediment analysis and the CFD analysis that were carried out to study the effect of sediment particles over the main components of SFH turbines when varying the guide vane opening. The simulation results were compared with actual site data.

4.1. Sediment Characterization

4.1.1. Particle Size

The analysis showed that 99.58% of the particles at the desilting basin were finer than 425 μm. In addition, the largest percentage (62.33%) of particles was finer than 75 μm. Figure 7 shows the particle size distribution of the sediment samples. The median grain size was determined as 9.28 μm. Furthermore, the characterization of the particle density was performed, obtaining a value of 2650 kg/m³.

4.1.2. Mineral Composition

Figure 8 shows the proportion of the mineral content in sediment samples. The results showed that plagioclase minerals represented the highest proportion of sediments in the samples. The hardness values for this mineral group lies between 6.5 and 7.5 in the Mohs scale. When comparing this mineral composition to other reports, a difference in the proportion of quartz and plagioclase sediments was observed. Quartz is typically the predominant mineral found in most rivers, while in this study, plagioclase minerals composed the greatest part of the sediments found in Pastaza River. Nevertheless, the hardness values did not seem to differ significantly from other studies [10,11,12,14].

4.1.3. Particle Shape

Figure 9 presents the shape of the sediments found in the samples. The analysis exhibited the sharp and slightly rounded edges of the sample particles, which are equivalent to angular and subangular particles based on IEC 62364 standards [33].

Comparing the shape of the sediments with the sample chart, the sphericity and roundness were found to be $S=0.7$ and $R=0.3$, respectively.

4.1.4. Sediment Concentration

Table 6. Particle concentration in Pastaza River (kg m⁻³).

Maximum	Minimum	Average
0.436	0.125	0.334

shows particle concentration values in the river under normal conditions. The average concentration value was used to determine the particle mass flow rate.

4.2. Flow Field Prediction

Since erosion is governed by the velocity, incidence angle, and concentration of the solid particles at the time of collision, the erosion prediction depends on the solutions of these parameters. Figure 10 presents the flow field in the turbine, where the main parameters that influence the sediment flow field are the inlet flow and the guide vane opening. In this context, the highest velocity of the flow was observed on the pressure side near the leading edge for the stay vanes. On the other hand, guide vanes and runner blades presented higher flow velocity on the suction side and trailing edge.

The efficiency of the turbine was calculated at four different operating points using the data from the flow numerical solution. These results were compared with the experimental data from SFH, as shown in Figure 11.

Satisfactory agreement between results was obtained, though a better prediction of the efficiency was obtained for higher flow rates.

4.3. Sediment Erosion Prediction

The sediment erosion patterns in critical turbine components are presented in Figure 12. As seen in Figure 12a, runner blades presented a higher erosion rate on the suction side near the trailing edge. In a similar manner, the guide vanes presented higher sediment erosion on the trailing edge of the suction side, as evidenced in Figure 12b. Comparing the flow field from Figure 10 with the eroded zones in Figure 12, it is evident that zones with the highest relative flow velocity experienced higher erosion rates, as expected.

4.3.1. Erosion Description: San Francisco Francis Turbines

This section compares actual site erosion damage on the turbine with the results from the CFD analysis. The images on the right of Figure 12 show the eroded components of a Francis turbine in SFH. Guide vanes presented a higher erosion rate near the clearance gap around the shaft and at both the leading and trailing edges of the vane since the inward flow accelerates near this region due to the decreasing net head pressure at the guide vane cascade. Regarding the runner, the most eroded areas were located at the leading edge and trailing edge of the blades due to the increase in particle velocity in these areas. Good agreement between numerical erosion results and site erosion was observed when comparing both components. The eroded areas of the runner blades coincided with the areas predicted by the numerical analysis. On the other hand, the eroded areas of the guide vanes did not match perfectly with the numerically predicted areas since the clearance gap near the hub and shroud that forces flow interaction with the shaft was not considered in this study.

4.3.2. Effect of Operating Condition on the Erosion Rate

Figure 13 shows the influence of the guide vane opening and flow velocity on the erosion rate at the stay vanes, guide vanes, and runner. The erosion rate on the walls is calculated as:

$$Er=\frac{1}{A}\int E{r}_{f}dA$$

(32)

where $E{r}_f$ is the facet value of the erosion rate and A is the cell area. The erosion increases when increasing guide vane opening since the velocity and the amount of particles impacting the walls is multiplied due to the rise in water flow. A dramatic increase in the erosion rate was observed for operating conditions with guide vane openings over 90%. The sudden increase in erosion rate observed past a certain operating point may be related to the increase in the intensity of turbulent vortices near the outlet of the blade. Previous works [34,35] have found a direct relation among turbulent flow, vortex formation, and accelerated erosion. Vortices and recirculation accelerate particles in the flow and change the impingement angle to critical values. Figure 14 shows the turbulence intensity of the flow surrounding the blade for the different operating conditions.

5. Conclusions

A CFD study replicating the operating conditions of the Francis turbines of San Francisco hydropower plant in Pastaza, Ecuador, was carried out. Flow conditions and erosion patterns were studied for different performance points, obtaining a detailed prediction of wear damage in different turbine components. From the results of the numerical analysis, the following can be concluded.

Erosion damage increases significantly for higher flow rates, when the opening of the guide vane exceeds an 85% aperture considering the closed position as a reference.

Operating the turbines at the previously mentioned conditions would result in unnecessary and accelerated erosion wear since the best performance point was obtained at a lower flow rate.

The operation of Francis turbines in sediment-laden rivers should be carried out with particular consideration of the effect that guide vane opening has on the formation of turbulent flow and vorticity. This situation can lead to accelerated erosion rates since vortices and recirculation can accelerate particles in the flow and change the impingement angle to critical values.

CFD is a powerful tool that can be used to prevent such occurrences and analyze the operating conditions in hydropower plants that best harness the available power without sacrificing mechanical integrity.

This study was conducted with the aim of contributing to the creation of a clear and cost-effective strategy to prevent and reduce erosion in existing hydropower plants and proposing an effective erosion-based operating procedure for Francis turbines in the Andean region.