r4.07.04 Couplage fluide-structure pour les structures tubulaires et les coques coaxiales#
Résumé :
Ce document décrit les différents modèles de couplage fluide-structure disponibles à partir de l’opérateur CALC_FLUI_STRU. Ces modèles permettent de simuler les forces de couplage fluide-élastique dans les configurations suivantes :
faisceaux de tubes sous écoulement transverse (essentiellement pour les tubes de GV),
passage tige de commande / plaque de logement (exclusivement pour les grappes de commande),
coques cylindriques coaxiales sous écoulement annulaire (espace cuve / enveloppe de cœur, …),
faisceaux de tubes sous écoulement axial (assemblages combustibles, …).
Pour chaque configuration, le modèle de forces fluide-élastiques est d’abord présenté. La résolution du problème modal est ensuite décrite. Les méthodes de résolution employées intègrent les spécificités des différents modèles de forces fluide-élastiques.
Excitation fluide-élastique agissant sur les faisceaux de tubes sous écoulement transverse (essentiellement pour les tubes de GV)#
Deux méthodes de simulation de l’excitation fluide-élastique sont disponibles dans Code_Aster .
La première remonte à la fin des années 70. Elle est très répandue dans la communauté scientifique, au sein de laquelle elle est connue sous la dénomination de « méthode de Connors ». Les résultats fournis par cette méthode dépendent en grande partie de la valeur que l’on attribue à l’un de ses principaux paramètres d’entrée : la « constante de Connors ». Des valeurs conservatives ont donc du être déterminées pour cette constante sur la base de nombreux essais réalisés dans le monde. La méthode de Connors est bien adaptée au dimensionnement des faisceaux tubulaires contre le risque vibratoire au stade de la conception. Elle est décrite ci-après dans le paragraphe § 2.5.
La seconde intègre davantage de physique que la méthode de Connors. Cependant, la modélisation complète des phénomènes étant trop complexe par rapport aux connaissances actuelles, cette deuxième méthode demeure basée sur un ensemble de corrélations expérimentales, dites corrélations fluide-élastiques. L’intégration de ce deuxième modèle d’excitation fluide-élastique dans le Code_Aster a été abordée dans la note de spécifications [bib.1]. La note de principe du logiciel FLUSTRU [bib. 3] constitue la documentation théorique de référence. Elle est rappelée dans ses grandes lignes dans les paragraphes § 2.1 à 2.4 ci-après.
Description de la configuration étudiée#
On considère un faisceau de tubes excité par un écoulement externe transverse. Les écoulements externes transverses ont tendance à déstabiliser le système mécanique lorsque la vitesse de l’écoulement augmente. Un cas industriel est celui des vibrations des tubes de GV. Sur ce composant, les écoulements transverses sont observés en entrée du faisceau de tubes (écoulement monophasique liquide), et dans la partie cintrée des tubes (écoulement diphasique) [Figure 2.1-a].
Figure 2.1-a : Schéma de générateur de vapeur
Du point de vue du couplage fluide-élastique, l’étude du comportement dynamique des différents tubes d’un faisceau soumis à un écoulement transverse est ramenée à l’étude d’un tube équivalent ; la définition du tube équivalent dépend de l’environnement du tube à traiter.
Lorsque le tube considéré possède des caractéristiques vibratoires sensiblement différentes de celles de ses voisins, ce tube peut être assimilé à un seul tube, vibrant au milieu d’un faisceau de tubes rigides.
Dans le cas contraire, le problème est plus complexe car on doit considérer un système mécanique avec couplage entre tubes du faisceau et comportant donc un grand nombre de degrés de liberté.
Pour traiter ce genre de configuration, un modèle a été développé au Département TTA, « le modèle global » [bib. 7] ; ce modèle permet la définition d’un système équivalent à un degré de liberté, qui représente le système couplé complet.
L’approche retenue pour conduire les calculs peut être résumée de la façon suivante [Figure 2.1-b] :
Compte tenu de la nature filaire des structures à étudier, le calcul du couplage fluide‑élastique dans le faisceau de tubes est réalisé en décrivant le tube par son abscisse curviligne.
Dans le calcul, l’environnement fluide du tube est caractérisé, à la fois par les propriétés physiques du fluide circulant à l’intérieur du tube (fluide primaire), et par celles du fluide circulant à l’extérieur du tube (fluide secondaire excitateur). Ces propriétés physiques, telle que la masse volumique, peuvent varier le long du tube, en fonction de l’abscisse curviligne.
La vitesse d’écoulement prise en compte pour le calcul de couplage fluide-élastique est la composante, normale au tube dans le plan du tube, de la vitesse du fluide secondaire. Cette vitesse peut varier le long du tube.
Afin de pouvoir prendre en compte les divers types possibles d’excitation, plusieurs zones d’excitation peuvent être définies le long de la structure. Dans le cas du générateur de vapeur, par exemple, on a intérêt à distinguer, d’une part les zones où l’excitation est exercée par un fluide à l’état monophasique, qui se situent en pied de tube, et d’autre part, la zone où l’excitation est diphasique à fort taux de vide, localisée dans la partie cintrée du tube.
Le calcul de couplage est réalisé à partir des caractéristiques mécaniques de la structure en « fluide au repos ». Les forces fluide-élastiques de couplage sont estimées à partir de corrélations adimensionnelles qui sont obtenues sur des expériences analytiques en similitude. Sur chaque zone d’excitation, on peut ainsi appliquer les corrélations adéquates ; les zones d’excitation doivent être disjointes.
Figure 2.1-b : représentation de la configuration à étudier
Pour cette configuration de couplage fluide-élastique, les notations suivantes seront utilisées :
\(L\) |
Longueur totale du tube |
\({L}_{k}\) |
Longueur de la zone \(k\) |
:math:`` |
Diamètre extérieur du tube |
\({d}_{i}\) |
Diamètre intérieur du tube |
\({\varphi}_{i}\) |
Déformée modale du mode \(i\) |
\({\rho}_{i}(x)\) |
Masse volumique du fluide externe à l’abscisse curviligne \(x\) |
\({\rho}_{i}(x)\) |
Masse volumique du fluide interne à l’abscisse curviligne \(x\) |
\({\rho}_{t}\) |
Masse volumique du tube (structure seule) |
\({\rho}_{\mathit{eq}}(x)\) |
Masse volumique équivalente à l’abscisse curviligne \(x\) |
\(U\) |
Vitesse du fluide externe spécifiée par l’utilisateur dans l’opérateur DEFI_FLUI_STRU |
\(V(x)\) |
Vitesse du fluide externe à l’abscisse curviligne \(x\) |
\({V}_{k}(x)\) |
Vitesse du fluide externe à l’abscisse curviligne \(x\) (zone d’excitation \(k\) ) définie par le produit de \(U\) et d’un profil de vitesse spécifiée par l’utilisateur dans l’opérateur DEFI_FLUI_STRU |
\({U}_{k}\) |
Vitesse moyenne du fluide externe calculée à partir de \({V}_{k}(x)\) pour la zone d’excitation \(k\) |
\(\overline{U}\) |
Moyenne des vitesses \({U}_{k}\) sur toutes les zones d’excitation |
Étapes du calcul#
La première étape du calcul consiste à calculer les caractéristiques de la structure en « fluide au repos ». On procède en considérant une masse équivalente du tube ; cette masse équivalente regroupe, d’une part la masse du tube seul, et d’autre part les masses ajoutées par les fluides interne et externe.
Une masse volumique équivalente est ainsi définie le long du tube en fonction de l’abscisse curviligne
par l’expression :
\({\rho}_{\mathit{eq}}(x)=\frac{1}{({d}_{e}^{2}-{d}_{i}^{2})}[{\rho}_{i}(x).{d}_{i}^{2}+{\rho}_{t}.({d}_{e}^{2}-{d}_{i}^{2})+{\rho}_{e}(x).{d}_{\mathit{eq}}^{2}]\) éq. 2.2- 1
avec
\({d}_{\mathit{eq}}^{2}=\frac{2.{C}_{m}.{d}_{e}^{2}}{\pi}\) éq. 2.2- 2
Dans l’équation [:ref:` éq. 2.2-1 < éq. 2.2-1 >`], le terme \({\rho}_{e}(x).{d}_{\mathit{eq}}^{2}\) représente la masse ajoutée par le fluide externe. Ce terme dépend, par l’intermédiaire du paramètre \({C}_{m}\) , de l’arrangement du faisceau de tubes (pas carré ou triangulaire), et du confinement du faisceau (pas réduit). Pour les calculs de couplage fluide-élastique des faisceaux de tubes soumis à un écoulement transverse, on utilise couramment, pour estimer le coefficient \({C}_{m}\) , des expressions analytiques déterminées à partir de résultats expérimentaux. L’ensemble des données nécessaires à l’estimation du coefficient \({C}_{m}\) est recueilli par l’opérateur DEFI_FLUI_STRU.
Connaissant la masse volumique équivalente du tube, les matrices élémentaires de masse et de raideur en eau au repos sont ensuite calculées au moyen du profil de masse volumique équivalente, par l’opérateur CALC_MATR_ELEM ; on utilise les options MASS_FLUI_STRU et RIGI_FLUI_STRU. L’opérateur CALC_MODES permet, après assemblage des matrices élémentaires, de calculer directement les modes en eau au repos de la structure étudiée.
Les forces fluide-élastiques de couplage sont calculées par l’opérateur CALC_FLUI_STRU à partir des corrélations adimensionnelles établies sur des maquettes analytiques en similitude. Ces forces de couplage, \([{B}_{ij}(U,s)]\) , dépendantes du mouvement de la structure sont ensuite prises en compte dans l’équation générale du mouvement [:ref:` éq. 1.2-1 < éq. 1.2-1 >`] pour calculer les caractéristiques du système couplé écoulement-structure pour une vitesse donnée d’écoulement.
Expression de la matrice de transfert des forces fluide-élastiques#
Dans le cas des faisceaux de tubes excités par un écoulement transverse, les forces fluide-élastiques de couplage sont des forces réparties le long de la structure. Elles sont caractérisées par des coefficients adimensionnels linéiques d’amortissement et de raideur ajoutés, dénommés respectivement \({C}_{d}\) et \({C}_{k}\) . L’expression des coefficients de la matrice de transfert des forces fluide-élastiques projetée sur la base modale de la structure en « fluide au repos » est alors la suivante :
\({B}_{ij}(U,s)=[(\underset{L}{\int}\text{}\frac{1}{2}{\rho}_{e}(x)V(x){d}_{e}{C}_{d}(x,{s}_{r}){\varphi}_{i}^{2}(x)\mathit{dx})s+\underset{L}{\int}\text{}\frac{1}{2}{\rho}_{e}(x){V}^{2}(x){C}_{k}(x,{s}_{r}){\varphi}_{i}^{2}(x)\mathit{dx}].{\delta}_{ij}\) éq. 2.3- 1
La dépendance des coefficients \({C}_{d}\) et \({C}_{k}\) vis-à-vis du mouvement de la structure et de la vitesse de l’écoulement du fluide est traduite par leur évolution en fonction de la fréquence réduite complexe , \({s}_{r}\) définie par:
\({s}_{r}=\frac{s.D}{U}\) éq. 2.3- 2
L’expression [:ref:` éq. 2.3-1 < éq. 2.3-1 >`] montre que l’on retient une matrice de transfert diagonale. Cela implique :
les différents modes propres de la structure sont assez éloignés les uns des autres pour que l’on puisse supposer qu’il n’y a pas couplage entre modes.
les déformées modales de la structure en « fluide au repos » ne sont pas perturbées par la mise en écoulement du fluide .
Ces deux hypothèses ont pu être vérifiées expérimentalement sur les faisceaux de tubes soumis à un écoulement transverse.
En pratique, compte tenu des différentes zones d’excitation prises en compte le long de la structure, les coefficients diagonaux de la matrice d’efforts fluide-élastiques projetée sur base modale s’écrivent :
\({B}_{ii}(U,s)=\sum_{k}[(\underset{{L}_{k}}{\int}\text{}\frac{1}{2}{\rho}_{e}(x){V}_{k}(x){d}_{e}{C}_{\mathit{dk}}(\frac{{\mathit{sd}}_{e}\overline{U}}{{\mathit{UU}}_{k}}){\varphi}_{i}^{2}(x)\mathit{dx})s+\underset{{L}_{k}}{\int}\text{}\frac{1}{2}{\rho}_{e}(x){V}_{k}^{2}(x){C}_{kk}(\frac{{\mathit{sd}}_{e}\overline{U}}{{\mathit{UU}}_{k}}){\varphi}_{i}^{2}(x)\mathit{dx}]\) éq. 2.3- 3
où \({C}_{\mathit{dk}}\) et \({C}_{kk}\) désignent respectivement les coefficients de couplage adimensionnels, d’amortissement et de raideur, retenus pour la zone d’excitation \(k\) . La vitesse fluide \(\frac{{\mathit{UU}}_{k}}{\overline{U}}\) intervenant dans la fréquence réduite complexe en argument des coefficients de couplage correspond à la vitesse moyenne sur la zone d’excitation \(k\) , après renormalisation du profil \({V}_{k}(x)\) , de sorte que sa moyenne sur toutes les zones d’excitation vaille \(U\) .
Il est par ailleurs très important de noter que chaque déformée modale prise en compte dans les équations 2.3-1, 2.3-3, etc. l’est en réalité uniquement par l’intermédiaire de sa composante en translation suivant la direction de la portance. Ceci est du au fait que les coefficients d’amortissement et de raideur ajoutés qui apparaissent dans ces équations ont été déterminés (expérimentalement) uniquement pour la direction de la portance. Cette remarque s’applique à toutes les méthodes de calcul d’instabilités fluide-élastiques de tubes de GV présentées dans ce document, y compris à la méthode de Connors présentée aux paragraphes § 2.5.1 et § 2.5.2. Il en résulte notamment que les matrices généralisées de masse, d’amortissement et de raideur qui apparaissent dans les équations associées aux calculs d’instabilité fluide-élastique des tubes de GV (comme par exemple l’équation 2.4-1) ne sont pas des matrices généralisées au sens habituel du terme, c’est à dire s’appuyant sur les trois composantes en translation et sur les trois composantes en rotation, mais des matrices généralisées que l’on peut qualifier d”«orientées suivant une direction privilégiée» dans la mesure où elles sont toutes calculées sur la base de la seule composante en translation des déformées modales selon la direction de la portance. Cette remarque s’applique uniquement à l’application «vibrations des tubes de GV» et, à l’intérieur de cette application, au calcul des instabilités fluide-élastiques.
Résolution du problème modal sous écoulement#
Dans la configuration « Faisceau de tubes soumis à un écoulement transverse », le problème est résolu sur la base modale caractérisant la structure en « fluide au repos ».
D’une manière générale, les caractéristiques du système couplé écoulement-structure sont obtenues en recherchant les solutions de l’équation :
\(\lbrace [{M}_{ii}]{s}^{2}+[{C}_{ii}]s+[{K}_{ii}]-[{B}_{ij}(U,s)]\rbrace (q)=0\) éq. 2.4- 1
où |
\([{M}_{ii}],[{C}_{ii}]\mathit{et}[{K}_{ii}]\) désignent respectivement les matrices diagonales de masse, d’amortissement et de raideur caractéristiques de la structure en “fluide au repos” ; |
\((q)\) désigne le vecteur des déplacements généralisés en “fluide au repos”. |
Comme la matrice d’efforts fluide-élastiques retenue est diagonale, et que les déformées modales sont supposées ne pas être modifiées sous écoulement, le problème de couplage fluide-élastique se ramène à la résolution de \(N\) problèmes scalaires, \(N\) désignant le nombre de modes pris en compte dans la base modale.
Pour chaque mode \(i\) et chaque vitesse d’écoulement \(U\) , le problème à résoudre s’écrit :
\({M}_{ii}{s}^{2}+[{C}_{ii}-\sum_{k}(\underset{{L}_{k}}{\int}\text{}\frac{1}{2}{\rho}_{e}(x){V}_{k}(x){d}_{e}{C}_{\mathit{dk}}(\frac{{\mathit{sd}}_{e}\overline{U}}{{\mathit{UU}}_{k}}){\varphi}_{i}^{2}(x)\mathit{dx})]s+{K}_{ii}-\sum_{k}(\underset{{L}_{k}}{\int}\text{}\frac{1}{2}{\rho}_{e}(x){V}_{k}^{2}(x){C}_{kk}(\frac{{\mathit{sd}}_{e}\overline{U}}{{\mathit{UU}}_{k}}){\varphi}_{i}^{2}(x)\mathit{dx})=0\)
éq. 2.4- 2
On notera que l’équation [:ref:` éq. 2.4-2 < éq. 2.4-2 >`] est non-linéaire en \(s\) ; ses solutions sont obtenues à l’aide d’une méthode itérative de type Broyden.
Pour chaque mode \(i\) , on obtient une solution \({s}_{i}\) de l’équation [:ref:` éq. 2.4-2 < éq. 2.4-2 >`]. On déduit ensuite de \({s}_{i}\) , pour ce mode, la pulsation \({\omega}_{i}\) et l’amortissement \({\xi}_{i}\) du système couplé écoulement-structure, en utilisant la relation :
\({s}_{i}=-{\xi}_{i}{\omega}_{i}+J{\omega}_{i}\sqrt{1-{\xi}_{i}^{2}}\) avec \({J}^{2}=-1\) éq. 2.4- 3
Le système couplé devient dynamiquement instable lorsque l’un des coefficients d’amortissement \({\xi}_{i}\) devient négatif ou s’annule.
Méthode de Connors#
Cas d’une zone unique d’excitation fluide#
En 1978, H.J. Connors propose de déterminer la vitesse critique \({V}_{\mathit{cn}}\) associée au mode d’ordre \(n\) d’un tube de Générateur de Vapeur(GV) suivant la relation[10] :
\(\frac{{V}_{\mathit{cn}}}{{f}_{n}{D}_{e}}=K\sqrt{\frac{\overline{m}{\delta}_{n}}{\overline{{\rho}_{s}}{D}_{e}^{2}}}\)
Dans cette relation:
désigne la vitesse critique inter-tubes d’instabilité pour le mode n,
désigne la fréquence propre d’ordre n du tube [1]
,
désigne le diamètre extérieur du tube,
désigne la constante de Connors,
désigne la masse linéique de référence du tube incluant les effets de masse ajoutée,
désigne le décrément logarithmique du mode n en fluide au repos, c’est à dire incluant l’amortissement de la structure et celui apporté par le fluide externe au repos, et
désigne la masse volumique de référence du fluide secondaire.
On rappelle que le décrément logarithmique
se définit comme :
Où
désigne l’amortissement modal réduit du mode n.
étant de l’ordre du pour cent, il est légitime de poser
, et donc l’approximation :
La constante d’instabilité K est déterminée expérimentalement à partir de résultats d’essais d’instabilité. Dans les études de dimensionnement vibratoire des faisceaux de tubes de GV, les valeurs adoptées usuellement pour cette constante sont :
K = 4 en cas d’écoulement transverse diphasique au niveau du chignon,
K = 2,9 en cas d’écoulement transverse monophasique au dessus de la plaque tubulaire.
En considérant comme masse linéique de référence du tube
sa masse linéique moyenne, on peut déterminer
sous la forme :
Où
désigne le diamètre intérieur du tube,
désigne sa longueur, s désigne l’abscisse curviligne le long du tube et
désigne la masse volumique équivalente du tube à l’abscisse s :
Où
désigne la masse volumique du tube supposée indépendante de l’abscisse curviligne,
et
désignent respectivement la masse volumique du fluide primaire et du fluide secondaire à l’abscisse curviligne \(s\) , et
est défini par:
, où \(\Delta\) désigne un diamètre équivalent donné par les relations :
pour un pas carré ( C vaut environ 2,0 pour les GV français)
pour un pas triangulaire ( C vaut environ 2,2 pour les GV français)
De même, en considérant comme masse volumique de référence du fluide secondaire
sa masse volumique moyenne, on peut déterminer
sous la forme :
Le mode n est instable si la vitesse critique
est inférieure à la vitesse efficace
associée au mode n, ainsi définie :
Où
désigne la déformée modale du mode \(n\) ,
désigne la vitesse de l’écoulement en fonctionnement (\(m/s\) ),
désignent la masse linéique du tube incluant les effets de masse ajoutée (\(\mathit{kg}/m\) ) supposée varier le long du tube, obtenue comme :
On définit le rapport d’instabilité pour le mode n au sens de Connors comme étant le rapport:
Cas de plusieurs zones d’excitation fluide#
La démarche d’application de la méthode de Connors mérite d’être précisée dans le cas où le tube est soumis à une excitation multiforme de la part du fluide, en particulier, si cette dernière est monophasique en bas de faisceau et diphasique dans le chignon. On rappelle qu’une telle situation est prise en compte dans le logiciel GEVIBUS [11].
Pour extrapoler la méthode de Connors à ce cas général, on procède en généralisant l’établissement de la démarche proposée par Connors [10].
Soit
l’énergie ajoutée par l’écoulement au cours d’un cycle à un tube vibrant dans son mode n:
Où, par rapport à l’exposé de Connors, la dépendance de la constante
à la zone d’excitation \(i\) est ajoutée,
désigne le nombre total de zones d’excitation, et
désigne la longueur de la \(i\) ‑ème zone d’excitation.
Soit
l’énergie dissipée au cours d’un cycle par le tube vibrant dans son mode n:
En égalant
et
, c’est à dire en se plaçant à l’instabilité, en introduisant comme Connors les variables de référence
et
(bien qu’elles ne paraissent pas indispensables), en posant
, où
désigne la constante de Connors associée à la \(i\) -ème zone d’excitation, et en cherchant à faire apparaître la vitesse efficace
telle que Connors la définit, on obtienttous calculs faits l’expression :
D’où:
On en déduit dans le cas d’une excitation multiforme l’expression de la vitesse critique
associée au mode d’ordre \(n\) :
On vérifie que, lorsque l’excitation est de même nature sur l’ensemble des zones excitées
, la relation ci-dessus retrouve la forme proposée par Connors:
Variante de la méthode#
Dans la méthode de Connors présentée dans les paragraphes ci-dessus (paragraphes § 2.5.1 et § 2.5.2), le calcul du rapport d’instabilité ne prend en compte qu’une seule composante des modes: la direction définie dans DEFI_FLUI_STRU. Une variante de cette méthode consiste à prendre en compte les trois composantes en translation des modes. La vitesse efficace et la vitesse critique s’écrivent alors:
\({V}_{\mathit{en}}=\sqrt{\frac{\underset{\mathit{tube}}{\int}(\frac{{\rho}_{s}(s)}{\overline{{\rho}_{s}}}{V}^{2}(s){\mathrm{\phi }}_{n}^{2}(s)\mathit{ds})}{\frac{{M}_{n}}{\overline{m}}}}\)
Où \({M}_{n}\) désigne la masse généralisée (non orientée suivant une direction privilégiée et prenant donc en compte à la fois les trois composantes en translation et les trois composantes en rotation) du mode \(n\) , et \({\varphi}_{n}\) désigne la déformée modale du mode \(n\) . Par \({\mathrm{\phi }}_{n}^{2}\) on entend ici la somme des carrés des trois composantes en translation de \({\varphi}_{n}\) . Les trois composantes en rotation ne sont pas prises en compte dans le calcul de \({\varphi}_{n}^{2}\) .
En posant:
\(r(s)=\frac{{\rho}_{s}(s)}{\overline{{\rho}_{s}}}\) et \(u(s)=\frac{V(s)}{{V}_{\mathit{moy}}}\)
où \({V}_{\mathit{moy}}\) est la vitesse moyenne., le rapport d’instabilité s’écrit :
\({R}_{n}=\frac{{V}_{\mathit{en}}}{{V}_{\mathit{cn}}}=\frac{{V}_{\mathit{moy}}}{{{f}_{n}{D}_{e}K[\frac{2\pi {\xi}_{n}{M}_{n}}{\overline{{\rho}_{s}}{D}_{e}^{2}\underset{\mathit{tube}}{\int}(r(s){u}^{2}(s){\mathrm{\phi }}_{n}^{2}(s)\mathit{ds})}]}^{1/2}}\)
Le calcul du rapport selon cette variante est systématiquement effectué par Code_Aster quand on demande la mise en œuvre de la méthode de Connors.
Le rapport d’instabilité calculé selon cette variante est fourni à côté du rapport d’instabilité calculé selon la méthode précisée aux paragraphes § 2.5.1 et § 2.5.2. La plupart du temps, les deux résultats sont identiques. S’il existe un écart, la raison de cet écart doit être recherchée dans la contribution des composantes en rotation du mode considéré, par exemple dans la contribution des rotations des parties droites des tubes autour de leur axe. On adoptera alors le résultat le plus pénalisant des deux.
Excitation fluide-élastique agissant sur la tige de commande au niveau de la plaque de logement (exclusivement pour les grappes de commande)#
Les forces fluide-élastiques agissant sur ce type de configuration ont été identifiées sur la maquette GRAPPE2 du département TTA. Les aspects théoriques de l’identification de ces sources sont développés en référence [bib. 4]. L’intégration du modèle GRAPPE2 dans Code_Aster est abordée dans la note de spécifications [bib. 2].
Description de la configuration étudiée#
La maquette GRAPPE2 représente la tige de commande, la partie supérieure du guide de grappe, et la manchette thermique d’un réacteur de type 900 ou 1300 MWe [Figure 3.1-a].
Figure 3.1-a: Schéma de principe de la maquette GRAPPE 2
Cette maquette est essentiellement constituée d’un tube cylindrique creux de faible épaisseur, fixé sur une âme centrale cylindrique pleine. Le tube creux est entièrement immergé dans de l’eau à température ambiante. Une plaque, représentant la plaque de logement, permet de reproduire le confinement annulaire. L’écoulement à travers la plaque peut être ascendant ou descendant. La tige de commande peut être centrée ou excentrée (50% du jeu moyen) au niveau de la plaque de logement.
Quatre configurations expérimentales sont donc possibles, en fonction du sens de l’écoulement et du centrage ou pas de la tige de commande. Les coefficients de forces fluide-élastiques ont été identifiés pour chacune de ces configurations et sont disponibles dans le Code_Aster .
La maquette GRAPPE2 a été dimensionnée en similitude géométrique, hydraulique et de fréquence réduite par rapport à la configuration réacteur. La seule donnée du diamètre de la tige de commande permet donc, en particulier, de déduire l’ensemble des autres grandeurs géométriques.
Étapes du calcul#
La première étape du calcul consiste à calculer la base modale de la structure en eau au repos, les effets de masse ajoutée induits localement au niveau du confinement de la plaque de logement étant négligés. Cette étape est réalisée par l’opérateur CALC_MODES.
Pour ce faire, une masse volumique équivalente homogène est affectée à l’ensemble de la structure, afin de prendre en compte la masse apparente ajoutée par le fluide, à l’exception de celle induite par les effets de confinement au niveau de l’espace annulaire. Cette masse volumique équivalente est définie par :
éq. 3.2- 1
où:
|
désigne un coefficient adimensionnel de confinement dépendant de la configuration étudiée ; \(\alpha =1\) est la valeur utilisée pour les calculs de grappes de commande. Elle correspond à un cylindre vibrant dans un domaine fluide illimité. |
|
désigne le rayon extérieur du tube, |
|
désigne l’aire de la section droite du tube, |
|
désigne la masse volumique du matériau constituant le tube vibrant. |
La seconde étape est la prise en compte du couplage avec l’écoulement fluide. Elle est réalisée à l’aide de l’opérateur CALC_FLUI_STRU.
Représentation de l’excitation fluide-élastique#
Soit \(x\) la direction de la fibre neutre du tube. L’excitation fluide-élastique identifiée sur la maquette GRAPPE2 est représentée par une force et un moment résultants, appliqués en un même point d’abscisse
, correspondant à la zone centrale du passage de la tige de commande à travers la plaque de logement. L’excitation est ainsi définie, dans la base physique, par la relation :
éq. 3.3- 1
où
désigne la dérivée par rapport à \(x\) de la distribution de Dirac
.
La force résultante,
, agit ainsi sous l’effet des déplacements transverses de la tige de commande ; et le moment résultant,
, agit sous l’effet de la rotation de cette dernière.
On note
le vecteur des déplacements transverses et
le vecteur des rotations associées, définis par :
éq. 3.3- 2
éq. 3.3- 3
Les relations suivantes sont utiliséespour calculer les forces et les moments fluide-élastiques résultants à partir des masses ajoutées
, des amortissements ajoutés
et des raideurs ajoutées
, coefficients adimensionnels identifiés sur la maquette GRAPPE2 :
éq. 3.3- 4
éq. 3.3- 5
Afin de simplifier l’écriture des équations, on note par la suite:
La vitesse réduite adimensionnelle
est ici définie à l’aide de la relation
, où
désigne la variable de Laplace.
Les expressions [:ref:` éq. 3.3-4 < éq. 3.3-4 >`] et [:ref:` éq. 3.3-5 < éq. 3.3-5 >`] font intervenir l’épaisseur
de la plaque de logement. Cette épaisseur se déduit de la valeur du diamètre de la tige de commande,
, du fait de la similitude géométrique avec la configuration réacteur. L’effort fluide-élastique
est ainsi complètement caractérisé par la donnée des grandeurs suivantes :
|
Masse volumique du fluide, |
|
Vitesse de l’écoulement moyen dans l’espace annulaire entre tige de commande et plaque de logement, |
|
Diamètre de la tige de commande, |
|
Coefficient de masse ajoutée associé au mouvement de translation, |
|
Coefficient d’amortissement ajouté associé au mouvement de translation, |
|
Coefficient de raideur ajoutée associé au mouvement de translation, |
|
Coefficient de masse ajoutée associé au mouvement de rotation, |
|
Coefficient d’amortissement ajouté associé au mouvement de rotation, |
|
Coefficient de raideur ajoutée associé au mouvement de rotation. |
Les coefficients adimensionnels de masse ajoutée, \({\mathit{Cm}}_{1}\) et \({\mathit{Cm}}_{2}\) , permettent la prise en compte des effets inertiels induits par le confinement local de la tige de commande au niveau de la plaque de logement. Ces effets sont estimés comme suit.
Soit
l’épaisseur de l’écoulement annulaire au niveau du confinement, déduite de
par similitude géométrique par rapport à la configuration réacteur ;
désigne le coefficient adimensionnel de confinement introduit par la relation [:ref:` éq. 3.2-1 < éq. 3.2-1 >`]. On obtient alors [bib. 4] :
On en déduit les valeurs de \({\mathit{Cm}}_{1}\) et \({\mathit{Cm}}_{2}\) par :
éq. 3.3- 6
éq. 3.3- 7
Les coefficients
et
sont directement déduits de la mesure et exprimés sous forme de corrélations adimensionnelles.
Projection sur base modale et expression des termes de la matrice de transfert d’effort fluide-élastique#
Décomposition du mouvement sur base modale
On note
la déformée modale du \(j\) ème mode de la structure. La décomposition du vecteur des déplacements dans la base modale s’exprime sous la forme :
éq. 3.4- 1
Où
,
et
correspondent aux trois composantes de translation caractérisant les déformées modales calculées à l’aide de Code_Aster .
Calcul de l’excitation généralisée associée au mode i
L’excitation généralisée
associée au mode
est définie par la relation :
éq. 3.4- 2
où
désigne la longueur de la structure sur laquelle on veut imposer les excitations GRAPPE2.
Les fonctions de transfert
et
étant définies à partir des relations [:ref:` éq. 3.3-4 < éq. 3.3-4 >`] et [:ref:` éq. 3.3-5 < éq. 3.3-5 >`], on en déduit, compte tenu des expressions [:ref:` éq. 3.3-1 < éq. 3.3-1 >`], [:ref:` éq. 3.3-4 < éq. 3.3-4 >`] et [:ref:` éq. 3.3-5 < éq. 3.3-5 >`] :
éq. 3.4- 3
D’où, après intégration :
éq. 3.4- 4
Remarque :
et
Résolution du problème modal sous écoulement#
Le problème modal est résolu en supposant, en première approximation, que les termes diagonaux de la matrice de transfert des efforts fluide-élastiques
sont prépondérants par rapport aux termes extra-diagonaux.
La matrice
étant ainsi réduite à sa diagonale, les déformées modales ne sont pas perturbées par la prise en compte du couplage fluide-élastique ; les seuls paramètres modifiés sont les fréquences propres et les amortissements réduits modaux.
Le problème modal sous écoulement se décompose alors en
problèmes scalaires indépendants, résolus par une méthode de type Broyden :
éq. 3.5- 1
où |
|
désigne la masse généralisée ajoutée par le fluide, |
|
désigne l’amortissement généralisé ajouté par le fluide, |
|
|
désigne la raideur généralisée ajoutée par le fluide. |
,
et
sont calculées à l’aide des relations :
éq. 3.5- 2
éq. 3.5- 3
éq. 3.5- 4
et
dépendent implicitement de
par l’intermédiaire de la vitesse réduite
.
Les trois grandeurs nécessaires pour dimensionner ces termes sont donc uniquement
et
étant déduites de
grâce à la propriété de similitude géométrique.
Comme cela a été indiqué précédemment, les coefficients adimensionnels
et
sont issus des corrélations empiriques identifiées expérimentalement sur la maquette GRAPPE2.
Excitation fluide-élastique agissant sur deux coques cylindriques coaxiales sous écoulement annulaire (exemple : espace cuve / enveloppe de cœur)#
L’intégration de ce modèle d’excitation fluide-élastique dans Code_Aster a été abordée dans la note de spécifications [bib. 2]. La note de principe du modèle MOCCA_COQUE [bib. 5] constitue la documentation théorique de référence.
Description de la configuration étudiée#
La configuration matérielle étudiée est composée de deux coques cylindriques coaxiales, séparées par un espace annulaire dans lequel s’écoule un fluide monophasique incompressible visqueux [Figure4.1‑a]. L’écoulement se fait dans la direction de l’axe de révolution des cylindres ; pour fixer les notations, on suppose dans la suite du document qu’il s’agit de l’axe
.
On note :
|
la longueur commune des deux coques cylindriques, |
|
le rayon intérieur de l’espace annulaire, |
|
le rayon extérieur de l’espace annulaire, |
|
le rayon moyen
|
|
le jeu annulaire |
|
les vecteurs de la base de coordonnées cylindriques. |
Figure 4.1-a : schéma de principe coques coaxiales
Étapes de calcul#
La première étape de calcul consiste à déterminer la base modale en air de la structure. Cette opération est réalisée par l’opérateur CALC_MODES. Ce calcul est nécessaire car la décomposition de la matrice de transfert des forces fluide-élastiques
est exprimée dans cette base.
La seconde étape concerne la prise en compte des forces fluide-élastiques. Elle intervient dans l’opérateur CALC_FLUI_STRU . Cette étape se décompose en huit sous-tâches :
Pré-traitements#
1°/ |
Détermination des grandeurs géométriques caractéristiques, à partir de la topologie du maillage : longueur commune des deux coques, rayon moyen, jeu annulaire moyen. |
2°/ |
Caractérisation des déformées modales en air : détermination des ordres de coque, des plans principaux, des nombres d’onde et des coefficients de déformées de poutre associées à chacun des modes de la structure, tant pour la coque interne que pour la coque externe. |
Résolution du problème modal en eau au repos#
3°/ |
Calcul de la matrice de masse ajoutée par le fluide
|
4°/ |
Calcul des caractéristiques modales de la structure en eau au repos en résolvant :
|
5°/ |
Calcul des déformées en eau au repos dans la base physique, par changement de base: |
Résolution du problème modal sous écoulement#
Pour chaque vitesse d’écoulement:
6°/ |
Calcul de
|
7°/ |
Calcul des forces fluide-élastiques induits par les effets d’amortissement et de raideur ajoutés, dans la base modale en eau au repos. |
8°/ |
Résolution du problème modal en négligeant les termes extra-diagonaux de cette dernière matrice, par la méthode de Broyden (boucle sur les sous-tâches 6° et 7°).
|
Fin de boucle sur les vitesses d’écoulement
Remarques :
Le calcul des termes de la matrice de transfert des forces fluide-élastiques nécessite la résolution du problème fluide instationnaire (sous-tâche 6°). Cette résolution n’est elle‑même possible que si l’on a préalablement déterminé certaines grandeurs géométriques caractéristiques de la configuration, ainsi que les coefficients des formes analytiques des déformées modales des structures (pré-traitements 1° et 2°).
Si l’utilisateur choisit de réaliser la première étape (calcul de la base modale par l’opérateur CALC_MODES) en prenant directement en compte les effets de masse ajoutée, ceux-ci ne doivent plus être pris en compte par l’opérateur CALC_FLUI_STRU. Pour cela, le mot-clé MASS_AJOU de la commande DEFI_FLUI_STRU doit être renseigné par “NON”. Les sous‑tâches 3° à 7° deviennent alors :
3°/ |
Calcul des effets de masse ajoutée par le fluide, dans la base modale de la structure en eau, afin de pouvoir retrancher ces effets de l’effort fluide-élastique global, puisque les termes de masse ajoutée sont déjà pris en compte. |
4°/ |
Sous-tâche supprimée. |
5°/ |
Sous-tâche supprimée. |
Pour chaque vitesse d’écoulement
6°/ |
Calcul de la matrice
|
7°/ |
Calcul des forces fluide-élastiques induites par les effets d’amortissement et de raideur ajoutés dans la base modale en eau : |
Les sous-tâches 1°, 2° et 8° ne sont pas modifiées.
Pour chaque vitesse d’écoulement
6°/ |
Calcul de la matrice
|
7°/ |
Calcul des forces fluide-élastiques induites par les effets d’amortissement et de raideur ajoutés dans la base modale en eau : |
Les sous-tâches 1°, 2° et 8° ne sont pas modifiées.
Résolution du problème fluide instationnaire#
Hypothèses simplificatrices#
Quelques hypothèses sur la nature de l’écoulement permettent de simplifier les équations de Navier-Stokes instationnaires, à la base du problème fluide-structure.
H1 |
On suppose que l’écoulement est la superposition d’un écoulement moyen stationnaire, obtenu lorsque les structures sont fixes, et d’un écoulement instationnaire induit par le mouvement des parois. |
H2 |
On suppose que les vibrations de structure sont de faible amplitude vis à vis de l’épaisseur de l’écoulement annulaire moyen. |
H3 |
On suppose que les perturbations de vitesse induites par les mouvements vibratoires sont, en moyenne sur un rayon, essentiellement dirigées dans les directions
|
H4 |
On suppose enfin que le champ de vitesse et de pression est uniforme, à l’ordre 1, dans la direction radiale. |
Ces hypothèses simplificatrices permettent de résoudre analytiquement le problème fluide. La matrice de transfert des forces fluide-élastiques
est déduite de l’écoulement instationnaire issu de cette résolution.
Analyse en perturbations#
Moyennant les hypothèses énoncées précédemment, l’analyse en perturbations du problème fluide conduit à rechercher l’écoulement instationnaire sous la forme :
éq. 4.3.2-1
éq. 4.3.2-2
éq. 4.3.2-3
éq. 4.3.2-4
avec :
éq. 4.3.2-5
éq. 4.3.2-6
On définit les variables
et
comme :
et
.
En limitant le développement des équations de Navier-Stokes au premier ordre, on obtient deux systèmes d’équations caractérisant la partie stationnaire et la partie perturbée de l’écoulement, le second système étant un système linéaire.
La résolution du problème fluide stationnaire conduit ainsi à :
constante et
éq. 4.3.2-7
Dans l’équation [:ref:` éq. 4.3.2-7 < éq. 4.3.2-7 >`],
désigne la masse volumique du fluide et
la partie stationnaire du coefficient de frottement à la paroi. Le fluide étant supposé incompressible, sa masse volumique n’est pas décomposée en partie stationnaire et partie fluctuante.
est déduit de la loi de Nikuradzé caractérisant les écoulements en conduite :
avec
éq. 4.3.2-8
où m désigne la valeur d’un exposant,
désigne la viscosité cinématique du fluide et
la rugosité des parois.
Il en découle:
avec
et
Le système différentiel linéaire d’ordre 1 caractérisant la partie instationnaire de l’écoulement induite par les mouvements de parois s’écrit dans le domaine de Laplace :
éq. 4.3.2-9
Trois conditions aux limites d’entrée-sortie permettent de résoudre ce système. La première de ces conditions est obtenue en supposant que l’écoulement est suffisamment régulier en amont de l’espace annulaire, pour que la composante tangentielle de la vitesse d’entrée puisse être négligée :
en
éq. 4.3.2-10
Les deux autres sont obtenues en appliquant l’équation de conservation de l’énergie cinétique, sous sa forme quasi-stationnaire, entre l’infini amont et l’entrée de l’espace annulaire, puis entre la sortie de l’espace annulaire et l’infini aval. On obtient alors respectivement, à l’ordre des perturbations :
éq. 4.3.2-11
Dans ces expressions,
et
représentent les parties stationnaires des coefficients de pertes de charge singulières d’entrée et de sortie. Ils prennent en compte la dissipation d’énergie induite, lorsque les parois sont fixes, par d’éventuelles brusques évolutions de la géométrie à l’entrée ou la sortie de l’espace annulaire. Dans la plupart des cas, ces coefficients peuvent être estimés simplement à l’aide de données de la littérature (Idel’cik par exemple). Lorsque la configuration géométrique d’entrée ou de sortie est très particulière, ces coefficients peuvent également être déterminés à l’aide d’un code de mécanique des fluides bidimensionnel adapté à l’étude des problèmes à parois fixes, de type N3S.
et
sont les parties instationnaires des coefficients de pertes de charge singulières. Ces coefficients prennent en compte les perturbations des lignes de décollement induites par les mouvements de structure. Ils peuvent être modélisés grâce à une approche quasi-stationnaire de même nature que celle introduite pour l’estimation du coefficient de frottement de paroi.
Le système [:ref:` éq. 4.3.2-9 < éq. 4.3.2-9 >`] est résolu analytiquement, à l’aide des conditions limites [:ref:` éq. 4.3.2-10 < éq. 4.3.2-10 >`] et [:ref:` éq. 4.3.2-11 < éq. 4.3.2-11 >`], en explicitant les fonctions
caractérisant le second membre.
Les perturbations
et
définissant le mouvement des parois, les parties perturbées du jeu annulaire et du rayon moyen sont alors définies, dans le domaine de Laplace, par :
éq. 4.3.2-12
éq. 4.3.2-13
Décomposition sur base modale#
Soit
le nombre de modes vibratoires de la structure dans la bande de fréquence étudiée. La décomposition sur base modale du mouvement des parois s’exprime de la manière suivante :
éq. 4.3.3-1
éq. 4.3.3-2
où |
et 
représentent les ordres de coque du i èmemode pour les mouvements respectifs des coques interne et externe, 
et 
permettent de caractériser les plans principaux de ces modes, 
et 
sont déduits des déformées de poutre des structures interne et externe associées au mode considéré, |
et |
représente le déplacement généralisé. |
Remarque :
Les fonctions
et
sont représentées, dans le cadre de la résolution analytique, sous forme de combinaisons linéaires de sinus, cosinus, sinus hyperbolique et cosinus hyperbolique :
éq. 4.3.3-3
éq. 4.3.3-4
avec
et
nombres d’onde du i èmemode pour les mouvements des coques interne et externe respectivement.
Les solutions du problème fluide
sont recherchées sous la forme de décompositions sur base modale déduites de celles de
explicitées par les relations [:ref:` éq. 4.3.3-1 < éq. 4.3.3-1 >`] et [éq4.3.3-2]. On obtient ainsi, dans le domaine de Laplace :
éq. 4.3.3-5
éq. 4.3.3-6
éq. 4.3.3-7
Expression des termes de la matrice de transfert des forces fluide‑élastiques#
L’effort fluide-élastique surfacique,
, est la résultante du champ de pression et des contraintes visqueuses et turbulentes exercées par l’écoulement sur les parois de la structure en mouvement.
éq. 4.3.4-1
L’effort fluide-élastique généralisé associé au \(i\) ème mode vibratoire de la structure,
, s’écrit ainsi :
éq. 4.3.4-2
Où
désigne la surface des parois de la structure mouillées par l’écoulement, et le vecteur
représente le i ème vecteur déformée modale dans cette expression. La représentation du champ de vitesses et de pression et la représentation sous forme d’une loi de paroi des contraintes visqueuses et turbulentes exercées sur la structure en mouvement permettent d’exprimer l’effort fluide-élastique généralisé
de la façon suivante :
éq. 4.3.4-3
avec
et
désignent respectivement les contributions des coques intérieure et extérieure. Ces contributions sont définies par :
éq. 4.3.4-4
et
éq. 4.3.4-5
Résolution du problème modal sous écoulement#
Comme on l’a expliqué au paragraphe [§ 4.2], on résout préalablement le problème modal en eau au repos, afin de prendre en compte le couplage inertiel entre modes. On estime ainsi la matrice de masse ajoutée par le fluide, en calculant
pour une vitesse moyenne de l’écoulement nulle. Les caractéristiques modales du système sous écoulement sont ensuite obtenues en perturbant les caractéristiques en eau au repos. On ne tient plus compte que de l’amortissement et de la raideur ajoutés : les termes de masse ajoutée précédemment calculés sont retranchés de la matrice
. Le couplage entre modes est alors négligé ; en conséquence, les déformées modales demeurent inchangées par rapport à celles en eau au repos. Les seuls paramètres perturbés par la mise en la fréquence et l’amortissement modal réduit. Ces paramètres sont calculés en résolvant
équations non linéaires mode par mode, par mise en œuvre d’une méthode de type Broyden.
Écoulement axial (exemple : assemblages combustibles)#
L’intégration de ce modèle d’excitation fluide-élastique dans Code_Aster a été abordée dans la note de spécifications [bib. 2]. La note de principe du modèle MEFISTEAU [bib. 6] constitue la documentation théorique de référence.
Description de la configuration étudiée#
On considère un faisceau de
cylindres circulaires mobiles en flexion et soumis à un écoulement incompressible de fluide visqueux, limité par une enceinte rigide cylindrique de section circulaire ou rectangulaire [Figure 5.1-a].
Figure 5.1-a: Faisceau sous écoulement axial
Les cylindres sont tous parallèles, dirigés suivant l’axe de l’enceinte. Ils ont une longueur commune, notée
. Pour simplifier les notations, on suppose par la suite que
est l’axe directeur. L’écoulement stationnaire est axial et supposé uniforme dans chaque section. La masse volumique du fluide pouvant être variable suivant l’axe
(gradients thermiques), la vitesse de l’écoulement stationnaire dépend aussi de la variable
.
Étapes de calcul#
La première étape concerne la détermination de la base modale en air du faisceau. Cette opération est réalisée par l’opérateur CALC_MODES. Cette étape est indispensable car les forces fluide-élastiques sont projetées sur cette base.
La deuxième étape concerne la prise en compte des forces fluide-élastiques par l’opérateur CALC_FLUI_STRU. Cette étape se décompose en 7 sous-tâches :
Pré-traitements#
1°/ |
Au moyen de la topologie du maillage, déduction des coordonnées des centres des cylindres du faisceau puis vérification de la bonne disposition des cylindres les uns par rapport aux autres (on vérifie notamment qu’il n’y a pas chevauchement entre deux cylindres) et par rapport à l’enceinte rigide. |
2°/ |
Détermination de la longueur d’excitation du fluide, commune à tous les cylindres, ainsi que d’une discrétisation associée suivant l’axe directeur. |
3°/ |
Constitution des tableaux donnant les déformées modales en air de chaque cylindre du faisceau, pour chacun des modes pris en compte pour le couplage fluide-structure. On interpole pour cela les déformées aux points de la discrétisation déterminée auparavant. |
Résolution du problème modal sous écoulement#
4°/ |
Résolution du problème fluide perturbé. La détermination du potentiel des vitesses perturbées nécessite l’inversion de systèmes linéaires d’ordres élevés appelant la mise en œuvre de la méthode de Crout. |
Pour chaque vitesse d’écoulement
5°/ |
Calcul des matrices de masse, d’amortissement et de raideur ajoutés donnant la matrice de transfert des forces fluide-élastiques dans la base modale en air :
|
6°/ |
Résolution du problème modal sous écoulement ; on résout le problème complet aux vecteurs et aux propres |
On ne néglige pas les termes extra-diagonaux de
. Après reformulation, la résolution est effectuée à l’aide de l’algorithme QR : obtention des masses, fréquences et amortissements modaux réduits sous écoulement
, déformées modales complexes
exprimées dans la base en air ; de ces dernières, on ne retient que la partie réelle après minimisation de la partie imaginaire (calcul d’un critère sur la partie imaginaire).
7°/ |
Restitution des déformées sous écoulement dans la base physique.
|
Fin de boucle sur les vitesses d’écoulement
Remarques :
La connaissance des coordonnées des centres des cylindres (pré-traitement 1°) est nécessaire à la résolution du problème fluide perturbé (sous-tâche 4°). Cette résolution conduit à l’estimation des termes de la matrice de transfert des forces fluide-élastiques (sous‑tâche 5°), qui font intervenir les perturbations de pression et de vitesse.
La détermination d’une longueur commune d’excitation et la création d’une discrétisation associée (pré‑traitement 2°) permettent de définir un domaine d’intégration sur les structures pour la projection des forces fluide-élastiques sur la base modale. L’interpolation des déformées modales aux mêmes points est donc requise (pré-traitement 3°).
Le comportement dynamique du faisceau sous écoulement peut également être étudié à l’aide d’une représentation simplifiée du faisceau (avec des tubes équivalents). Les étapes de calcul pour la prise en compte du couplage fluide-structure sont alors identiques à celles décrites précédemment, les seules différences apparaissant dans les pré-traitements. Cette deuxième approche est décrite plus précisément dans la note [bib. 2]. Dans l’étape 1° des pré‑traitements, les coordonnées des centres des cylindres du faisceau sont alors spécifiées directement par l’utilisateur, qui établit également la correspondance entre les cylindres du faisceau et les poutres de la représentation simplifiée donnée par le maillage. Dans l’étape 3° des pré-traitements, cette correspondance permet d’affecter aux cylindres du faisceau, aux points de discrétisation déterminés dans l’étape 2°, les déformées modales des poutres de la représentation simplifiée.
Résolution du problème fluide instationnaire#
Hypothèse simplificatrice#
H1 |
Le champ des vitesses fluides instationnaires est déterminé analytiquement en supposant que l’écoulement perturbé
|
éq. 5.3.1-1
Un tel champ de vitesses admet un glissement sur les parois des cylindres qui permettra de calculer la contrainte visqueuse par une loi de frottement.
H2 |
Le mouvement des cylindres ne génère des perturbations de vitesse
|
H3 |
Le champ de pression est décomposé en parties stationnaire et perturbée selon
|
éq. 5.3.1-2
où |
désigne le diamètre hydraulique du faisceau, 
désigne le coefficient de frottement local pour la vitesse stationnaire 
. Il dépend du nombre de Reynolds, calculé à l’aide de la vitesse stationnaire 
, du diamètre hydraulique du faisceau et de la rugosité des parois. Ce coefficient est déduit de la loi de Nikuradzé caractérisant les écoulements en conduite ; |
désigne le champ de pesanteur. Son action sur le champ de pression stationnaire dépend de l’inclinaison du faisceau 
. |
Détermination du potentiel des vitesses perturbées#
On recherche une solution analytique pour
sous la forme d’une superposition de singularités élémentaires qui s’écrivent :
éq. 5.3.2-1
au centre de chaque cylindre \(k\) et :
éq. 5.3.2-2
au centre de l’enceinte rigide quand celle-ci est circulaire où:
|
désigne l’ordre de troncature de la série de Laurent
|
|
désignent les coordonnées polaires dans un plan perpendiculaire à l’axe
|
|
désignent les coordonnées polaires dans un plan perpendiculaire à l’axe
|
Les coefficients
et
des expressions [:ref:` éq. 5.3.2-1 < éq. 5.3.2-1 >`] et [:ref:` éq. 5.3.2‑2 < éq. 5.3.2‑2 >`] sont déterminés en appliquant la condition aux limites de non pénétration :
sur le contour de chaque cylindre mobile
, cette condition s’écrit :
où |
et 
désignent les composantes du déplacement de la fibre neutre du cylindre k à l’abscisse 
dans le repère 
, |
désignent les coordonnées polaires dans le repère 
dont l’origine est prise au centre du cylindre 
, |
|
désigne le rayon du cylindre 
, 
. |
sur le contour d’une enceinte rigide circulaire, elle s’écrit :
où
désigne le rayon de l’enceinte.
Dans le cas d’une enceinte rigide rectangulaire, cette condition est prise en compte par une méthode dérivée de la méthode des « images » [bib. 6] ; le problème fluide confiné par l’enceinte rectangulaire est rendu équivalent au problème en milieu infini en créant des images des cylindres mobiles du faisceau par rapport aux côtés de l’enceinte. Cette méthode conduit à introduire de nouvelles singularités de la forme [:ref:` éq. 5.3.2-1 < éq. 5.3.2-1 >`], placées au centre des cylindres «images», dans l’expression de
. Elle n’ajoute cependant pas d’inconnue au problème puisque les coefficients pour ces nouvelles singularités sont dérivés de ceux des cylindres mobiles du faisceau par le jeu des images.
Finalement, le potentiel des vitesses perturbées s’écrit :
éq. 5.3.2-3
Où
désigne le nombre de cylindres mobiles du faisceau. Les fonctions
sont des combinaisons linéaires de
et
dont les coefficients sont déterminés par les conditions aux limites précédentes. Cela nécessite la résolution de systèmes linéaires d’ordres élevés et à matrices pleines. Les inversions sont réalisées en mettant en œuvre la méthode de Crout.
Modélisation des forces fluides#
On retient d’abord les forces dues aux perturbations de pression
, reliées au potentiel des vitesses perturbées par :
éq. 5.3.3-1
La résultante du champ de pression perturbée autour de chaque cylindre mobile est une force linéique
agissant suivant
. Cette force dépend linéairement des
et
, engendrant ainsi des termes de masse, d’amortissement et de raideur ajoutés.
On prend ensuite en compte les forces liées à la viscosité du fluide.
Dans une approche quasi-statique, on considère l’action du champ de vitesse fluide
autour d’un cylindre à l’instant
: dans le repère lié au cylindre, l’écoulement, de vitesse
à l’ordre 0, présente une incidence par rapport au cylindre qui est fonction des perturbations de vitesse et du mouvement du cylindre lui‑même. Il en résulte une force de traînée et une force de portance. On montre que les composantes suivant
de la force linéique résultante
s’écrivent, pour le cylindre
:
éq. 5.3.3-2
éq. 5.3.3-3
où |
désigne la pente à incidence nulle du coefficient de portance autour d’un cylindre très faiblement incliné ( 
= 0,08). |
et 
désignent les moyennes des perturbations de vitesse suivant les axes 
et 
autour des cylindres, qui dépendent linéairement des 
et 
(cf. [:ref:` éq. 5.3.2-3 < éq. 5.3.2-3 >`]). |
Ces forces engendrent des termes d’amortissement et de raideur ajoutés.
On prend enfin en compte l’action du champ de pression stationnaire sur les structures mobiles déformées. On montre que la force linéique résultante
sur le cylindre
a pour composantes, à l’ordre 1 :
éq. 5.3.3-4
éq. 5.3.3-5
Ces forces n’engendrent que des termes de raideur ajoutée et aucun couplage entre cylindres.
Les expressions [:ref:` éq. 5.3.3-1 < éq. 5.3.3-1 >`], [:ref:` éq. 5.3.3-2 < éq. 5.3.3-2 >`] et [:ref:` éq. 5.3.3-3 < éq. 5.3.3-3 >`] mettent en évidence la nécessité de résoudre le problème fluide perturbé avant d’estimer les forces fluide-élastiques.
Expression des termes de la matrice de transfert des forces fluide‑élastiques#
Récapitulatif des forces linéiques
Pour chaque cylindre
, les forces fluide-élastiques s’écrivent suivant
et
:
éq. 5.3.4-1
et sont combinaisons linéaires des :
Décomposition du mouvement sur base modale
Le mouvement du faisceau de cylindres est décomposé suivant
modes de vibration en air. On note
(
et
) les déformées suivant
du cylindre
correspondant au j ème mode du faisceau. Les composantes du déplacement de la fibre neutre du cylindre
à l’abscisse
peuvent alors s’écrire:
éq. 5.3.4-2
éq. 5.3.4-3
où
est le vecteur des déplacements généralisés.
Projection des forces sur base modale
On note
la projection des forces fluide-élastiques suivant le i ème mode du faisceau.
éq. 5.3.4-4
est une combinaison linéaire des
On note
le vecteur des forces fluide-élastiques modales :
qui s’écrit :
éq. 5.3.4-5
Où |
désigne la matrice des termes de masse ajoutée par le fluide, |
désigne la matrice des termes d’amortissement ajouté par le fluide, |
|
désigne la matrice des termes de raideur ajoutée par le fluide. |
Ces matrices sont carrées réelles d’ordre
et leurs termes sont indépendants du mouvement des structures. La matrice
est symétrique ; les matrices
et
ne le sont pas nécessairement.
La projection des équations du mouvement sur base modale fournit :
éq. 5.3.4-6
où |
désignent les matrices de masses, d’amortissements et de raideurs de structure en air ; ces matrices sont d’ordre 
et diagonales. |
Dans le domaine de Laplace, la relation [:ref:` éq. 5.3.4-6 < éq. 5.3.4-6 >`] devient :
éq. 5.3.4-7
On introduit ensuite la matrice de transfert des forces fluide-élastiques
définie par :
éq. 5.3.4-8
Et l’on retrouve la relation [:ref:` éq. 1.2-1 < éq. 1.2-1 >`] du paragraphe [§ 1.2] :
Résolution du problème modal sous écoulement#
Le problème modal sous écoulement est formulé par la relation [:ref:` éq. 5.3.4-7 < éq. 5.3.4-7 >`] du paragraphe précédent.
Ce problème est résolu après réécriture sous la forme d’un problème standard aux vecteurs et aux valeurs propres du type
.
La nouvelle formulation est la suivante :
éq. 5.4- 1
Remarques :
On double la dimension du problème par rapport à celle du problème initial.
Les propriétés des matrices
et
permettent l’inversion.
La résolution de ce problème se fait au moyen de l’algorithme QR. Les modules mis en œuvre par l’opérateur CALC_FLUI_STRU sont les mêmes que ceux utilisés par CALC_MODES.
Le problème aux éléments propres que l’on résout est un problème complexe. On obtient donc un nombre pair de valeurs propres complexes conjuguées deux à deux. On ne conserve que celles dont la partie imaginaire est positive ou nulle.
Les vecteurs propres sont complexes, définis à une constante complexe multiplicative près. Comme l’on ne prend en compte que des modes réels, il s’agit en premier lieu de déterminer, pour chaque vecteur propre, la constante qui minimise la partie imaginaire du vecteur par rapport à sa partie réelle, au sens de la norme euclidienne. Les vecteurs propres sont alors redéfinis par rapport à cette norme. Compte tenu de la normalisation utilisée, il est alors possible de ne conserver dans le concept mode_meca que la partie réelle des vecteurs propres. On restitue cependant, dans le fichier MESSAGE, des indicateurs sur les rapports entre partie imaginaire et partie réelle des vecteurs propres ainsi normés, afin que l’utilisateur puisse juger du biais introduit par la non prise en compte de la partie imaginaire des vecteurs normés.
Prise en compte de la présence des grilles du faisceau de tubes#
La modélisation décrite précédemment, des forces induites par un écoulement axial sur un faisceau de cylindres, ne prend pas en compte la présence des grilles du faisceau (par exemple, les grilles de mélange et de maintien des assemblages combustibles). Une comparaison de ce modèle avec des essais effectués sur la maquette CHAISE (dans la configuration d’un faisceau de neuf tubes flexibles comportant une grille) est présentée dans une note de synthèse [bib. 8]: on constate que le couplage fluide-élastique entre la grille et l’écoulement axial n’est pas négligeable et qu’il engendre un accroissement de l’amortissement modal réduit des tubes. L’objet de ce paragraphe est la description des effets supplémentaires dus aux grilles et de leur prise en compte dans le modèle MEFISTEAU.
Description de la configuration des grilles#
On restreint ici l’étude à deux types de grilles :
les grilles de maintien qui se situent aux extrémités du faisceau,
les grilles de mélange qui sont réparties entre les grilles de maintien.
Les grilles sont toutes positionnées perpendiculairement au faisceau de cylindres et se présentent sous la forme d’un réseau prismatique à base carrée de côté
et de hauteur
(suivant l’axe \(x\) des cylindres). Les grilles d’un même type sont caractérisées par des dimensions identiques.
Étapes de calcul supplémentaires#
La première étape supplémentaire concerne la spécification du type de configuration des grilles par l’opérateur DEFI_FLUI_STRU , puis la vérification de la bonne disposition des grilles les unes par rapport aux autres, et par rapport aux extrémités du faisceau.
La deuxième étape concerne la résolution du problème modal sous écoulement. Dans la boucle sur les vitesses d’écoulement, la matrice de transfert des forces fluide-élastiques dans la base modale en air est complétée par le calcul d’une matrice d’amortissement ajouté et d’une matrice de raideur ajoutée, liées aux grilles.
Modélisation des forces fluides exercées sur les grilles#
Calcul du saut de pression
Tout d’abord, la présence de grilles perturbe le champ de pression stationnaire
; on considère chaque grille comme une singularité entraînant un saut de pression, dont l’expression se met sous la forme:
éq. 5.5.3-1
où |
désigne le coefficient de perte de charge due à la grille, |
désigne la vitesse stationnaire de l’écoulement au niveau de la grille, |
|
désigne la masse volumique de l’écoulement au niveau de la grille, |
|
désigne la position axiale du milieu de la grille le long du faisceau. |
La masse volumique
est calculée par interpolation linéaire du profil de masse volumique
de l’écoulement en l’absence de grille. La vitesse stationnaire
est calculée en application de la conservation du débit massique, qui se traduit par l’équation suivante :
Calcul des forces fluides ponctuelles exercées sur chaque grille
Suivant la même approche quasi-statique que celle effectuée au paragraphe [§5.3.3], on montre que l’action du champ de vitesse fluide
autour d’une grille implique une force de traînée et une force de portance, en fonction de l’incidence de l’écoulement par rapport à la grille. Les composantes y et z de la force ponctuelle résultante f
s’écrivent donc, pour chaque cellule élémentaire
d’une
où |
et 
désignent respectivement le profil de masse volumique et de vitesse stationnaire de l’écoulement en pied de faisceau, |
désigne la section fluide du faisceau en l’absence de grille, |
|
désigne la section fluide du faisceau au niveau de la grille: 
avec 
la section solide de la grille. |
On en déduit l’expression:
Le coefficient de perte de charge
est calculé à partir de l’expression de la force hydrodynamique totale qui s’applique sur la grille, et nous obtenons:
éq. 5.5.3-2
Le 1er terme (en
) provient de l’effort de traînée;
est le coefficient de traînée de la grille. Le 2 ème terme (en
) est un terme correcteur de l’effort de frottement appliqué sur le faisceau seul à l’altitude de la grille (
est le périmètre mouillé du faisceau en l’absence de grille).
En introduisant l’expression [:ref:` éq. 5.5.3-2 < éq. 5.5.3-2 >`] dans la relation [:ref:` éq. 5.5.3-1 < éq. 5.5.3-1 >`], on obtient donc l’expression du saut de pression
pour chaque grille d’altitude
. Ce saut de pression est pris en compte au niveau du calcul du champ de pression stationnaire
, de la manière suivante:
grille:
où |
désigne la pente à incidence nulle du coefficient de portance autour d’une grille très faiblement inclinée. |
désigne la section solide de la cellule élémentaire 
de la grille (qui en comprend 
). |
Ces forces vont donc engendrer des termes supplémentaires d’amortissement et de raideur ajoutés, que l’on obtient après décomposition modale du mouvement et projection de ces forces sur la base modale.
Prise en compte de l’amortissement en fluide au repos#
Jusqu’à présent, l’amortissement apporté à un faisceau de tubes par la présence d’un fluide au repos n’a pas été pris en compte dans la modélisation. On propose donc ici un modèle d’amortissement en fluide au repos, dont l’annexe 1 de la note de synthèse des essais CHAISE [bib8] constitue la documentation de référence.
Modélisation de la force fluide au repos exercé sur un faisceau de tubes#
La méthode de calcul de l’amortissement en fluide au repos qui est mise en œuvre ici, est une généralisation de la méthode de CHEN [bib. 9].
Il s’agit de calculer la force résultante sur chaque tube des contraintes dues au cisaillement dans la couche limite. C’est un problème non linéaire car le coefficient d’amortissement fluide dépend de la fréquence On introduit donc les simplifications suivantes :
le problème est écrit à l’aide des fréquences en eau au repos calculées sans prendre en compte l’amortissement fluide,
on néglige le couplage entre modes.
La force linéique
s’exerçant sur le tube
soumis à un mouvement harmonique du faisceau suivant le mode
à la fréquence
est donnée par la relationsuivante:
éq. 5.6.1-1
où |
|
désigne la vitesse de glissement entre le tube
|
|
désigne le coefficient de traînée d’un cylindre de rayon
|
La relation obtenue en remplaçant [:ref:` éq. 5.6.1-2 < éq. 5.6.1-2 >`] et [:ref:` éq. 5.6.1-3 < éq. 5.6.1-3 >`] dans l’équation [:ref:` 5.6.1-1 < 5.6.1-1 >`] est linéarisée par un développement en série de Fourier (du terme
) dont on ne retient que le premier terme; il vient:
Projection sur base modale
Par projection sur base modale et en négligeant le couplage entre modes, on obtient la force généralisée s’exerçant sur le faisceau de tube suivant le mode
:
est donc proportionnelle à
et le vecteur de force modale associé
se met sous la forme :
où
désigne la matrice d’amortissement ajouté par le fluide au repos.
Bibliographie#
GAY, T. FRIOU : Résorption du logiciel FLUSTRU dans ASTER. HT32/93/002/B
PEROTIN, M. LAINET : Intégration de différents modèles d’excitations fluide-élastiques dans le Code_Aster ® : spécifications. HT-32/96/014/A
GRANGER, N. GAY : Logiciel FLUSTRU Version 3. Note de principe. HT32/93/013/B
GRANGER : Compléments théoriques pour l’interprétation des essais GRAPPE2 sous écoulement. HT32/92/025/A
PEROTIN : Note de principe du modèle MOCCA_COQUE. HT32/95/021/A
BEAUD : Note de principe du modèle MEFISTEAU. HT-32/96/005/A
GRANGER : « A Global Model For Flow-Induced Vibration Of Tube Bundles In Cross-Flow »ASME Journal of Pressure Vessel Technology, 1991, Vol.113, pp.446-458.
J-L. SAGE, F. BEAUD, P. MANDROU : Synthèse des essais CHAISE en écoulement axial et interprétation avec le modèle MEFISTEAU. HT-32/99/003/A
R.D. BLEVINS : « Flow-Induced Vibrations », Krieger Publishing Company, 1994, pp308-310.
Connors H.J., «Fluid-elastic vibration of heat exchanger tube arrays», Journal of Mechanical Design, vol. 100, april 1978, [Galaxie 79H513944]
«Stabilité vibratoire des tubes de générateur de vapeur : impact sur les marges de la prise en compte de la fonctionnalité ‘Zones multiples d’excitation’ dans le logiciel GEVIBUS », ADOBES A., DUBRUQUE D., note EDF-R&D HI-86-02-009-A
Description des versions du document#
Version Aster |
Auteur(s) Organisme(s) |
Description des modifications |
8.4 |
T.Kestens EDF-R&D/MFEE |
Texte initial |
9.2 |
A.Adobes, E. Longatte EDF-R&D/MFEE |
Ajout de précisions concernant le mode de calcul de l’instabilité fluide-élastique par la méthode de Connors et par la méthode des corrélations fluide-élastiques. |