Accélération De Convergence Des Algorithmes De Résolution .

2y ago
26 Views
5 Downloads
1.61 MB
36 Pages
Last View : 1m ago
Last Download : 3m ago
Upload by : Macey Ridenour
Transcription

RAPPORT de D.E.A.Centre de Mathématique et d’Informatique de MarseilleStage du 2 avril au 30 août 2002Soutenance : 9 Septembre 2002Accélération de convergence des algorithmes derésolution des problèmes d’aérodynamique stationnaireM ICKAËL FONTAINERéférence : WN/CFD/02/85Encadrement :S. Champagneux et G. Chevalier, CERFACS-CFD, 42 av. G. Coriolis, 31500 Toulouse

RemerciementsJe tiens tout d’abord à exprimer ma reconnaissance à Thierry POINSOT, chef de l’équipe CFD, pourm’avoir accepter au sein de son équipe.Je remercie Guilhem CHEVALIER de m’avoir proposé ce stage, Steeve CHAMPAGNIEUX de m’avoirencadré et conseillé dans toutes les phases du stage, Marc MONTAGNAC pour son aide précieuse sur elsA.Un merci aussi à ceux qui oeuvrent pour que les machines répondent toujours présent : le groupe CSG,avec en particulier Isabelle et Gérard .Je tiens également à remercier les seniors et ingénieurs , en particulier Benedicte et Jean Christophe ,pour nos nombreuses discussions.Merci enfin à la documentaliste, Jenny, à la secrétaire de l’équipe, Marie, ainsi qu’à tous les stagiaireset thésards qui contribuent largement à la bonne ambiance qui règne au sein du Cerfacs, grâce à eux j’aipassé cinq mois très agréables.1

RésuméDiminuer les temps de calcul sur un code de simulation numérique, tout en conservant une bonne précision des résultats est un objectif constant dans le milieu de la mécanique des fluides numérique. L’objectifde ce travail est de trouver des moyens d’atteindre ce but pour le code de calcul elsA.Des tests ont été effectués pour evaluer l’influence de paramètres géométriques et topologiques du maillagesur les temps de calcul sans donner de résultats probants.L’utilisation d’un procédé s’inspirant des méthodes multizones, a par contre, donné de meilleures résultatsen ce qui concerne la diminution des temps de calculs.AbstractOne of the main purpose in cfd is to decrease the calculation time, whyle keeping a great accuracy ofthe results. The aim of this work is to find methods in order to reach that objective with the elsA software.Some tests have been done to evaluate the influence of different block orientation and numerotation of thegrid on the calculation time, but no real effects was noticed.By using a kind of zonal approach, an interesting decreasing of calculation time was obtained.2

Table des matièresRemerciements1Introduction51 Analyse du problème et des algorithmes de résolution1.1 Modèle Physique . . . . . . . . . . . . . . . . . .1.2 Modélisation numérique . . . . . . . . . . . . . .1.2.1 Formulation Volumes Finis . . . . . . . . .1.2.2 Discrétisation spatiale . . . . . . . . . . .1.2.3 Intégration en temps . . . . . . . . . . . .1.2.4 Résolution du système d’équations . . . . .778889102 Stratégies liées aux algorithmes de résolution2.1 Préliminaires . . . . . . . . . . . . . . .2.2 Tests . . . . . . . . . . . . . . . . . . . .2.2.1 Orientation des blocs . . . . . . .2.2.2 Numérotation des blocs . . . . . .13131414153 Approche zonale3.1 Principe . . . . . . . . . . . . . . . . .3.2 Adaptation au code elsA . . . . . . . .3.2.1 Structure du code elsA . . . . .3.2.2 Modifications apportées au code3.2.3 Résultats . . . . . . . . . . . .3.3 Premiers tests sur un cas 3D . . . . . .3.3.1 Résultats . . . . . . . . . . . esMaillage 30 blocs du profil 2D RAE2822 . . . . .Orientation des blocs . . . . . . . . . . . . . . . .Numérotation des blocs . . . . . . . . . . . . . . .Exemple de résultats sur la configuration init . . .29293031333

Table des figures2.12.22.32.42.5Profil du RAE2822 . . . . . . . . . . . . . . . . . . . . . .Découpage du domaine de calcul en blocs . . . . . . . . . .Visualisation des blocs du domaines de calcul près du profilOrientation initiale des blocs ( init ) . . . . . . . . . . . . .Numérotation initiale des blocs . . . . . . . . . . . . . . . .13131415163.13.23.33.4Schéma de principe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Structure en “couches” du code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .comparaison des Cp obtenus sans et avec utilisation du principe zonal sur le cas test init . .comparaison des composantes suivant x, y et z du vecteur frottement a la paroi obtenus avecet sans utilisation du principe zonal . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Maillage de la zone de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Maillage sur l’aile . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Coupe du domaine à y constant . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Maillage de la zone de calcul . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Visualisation du maillage près du profil . . . . . . . . . . . . . . . . . . . . . . . . . . .orient1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .orient2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .orient3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .orient4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .orient5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .num1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .num2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .num3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .num4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .num5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .num6 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Variation des résidus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Variations des efforts . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .champs de vitesse autour du profil ( nombre de Mach ) . . . . . . . . . . . . . . . . . . .champs d’énergie turbulente autour du profil avec visualisation des blocs . . . . . . . . . 30303031313132323233333434

IntroductionL’ étude de phénomènes aérodynamiques autour de composantes ( aile, fuselage, nacelles . ) d’avionsnécessite l’utilisation de codes de simulations numériques. Ces derniers sont nombreux et adaptés souvent aux types d’études réalisées, mais pour plus de facilité d’utilisation, trois organismes, Airbus France,l’Onera et le Cerfacs, se sont regroupés afin de mettre au point un code généraliste commun en CFD 1 , celogiciel, elsA 2 est actuellement en cours de développement dans les laboratoires, et une première versionsera mise en production au niveau industriel en fin d’année 2002.Mon stage s’incrit donc dans cette phase de développement du code, notamment dans l’ optimisation de lavitesse de résolution des calculs. En effet, comme avec tous les codes de CFD, les calculs peuvent devenirvite très longs en instationnaire ou sur des configurations tridimensionnelle par exemple , par conséquent ilest nécessaire de trouver des moyens pour réduire ces temps d’exécution, ce qui a été le but du stage.Mais afin de décrire les procédés utilisés pour accélérer la résolution des problèmes d’ aérodynamique, ilfaut dans un premier temps analyser les algorithmes mis en jeu dans le code afin d’en comprendre les mécanismes, puis déterminer les gains potentiels réalisables.Ces différents points seront développés dans la suite du rapport, mais avant tout une présentation du CERFACS doit être faite pour situer la cadre dans lequel s’est déroulé le stage.Le CERFACS ( Centre Européen de Recherche et de Formation Avancée en Calcul Scientifique ), fondé en1987, est un laboratoire de réputation internationale dans le calcul scientifique à haute performance. Uneéquipe de plus de 100 personnes, chercheurs, post-doctorants, doctorants, consultants et stagiaires représentant de nombreuses nations y travaillent.Ce groupe de personnes est divisé en quatre équipes représentant les axes de recherche du laboratoire,– Algorithmique et Calcul parallèle (développement d’algorithmes de résolution numérique optimiséspour des calculateurs à hautes performances).– Climatologie (Modélisation des interactions océan-atmosphère, modèles climatiques, assimilation dedonnées et prévisions climatiques ).– Electromagnétisme (Simulation de propagation d’ondes, application à la compatibilité électromagnétique).– Mécanique des fluides numérique (CFD)L’ équipe CFD regroupe deux domaines de recherche, l’ aérodynamique (Simulation d’écoulements autourd’un avion et phénomènes liés, par exemple : sillage, intéraction choc/couche limite, etc.) et la combustion(Combustion interne aux moteurs pour l’automobile et les turbines, études pour la réduction des émissionspolluantes etc . ).1 Computational2 ensembleFluid Dynamiclogiciel pour la simulation en Aérodynamique5

6

CHAPITRE 1. ANALYSE DU PROBLÈME ET DES ALGORITHMES DE RÉSOLUTIONChapitre 1Analyse du problème et des algorithmesde résolutionLe logiciel elsA est un outil de simulation numérique en Aérodynamique : il résout les équations deNavier-Stokes compressibles 3D stationnaires et instationnaires discrétisées via une approche volume finissur des maillages de type multi-blocs structurés. Le problème résolu est de type système d’EDP soumis àdes conditions initiales et aux limites.Dans une première partie, il sera question de la modélisation physique des problèmes étudiés, puis dans uneseconde, l’accent sera mis sur l’aspect numérique .1.1 Modèle PhysiqueLe cadre de l’étude sera le suivant : placé en régime continu, les fluides étudiés seront considérés commeconstitué d’une seule espèce et avec une variation de la masse volumique ρ suffisamment faible pour queles effets de gravité soient négligés.Il est admis que cette description relève des équations de Navier-Stokes compressibles complétées par deslois de comportement et des lois d’états , ces équations n’étant autre que des équations de bilan sur lesvariables conservatives : la masse volumique, ρ, la quantité de mouvement par unité de volume, ρU etl’énergie totale par unité de volume, ρE : ρbilan de massediv ρU0(1.1) t bilan de quantit é de mouvement ρE tbilan d énergie ρU tdiv ρU pUdiv ρE UUpI τ 0 τ U q 0 (1.2)(1.3)où U est le champ de vitesse exprimé dans un repère absolu, p, τ et q, respectivement la pression, letenseur des contraintes dues à la viscosité et le vecteur du flux de chaleur dû à la conductivité thermique.Pour fermer le problème, il faut préciser les lois de comportement et d’état du système.Le modèle du gaz parfait, à chaleur spécifique constante sera alors retenu pour le fluide, ce qui signifie quel’energie interne e et la pression p seront donnés par les relations,e Cv Tp rgaz ρT(1.4)où rgaz représente le rapport de la constante universelle des gaz parfaits à la masse molaire du gazconsidéré, et Cv , la chaleur spécifique à volume constant.Par ailleurs pour un fluide newtonien, le tenseur des contraintes est donné par la loi,τ λ divU I2µD(1.5)7

1.2. MODÉLISATION NUMÉRIQUEoù D est le tenseur de taux de déformation, et λ et µ les deux coefficients de viscosité du fluide, vérifiantl’hypothèse de Stokes, 3µ 2λ 0.Le vecteur du flux de chaleur est quant à lui donné par la loi de Fourier,qKT gradT(1.6)où KT désigne le coefficient de conductivité thermique.1.2 Modélisation numérique1.2.1 Formulation Volumes FinisLe code elsA utilise une formulation en volumes finis sur un maillage constitué de cellules hexaédriquesoù les inconnues sont calculées au centre de ces dernières ( approche “cell centered”). Sur ces cellules serontdonc déterminés des bilans de flux établis à partir des équations de Navier-Stokes ( 1.1, 1.2 et 1.3) écritessous la forme,ddt Fc W s WdΩΩt Fd W gradW Ω t ndΣ T W gradW dΩ Ωt(1.7)où W est le vecteur des variables conservatives, Fc, celui du flux convectif,Fd, le flux diffusif et T, leterme source volumique.W ρρUρE (1.8)1.2.2 Discrétisation spatiale V Ω T̃ La discrétisation spatiale de l’équation (1.7) sur une cellule élémentaire Ω t du domaine , donne, aprèsque les termes aient été moyennés 1 sur chaque Ω t de volume V Ω t : dV Ω W̃ F dt 6n ΣiΩi 1(1.9)Ω La formulation précédente est une formulation exacte, mais à présent le flux exact F n Σi va être approché, par un flux numérique F W Wi de même que W̃ et T̃, ce qui donne alors, avec l’hypothèse d’unmaillage fixe 2 ,dWΩdt F W W N 1V Ω1 RV Ω 6ΩΩiV Ω TΩΣii 1où RΩ est le résidu numérique de modélisation. Ω(1.10)Le flux F W Wi est constitué d’un flux convectif Fc et d’un flux diffusif Fd, ces flux peuvent être alorsapprochés par plusieurs schémas, centrés ou décentrés, mais pour la suite , il sera retenu une approche centrée de Jameson stabilisée par un flux de dissipation artificielle pour le flux convectif et des schémas centréspour la détermination du flux diffusif.Le terme source quant à lui peut comporter plusieurs composants, production, dissipation, bas Reynolds.Ces termes s’écrivent entre autres en fonction de gradU et de grad gradU qu’ il faudra discrétiser. gradUle sera à la manière décrite en 1.13 et pour grad gradU une extrapolation de gradU à l’ordre un serautilisée.1 par8!" !V1ex : W̃Ω tV tdit, V Ω t2 autrementΩtW t dΩΩ

CHAPITRE 1. ANALYSE DU PROBLÈME ET DES ALGORITHMES DE RÉSOLUTIONSchéma de JamesonPour ce schéma, le flux de Jameson, approximation du flux convectif, s’écrit sur les faces Σ l du volumede contrôle Ω,# % & ' 12 ( Fc # W % )FJam WΩ WΩl NΣl#Fc WΩlΩ% * &N ΣlDl(1.11)oú W et Wl sont respectivement les inconnues conservatives dans la cellule Ω et dans la cellule voisineΩl , adjacente à la face Σl , et où Dl désigne le flux de dissipation artificielle.Le schéma de Jameson s’écrit alors :' V # Ω1 % , - 21 ( Fc # W % ) Fc# W % * & N 6dWΩdt / /ΩDiΣlΩll 1 Dj Dk V # Ω% T .Ω(1.12) /où Di , D j et Dk représentent l’opérateur de dissipation artificielle, appliqué dans la cellule i jk, dans lesdirections i , j et k du maillage structuré. Cet opérateur comporte pour chaque direction, un terme dedissipation non linéaire du second ordre pour capturer les discontinuités de l’écoulement et un terme dedissipation linéaire du quatrième ordre qui avec celui du second ordre sert à assurer la stabilite de l’approximation. Ces opérateurs de dissipations comportent même des termes appelés senseurs ( combinaison d’unedifférence seconde de la pression statique et de la vitesse ) dont le rôle est de détecter les ondes de chocainsi que les discontinuités de contact.Détermination du flux diffusifLe flux diffusif comporte des termes de gradients de vitesse, de température, de quantités turbulentes.La discrétisation de ce terme fait donc intervenir deux étapes, une d’évaluation numérique des gradients,une autre de discrétisation du flux diffusif.L’évaluation des gradients se fait grâce à une formule de la moyenne sur la cellule Ω de centre M et duthéorème de Green-Ostrogradski, par exemple pour gradU,gradUM'' - V #1Ω%13 0 U4 2 5 ndΣ67 8 9 :;61V Ω# %0ΩgradUdΩΣii 112UΩ UΩ i(1.13)N Σila densité de flux elle, se déterminera pour les faces intérieures Σi , par la demi somme des densités deflux dans la cellule Ω et dans la voisine Ωi ,0Σi&Fd ndΣ' 12 # Fd )Ω%&FdΩi NΣi(1.14)1.2.3 Intégration en tempsEn notant RΩ , le résidu de modélisation ( cf. 1.10 ), l’équation différentielle ordinaire à résoudre devient,dWΩdt' V 1# Ω% RΩ(1.15)La méthode utilisée pour déterminer la solution de 1.15 est du type “Alternative Explicite/Implicite”,puisque l’intégration de cette équation se fait de manière explicite, le schéma le plus fréquemment utilisépour la discrétisation de cette expression dans le code elsA étant le schéma de Runge Kutta à 4 pas ( schémaexplicite, du second ordre en temps, de domaine de stabilité linéaire assez étendu ), alors que par la suite, larésolution du système se fait de manière implicite ( cf. 1.2.4 ).9

1.2. MODÉLISATION NUMÉRIQUE ?@où α114,α2@13,α3@12@W0.Wn Wi@BAWi@W0@W4.Wnet α4L@1C D CE DF W Gαi V tΩ RΩi 1i HJI 1 G 4Ki(1.16)1.Dans le cas d’un écoulement stationnaire, cas étudié par la suite ( Chap. 2 ), la méthode sera dite “pseudoinstationnaire” puisque qu’elle n’est pas nécessairement consistante avec une évolution instationnaire ayantun sens physique, le temps étant alors un paramètre itératif permettant de converger vers une solution stationnaire.Dans l’expression du résidu explicite, le flux convectif est recalculé à chaque pas Runge-Kutta alors que leflux diffusif ainsi que le terme de dissipation numérique ( cf. section 1.2.2 ) et le terme source sont figés àla première étape pour ne pas rendre le coût de l’intégration en temps trop prohibitif.Détermination du pas de tempsLa détermination du pas de temps nécessite la connaissance des valeurs propres de la matrice d’amplification du schéma numérique. Dans le cas pseudo-instationnaire une convergence rapide vers un étatstationnaire est recherchée, cela implique l’utilisation d’un pas de temps local ( différent d’une maille àl’autre ) calculé comme le minimum de deux pas de temps locaux associés aux opérateurs hyperboliques etparaboliques des équations de Navier-Stokes.Le critère de stabilité de l’opérateur hyperbolique 3 , donne la condition suivante sur le pas de temps ,@ tCCFLM U hM F a(1.17)avec CFL, le nombre de Courant-Friedrich-Lewy 4 introduit pour assurer la stabilité du schéma numérique , h une dimension caractéristique du volume de controle et a la vitesse du son.2Le critère de stabilité de l’opérateur parabolique donne en monodimensionnel la condition, t D h2κ .NLe pas de temps local sera donc le minimum des deux pas de temps, t@O G Pmin tC tD(1.18)1.2.4 Résolution du système d’équationsLa résolution du système se fait par une méthode d’inégration implicite, ainsi à partir de notre système(1.15), le problème peut se ramener à la résolution deLoù RnL AL AWntnV111@WntnRnL1(1.19)sera déterminé par linéarisation à l’ordre 1 en temps, au voisinage de R n sur le domaine Ω,L QnRΩ1 tC xRnΩW XXFS RΩ WΩ WF6 RΩ WΩ WΩiTRUi 1iV(1.20)J1monodimensionnel ce critère s’écritλles équation de Navier-Stokes, il n’est pas possible de déterminer ce nombre par une analyse théorique de la stabilité numérique, il est donc pris voisin d’un nombre calculé moyennant de nombreuses approximations dans les équations3 en4 pour10

CHAPITRE 1. ANALYSE DU PROBLÈME ET DES ALGORITHMES DE RÉSOLUTIONYLe système (1.15) peut alors se mettre sous la forme, hI J W tZ [\ ]R(1.21)Méthode ADILa technique des directions alternées ( ADI5 ) consiste à substituer à l’opérateur implicite, un opérateurfactorisé suivant les directions du maillage.En effet, une méthode implicite classique demande la résolution d’un grand nombre d’équations, M N àla fois et à chaque pas de temps, la méthode ADI quant à elle se base sur le principe d’être implicite dansune direction x ou y à un pas de temps, et implicite dans l’autre direction au pas de temps suivant. Au coursd’un cycle complet, il sera donc résolu M puis N

RAPPORT de D.E.A. Centre de Mathématique et d’Informatique de Marseille Stage du 2 avril au 30 août 2002 Soutenance : 9 Septembre 2002 Accélération de convergence des algorithmes de résolution des problèmes d’

Related Documents:

BA-XS-CP-880 G4 Acc. to IEC 62471/Acc. to IEC 62778/Acc. to IEC 61347-1/Acc. to EN 62031/Acc. to IEC 60598-1/Acc. to EN 55015/Acc. to EN 61547 BA-XS-CP-865 G4 Acc. to IEC 62471/Acc. to IEC 62778/Acc. to IEC 61347-1/Acc. to EN 62031/Acc. to IEC 60598-1/Acc. to EN 55015/Acc. to EN

Jul 17, 2003 · ACC 511 Financial Accounting 5110 ACC 512 Managerial Accounting Systems 5120 ACC 515 Federal Income Taxation 5150 ACC 518 Intro AIS & Databases 5180 ACC 521 Federal Income Tax II 5210 ACC 522 Auditing 5220 ACC 524 Government/Not-for-Profit Acct 5240 ACC 526 AIS: Audit & Control 5260 AC

ACC 111 Principles of Accounting I 4.00 II Lehman ACC 171 Principles of Accounting I ACC 111 Principles of Accounting I 4.00 II Lehman ACC 1000 100-Level Elective ACC 112 Principles of Accounting II 4.00 II Lehman ACC 272 Principles of Accounting II ACC 113 Princ of Intermdte Accountg 4.00 II Lehman ACC 3

55015/Acc. to EN 61547/Acc. to EN 61000-3-2/Acc. to EN 62384/Acc. to EN 62386/Acc. to IEC 62386-101:Ed2/Acc. to IEC 62386-102:Ed2/Acc. to IEC 62386-207:Ed1 Protection class II Type of protection IP20 Logistical data Commodity code 850440829000 Download Data File User instruction OPTOTRONIC L

ACC 102 Accounting Principles II Note: MCC ACC 101 ACC 102 UR ACC 201 ACC 201 Accounting Using QuickBooks ACC 110 Fundamentals of Accounting I Must complete ACC . ANT 205 Archaeology Field School ANT 216 Special Topics in Anthropology ART ART 110 Comics and Sequential Art Accepted as Intro Studio course

1443 Series Accelerometers Specifications Catalog Numbers 1443-ACC-GP series, 1443-ACC-VO series, 1443-ACC-IS series, 1443-ACC-AT series, 1443-ACC-LF-T, 1443-ACC-HF-T Summary of Changes . Portable Data Collector an d Permanently Installed Accelerometers Figure 3 - Portable Data Collector and Accelerometer

dept crs instructor title author isbn ed. new used acc 201 depusoir financial & managerial acc miller-nobles 9780134491554 6th 260.00 acc 202 tba financial & managerial acc miller-nobles 9780134491554 6th 260.00 acc 301 tba intermediate accounting i keiso 9781119503583 17th 180.00 acc 303 tba intermediate accounting iii kieso 9781119503583 17th 180.00

the ACC website - contact our provider relationship team for details of the Engagement and Performance Manager in your region. Provider Contact Centre 0800 222 070 providerhelp@acc.co.nz. Provider Registration 04 560 5211 registrations@acc.co.nz ACC Portfolio Manager or Advisor vrs@acc.co.nz .