# Élaboration d'un modèle d'usure par érosion de cavitation à partir de simulations d'écoulements cavitants

Laurent Krumenacker, Regiane Fortes Patella, Antoine Archer

► To cite this version:

Laurent Krumenacker, Regiane Fortes Patella, Antoine Archer. Élaboration d'un modèle d'usure par érosion de cavitation à partir de simulations d'écoulements cavitants. 21e Congrès Français de Mécanique, 2013, Bordeaux, France. hal-02510597

HAL Id: hal-02510597

https://hal.science/hal-02510597v1

Submitted on 19 Mar 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Élaboration d'un modèle d'usure par érosion de cavitation à partir de simulations d'écoulements cavitants

L. Krumenacker<sup>a</sup>, R. Fortes-Patella<sup>a</sup>, A. Archer<sup>b</sup>

a. Laboratoire des Ecoulements Géophysiques et Industriels

b. EDF R&D, Chatou

## Résumé :

Cette étude présente une nouvelle méthode de prédiction de l'érosion de cavitation par voie numérique. Celle-ci repose sur un modèle (sans couplage retour avec le comportement du matériau) appliqué lors du post-traitement de simulations RANS instationnaire en régime cavitant. A partir des grandeurs du calcul numérique, une répartition de bulles est considérée dans chaque maille située dans les zones de cavitation en sortie de poche et des équations de type Rayleigh-Plesset sont ensuite utilisées pour calculer les ondes de surpression générées par leur implosion. Les calculs sont réalisés sur des géométries d'hydrofoil et comparés à la fois d'un point de vue écoulement (taille de poche, pression en entrée, vitesse de l'écoulement) mais aussi marquages (zone d'érosion, taille et formes des impacts) à des résultats expérimentaux réalisés au LMH-EPFL (Ecole Polytechnique de Lausanne) [4].

## Abstract :

This study presents a new method to predict cavitation erosion. It is based on the post-treatment of URANS simulations in cavitant flows. Bubbles are initialized in the cavitation area from CFD data and the equations of Rayleigh-Plesset are used to calculate the magnitude of pressure waves generated during the collapse. The calculations are performed on hydrofoil geometries and compared with experimental works from LMH-EPFL concerning flows characteristics (size of pocket, pressure, velocity) as well as pitting (erosion area, size and forms of the pits).

**Mots clefs :** Cavitation, Erosion de cavitation, URANS

## 1 Introduction

D'un point de vue industriel, la cavitation peut paraître au sein des circuits hydrauliques au niveau des diaphragmes, vannes ou venturi mais aussi des turbomachines, sur les aubes des inducteurs, des pompes, au niveau de l'entrefer et sur les hélices des navires. Elle donne lieu à des phénomènes de grossissement et d'implosions de structures de vapeur dont la taille et les vitesses caractéristiques varient très rapidement sur des intervalles de temps très réduits. Les ondes de surpressions générées alors peuvent endommager la surface des matériaux.

Cette étude a pour objectif de proposer une nouvelle méthode de prédiction de l'érosion de cavitation à partir d'un post-traitement de calculs numériques d'écoulements cavitants. La méthode prend en compte la tension de surface et la quantité d'air dissous dans le fluide en considérant des bulles sphériques dans l'écoulement.
## 2 Présentation de la méthode

La méthode de post-traitement repose sur une modélisation d'un champ de bulles sur la base des valeurs macroscopiques obtenues numériquement dans chaque maille du domaine de calcul. Ces bulles seront ensuite soumises à une fluctuation de pression importante en sortie de poche de cavitation et engendreront des ondes de surpression avant de se condenser totalement (voir figure 1).

<figure><img src="image_2.png" /><figcaption>FIGURE 1 – Présentation de l'approche considérée</figcaption></figure>

Pour quantifier ces ondes de surpressions des équations de type Rayleigh-Plesset [1] (voir équation 1) seront utilisées.

$$R \frac{d^2 R}{dt^2} + \frac{3}{2} \left( \frac{dR}{dt} \right)^2 = \frac{p_R - p_\infty}{\rho_l} \quad \text{et} \quad p_R = p_{vap} + p_{air} - \frac{2\sigma}{R} \qquad (1)$$

- R est le rayon de la bulle.

- $p_{vap}$ représente la pression de la vapeur, $p_{air}$ celles des autres gaz occlus dans la bulle (dioxygène, diazote,...) considérés comme étant de l'air et $p_\infty$ est la pression dans le liquide située à l'infini.

- $\frac{2\sigma}{R}$ représente le terme de tension superficielle, avec $\sigma$ la tension superficielle air/eau.

### 2.1 Code homogène

Pour cette étude, le code de calcul "IZ", développé pour le compte du CNES au sein du LEGI est employé sur des géométries de type hydrofoil [4]. L'approche homogène dans les calculs cavitants permet de calculer l'écoulement en considérant non plus deux phases, une liquide et une vapeur mais un seul fluide. Le taux de vide $\alpha$ est le rapport entre le volume de la phase gazeuse et le volume total.

$$p_m = \alpha p_{gaz} + (1 - \alpha) p_{liq} \quad , \quad \rho_m = \alpha \rho_{gaz} + (1 - \alpha) \rho_{liq} \quad \text{et} \quad u_{gaz} = u_{liq} = u_m \qquad (2)$$

La conservation de l'équation de quantité de mouvement moyennée pour un modèle homogène s'écrit d'après Ishii [6] :

$$\frac{\partial \rho_m u_{m,i}}{\partial t} + \frac{\partial \rho_m u_{m,i} u_{m,l}}{\partial x_l} = - \frac{\partial p_m}{\partial x_i} + \frac{\partial \sigma_{m,il} - \tau_{m,il}}{\partial x_l} + \sum_k M_{k,i}$$
En considérant que la phase vapeur est constituée uniquement de bulles sphériques, le terme d'interface peut alors s'écrire $\sum_k M_{k,i} = \frac{\partial}{\partial x_i} \left( \alpha \frac{2\sigma}{R} \right)$. Le terme de tension superficielle peut être considéré dans l'équation de quantité de mouvement pour le modèle homogène en notant $\hat{p}_m = p_m - \alpha \frac{2\sigma}{R}$:

$$ \frac{\partial \rho_m u_{m,i}}{\partial t} + \frac{\partial \rho_m u_{m,i} u_{m,l}}{\partial x_l} = - \frac{\partial \hat{p}_m}{\partial x_i} + \frac{\partial \sigma_{m,il} - \tau_{m,il}}{\partial x_l} \quad (3) $$

Enfin pour fermer les équations, on utilise une loi adaptée de la loi barotrope proposée par [2].

$$ \alpha = \frac{1}{2} \left[ 1 - \sin \left( \frac{\hat{p}_m - p_{sat}}{c_{min}^2} \frac{2}{\rho_{liq} - \rho_{gaz0}} \right) \right] \quad (4) $$

Où $\rho_{gaz0}$ est la masse volumique pour la phase pure gazeuse.

## 2.2 Modélisation de la pression des phases

L'eau industrielle n'est pas pure et contient des gaz dissous. Or certains auteurs s'accordent à penser que ces gaz ont un rôle dans le mécanisme d'érosion de cavitation. Pour faciliter l'étude, on considère que ces gaz dissous peuvent être représentés par de l'air (mélange de diazote et d'oxygène dans les mêmes proportions que dans les conditions atmosphériques).

De plus, on suppose que la loi des gaz parfaits est applicable et que la loi de Henry [5] est respectée entre les phases liquide et gazeuse, $C_{air} = H_{air}p_{air}$ où $H_{air}$ est la constante de Henry. En notant $C_{air0}$ la masse volumique d'air total contenue dans la phase liquide et dans la phase gazeuse, on obtient alors :

$$ p_{air} = \frac{C_{air0}}{\left( \frac{\alpha}{r_{air}T} + H_{air} (1 - \alpha) \right)} \quad (5) $$

La masse volumique d'air total est supposée constante dans l'ensemble du calcul.

Comme dans l'équation de Rayleigh-Plesset, la condensation/évaporation à l'interface de la bulle est considérée comme instantanée. La pression de vapeur est donc toujours égale à la pression de vapeur saturante $p_{vap} = p_{sat}$.

Enfin, les bulles sont supposées à l'équilibre mécanique dans la poche de cavitation donnant la relation :

$$ p_{liq} = p_{vap} + p_{air} - \frac{2\sigma}{R} \quad (6) $$

Ainsi, à partir des équation 2 et 6, la pression de la phase liquide est obtenue à partie de la pression du mélange et est considérée comme étant la pression à l'infini dans les calculs de type Rayleigh-Plesset.

$$ p_{liq} = \hat{p}_m = p_{\infty} \quad (7) $$

## 2.3 Implosion des bulles et endommagement du matériau

Avec ces hypothèses, les résultats des calculs numériques nous permettent d'initialiser dans chaque maille du calcul une distribution de bulles de rayon $R = \frac{2\sigma}{p_{liq} - p_{vap} - p_{air}}$. Chaque maille contient $N_{bulle} = \frac{\alpha V_{maille}}{\frac{4}{3}\pi R^3}$.

Ces bulles se déplacent à la même vitesse que l'écoulement et subissent alors une dérivée temporelle de la pression, approximée par sa composante convective, moteur de l'équation de Rayleigh-Plesset. Seuls
les gradients de pression permettant une forte l'implosion des bulles sont considérés, limitant l'étude à la périphérie des structures de vapeur. Par souci de simplicité, on considère que malgré les fortes pressions atteintes, la loi des gaz parfaits est toujours applicable. De même, le temps d'implosion des bulles est supposé suffisamment court pour pouvoir négliger les échanges à l'interface.

$$ \frac{dp_{liq}}{dt} = \frac{\partial p_{liq}}{\partial t} + u_{liq,i} \frac{\partial p_{liq}}{\partial x_i} \approx u_{liq,i} \frac{\partial p_{liq}}{\partial x_i} $$

Dans le modèle actuel, chaque bulle implose et marque la paroi de façon individuelle mais un comportement d'ensemble, considérant le nuage de bulles, est en cours d'étude. On note $P_{max}$ et $\delta t$, les pressions maximales et les temps de passage à mi-amplitude de l'onde de surpression résultant de l'implosion des bulles modélisées. En notant r, leur distance entre leur centre et la paroi, l'énergie de l'onde et le volume de l'impact sont obtenus par [3]. L'endommagement du matériau est alors considéré comme proportionnel à l'énergie de l'onde, via un coefficient $\beta$ déterminé expérimentalement par des essais de références [3]. Nous définissons ainsi une relation entre l'écoulement fluide et la sollicitation du matériau :

$$ E_{onde} = \frac{4\pi r^2}{\rho c_{\infty}} P_{max}^2 \delta t = \beta V_{trou} \quad (8) $$

Avec $c_{\infty} = 1480m.s^{-1}$, la vitesse du son dans la phase liquide et $\beta_{aluminium} = 4J/mm^3$.

En notant $Nb_s$, le nombre de bulles par seconde sortant d'une cellule CFD et générant ainsi une onde de surpression, le volume des empreintes créées sur la paroi la plus proche peut être additionné. Le taux d'endommagement surfacique par seconde $V_d$ ($\mu m^3/mm^2/s$) est alors :

$$ V_d = \sum \frac{Nb_s V_{trou}}{A_{cellule_{paroi}}} \quad (9) $$

## 3 Simulation et essais

L'hydrofoil considéré est un profil NACA 65012 en acier inoxidable, la zone de marquage est quant à elle, en aluminium. La corde est de 100 mm et l'envergure de 150 mm, le profil est placé avec une incidence de 6 degrés. La veine d'essai possède une section carrée de 150x150 mm² et une longueur de 750 mm. Les simulations reprennent le profil 2D de cette géométrie. Deux incidences ont été étudiées, une à 4 degrés, l'autre à 6 degrés.

<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Maillage NACA 65012 - incidence 4 degrés</figcaption></figure>

La campagne d'essais expérimentaux a été menée par Pereira [4] à l'École Polytechnique Fédérale de Lausanne au laboratoire de machines hydrauliques. Outre l'incidence, les deux paramètres d'essais sont la vitesse $v$ moyenne dans la veine et le paramètre de cavitation d'entrée.

$$ \sigma_{in} = \frac{p_{in} - p_{sat}}{\frac{1}{2} \rho_l v^2} $$
**3.1 Incidence 4 degrés**

Avec une incidence à 4 degrés, et une vitesse de $v = 20 \ m.s^{-1}$ et $\sigma = 1.1$, la poche de cavitation n'est pas stable comme le montre la figure 3. La zone marquée par la cavitation obtenue numériquement, centrée sur 50 % de la corde, est légèrement décalée par rapport aux résultats expérimentaux, centrée sur 40 %. L'amplitude du marquage a également plusieurs ordres d'écart par rapport aux essais. (figure 4)

<figure><img src="image_4.png" /><figcaption>FIGURE 3 – taux de vide maximale sur la corde en fonction du temps pour la configuration 4 degrés $v = 20m.s^{-1}$, $\sigma_{in} = 1.1$.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIGURE 4 – Résultats d'expériences pour le profil à 4 degrés d'incidence pour une vitesse $v=20 \ m.s^{-1}$ et un paramètre de cavitation $\sigma_{in} = 1,1$.</figcaption></figure>

**3.2 Incidence 6 degrés**

Les résultats de [4] pour le profil à 6° avec une vitesse de $15 \ m.s^{-1}$ et un paramètre de cavitation de 1.86 présentent une zone d'érosion autour de 20% de corde et autour de 50% pour un paramètre de cavitation de 1,55 (voir figure 5 et 6). L'écoulement présentant plus de lâchers, la résolution numérique est plus délicate à obtenir mais l'on peut noter une bonne approximation de la zone d'érosion.

**4 Conclusions**

La méthode de prédiction d'endommagement proposée présente l'intérêt de prendre en compte le terme de tension superficielle, ainsi que de considérer les gaz dissous contenus dans l'eau. Les résultats qualitatifs concernant notamment les zones érodées obtenues numériquement coïncident avec celles obtenues expérimentalement par Pereira [4]. Néanmoins, avec les paramètres du modèle considérés lors de cette première phase de l'étude, les amplitudes du taux d'endommagement surfacique sont fortement sous estimées si comparées aux résultats expérimentaux. On souligne le fait que les incertitudes de mesure n'ont pas pu être quantifiées dans cette étude et que des travaux expérimentaux sont développés actuellement dans l'équipe pour obtenir une base de données expérimentale plus fiable afin de caler et valider les modèles physiques proposés.
<figure><img src="image_6.png" /><figcaption>FIGURE 5 – Résultats d'expériences pour le profil à 6 degrés d'incidence pour $v=15 \ m.s^{-1}$ et $\sigma_{in} = 1,86$</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>FIGURE 6 – Résultats d'expériences pour le profil à 6 degrés d'incidence pour $v=15 \ m.s^{-1}$ et $\sigma_{in} = 1,5$</figcaption></figure>

L'évolution du modèle en vue de prendre en compte une éventuelle amplification d'impacts liée à l'effet d'ensemble de collapsus de bulles est envisagée et devrait permettre d'obtenir d'ondes de surpressions plus fortes capables d'endommager suffisamment le matériau.

## Remerciements

Les auteurs tiennent à remercier EDF R&D et le CETIM pour avoir financé cette étude.

## Références

[1] Blake and al. (auteur original Prosperetti) 1994 *Bubble Dynamics* : Some things we did not know 10 years ago *Bubble Dynamics and interface phenomena*

[2] Coutier-Delgosha O., Reboud J-L., Delannoy Y. 2003 Numerical simulation of the unsteady behaviour of cavitating flows *Int. J. for Numerical Meth. In fluids* **42** 527-548

[3] Fortes-Patella and al. 2013 Energy Balance in Cavitation Erosion : From Bubble Collapse to Indentation of Material Surface *Journal of Fluids Engineering* **135**

[4] Pereira, F., Avellan F., Dupont J.M. 1998 Prediction of Cavitation Erosion : an Energy Approach *Journal of Fluid Engineering, Transactions of the ASME* **120**

[5] Sander 1999 Compilation of Henry's Law Constants for Inorganic and Organic Species of Potential Importance in Environmental Chemistry http://www.mpch-mainz.mpg.de/sander/res/henry.html

[6] Yoon and Ishii and Revankar 2006 Choking flow modeling with mechanical and thermal nonequilibrium *Int. Journal of Heat and Mass Transfer* **49** 171-186