Un autre calcul des coefficients de relaxation 

Dans ce paragraphe, nous présentons une variante du calcul des coefficients de relaxation. Il doit être noté que, en effet, le développement de Chapman-Enskog ne donne qu’une indication sur la stabilité de la méthode : elle n’assure que la stabilité du problème quasi-linéaire associé. Nous allons vérifier qu’une approche différente permet de diminuer la diffusion numérique. Par ailleurs, il doit être souligné que l’approche présentée dans l’article précédent est relativement complexe. En effet, si il est difficile d’assurer la positivité des valeurs propres d’une matrice d’ordre 3, il sera humainement impossible d’assurer la même condition pour une matrice d’ordre plus élevé. Or, comme nous allons le vérifier dans la deuxième partie de cette thèse, nous avons pour objectif final de traiter des système d’équations aux dérivées partielles de plus grande taille (par exemple le système à K constituants ou bien le système triphasique). Dans ce cadre, la méthode actuelle n’est plus envisageable et une alternative doit être trouvée. Le plan de ce paragraphe est le suivant. Nous commençons par analyser la structure du système de relaxation en coordonnées Lagrangiennes. Puis nous faisons des liens entre le calcul des coefficients de relaxation actuels et des systèmes de lois de conservations plus simples (le p-système et une équation de conservation scalaire). Nous présentons alors la nouvelle manière de calculer les coefficients de relaxation. Cette nouvelle stratégie est soutenue par une expérience numérique spécialement choisie. Reprenons le système d’équilibre gaz-liquide en coordonnées Eulériennes (2.16) avec S = 0 :    ∂tρ + ∂x(ρv) = 0, ∂t(ρv) + ∂x(ρv2 + P(u)) = 0, ∂t(ρY ) + ∂x(ρY v − σ(u)) = 0. (4.1) 4.1.1 Système relaxé en coordonnées Lagrangiennes Afin de bien mettre en évidence les véritables non-linéarités du système, plaçons-nous en coordonnées Lagrangiennes. En introduisant le co-volume τ = 1 ρ et la variable Lagrangienne y définie par on aboutit à    ∂tτ − ∂yv = 0, ∂tv + ∂yP = 0, ∂tY − ∂yσ = 0. (4.3) Par abus de notation, P et σ doivent être maintenant considérés comme des fonctions de (τ, v, Y ). On voit tout de suite que les difficultés théoriques sont concentrées dans P et σ, pression et impulsion apparentes. Ce sont par conséquent ces termes-là qu’il faut relaxer. 

Remarque 

 On considère ici le système (4.3). Introduisons une variable supplémentaire Π censée être égale à P et considèrons le système 4 × 4 .Le système est hyperbolique et ses valeurs propres sont (0, ±a, −σY ). La relaxation de l’équation de quantité de mouvement introduisant le coefficient a peut donc être vue comme un moyen de forcer l’hyperbolicité du système et d’assurer la positivité de la densité, comme nous allons le voir par la suite. Il apparaˆıt par ailleurs que la dernière équation de (4.4) est non-linéaire. La relaxation de cette équation peut donc être vue comme un moyen de forcer la propriété physique Y ∈ [0, 1] grˆace au coefficient b, comme nous allons le voir par la suite. 

Ancienne technique de calcul des coefficients de relaxation

 Dans ce paragraphe, nous présentons la technique initialement envisagée pour le calcul des coefficients de relaxation ainsi que les raisons qui nous ont fait abandonner cette technique. Dans cette technique, le système de relaxation était réduit par le développement de Chapman-Enskog à la forme quasi-linéaire : ∂tw + Γ(w)∂yw = λ −1 ∂y [∆(w)∂yw] , (4.9) avec w = (τ, v, Π, Y, Σ)T et Γ la Jacobienne du système d’équilibre en coordonnées Lagrangiennes.Nous voudrions ici mettre en lumière la structure très singulière de ces coefficients. Le coefficient −Pτ +P 2 v dans A est un terme intrinsèquement lié au sous-système (A) : il contient les dérivées de P par rapport aux variables d’équilibre régies par les équations d’évolution de ce soussystème, c’est-à-dire τ et v. Le coefficient σ 2 Y dans B est un terme intrinsèquement lié au sous-système (B) : il contient la dérivée de σ par rapport à la seule variable d’équilibre régie par une équation d’évolution dans ce sous-système, c’est-à-dire Y . Tous les autres termes sont en fait des termes de couplages entre sous-systèmes. Notons par ailleurs que l’étude (en annexe) du système de relaxation avec glissement nul produit un coefficient A = −pτ . D’autre part, si l’équation d’évolution en Y était une simple équation scalaire, découplée des autres variables, le coefficient sortant du calcul serait B = σ 2 Y . D’un tout autre point de vue, certains cas-tests exhibent des coefficients a et b tout à fait surestimés, lorsqu’on considère l’approche du rapport antérieur. Cela se traduit par une diffusion trop importante. Un cas-test de ce type est présenté à la fin de ce chapitre. C’est cet ensemble de remarques qui nous a conduit à formuler le développement de ChapmanEnskog suivant. 

Traitement aux limites

Nous présentons ici la technique utilisée pour calculer les points fictifs amont et aval. Le traitement des conditions aux limites est à la fois essentiel et extrêment délicat. C’est un problème peu connu d’un point de vue théorique (pour les systèmes d’équations aux dérivées partielle) et difficile d’un point de vue numérique dans le cadre des écoulements en conduite. Plus spécifiquement, comme toutes les informations qui entrent dans la conduite pétrolière entrent par l’amont, le traitement de la condition limite amont est très important. Autrement dit, il est inutile de construire un schéma numérique intérieur précis si la condition à la limite amont n’est pas traitée de façon adéquate. Or, il y a une difficulté à prendre en compte les conditions aux limites lorsqu’on considère un schéma semi-implicite. En effet, lorsqu’on considère un schéma explicite, une stratégie qui consiste à contrˆoler la solution du problème de Riemann sur l’arête limite est suffisante (c’était notre premier choix, voir [10]). En fait, si on utilise un schéma explicite, une méthode qui se contente de fournir le flux numérique sur l’arête limite (comme dans [27]) fonctionne bien. Mais dans le cadre d’un schéma semi-implicite, dans la phase mathématique du schéma, les états fictifs u0 = (ρ0, ρ0Y0, ρ0v0) et uI+1 = (ρI+1, ρI+1YI+1, ρI+1vI+1) sont utilisés dans le calcul des matrices de diffusion ARoe(u0, u1) et ARoe(uI , uI+1). Ainsi, notre schéma semi-implicite nous oblige à connaitre explicitement les états fictifs, ce qui n’est pas sans poser de problèmes. Nous décidons de nous inspirer, pour la prise en compte des conditions aux limites, de la stratégie qui a fonctionné pour le schéma Tacite. Cette stratégie est globalement fondée sur le signe des valeurs propres. De plus, nous imposons que les variables de relaxation correspondant aux états fictifs doivent être à l’équilibre. Nous allons voir que, concernant la prise en compte de la condition limite amont, la méthode implique que l’état fictif amont est solution d’un problème non-linéaire, ce qui oblige à utiliser une méthode de résolution approchée (par exemple une méthode de Newton). A l’heure actuelle, il n’existe pas de théorème d’existence de l’état fictif amont dans le cadre des écoulements diphasiques et la méthode de Newton n’est pas garantie de converger. En ce qui concerne la condition limite aval, nous imposons une condition dite de “non-retour de liquide” qui n’est pas usuelle. Cette condition, difficile à prendre en compte numériquement, est le sujet de résultats théoriques présentés par Gallou¨et dans [29] et qui sont l’application de résultats dus ˆ à Vovelle [62]. Ils montrent en particulier que la stratégie utilisée ici fonctionne correctement, en particulier parceque la fonction flux numérique reste continue en fonction de l’état fictif.