Analyse mathématique et numérique d’écoulements de fluides à seuil

Les lois de conservation

L’idée d’une certaine conservation de la masse est apparue très tôt dans l’Antiquité, tandis que la conservation du mouvement est apparue beaucoup plus tard (voir par exemple le premier chapitre de [121]), notamment avec la seconde loi de Newton. Cette seconde loi de Newton a été difficile à transposer de la mécanique du point à celle des fluides. Les concepts de pression, d’abord hydrostatique, (cf. Pascal) puis dynamique (cf. Bernoulli) permirent cette généralisation. Ces lois de conservation ont été présentées dans de nombreux livres de synthèse (par exemple [46, 75]). On se place dans toute la suite sur un domaine Ω ∈ R d . — Conservation de la masse : Il faut attendre Euler pour avoir une formulation moderne de cette relation : ∂tρ + div(ρu) = 0, sur Ω (1.1) avec ρ la densité du matériau et u le champ de vitesse. Dans toute la suite, on travaillera avec des fluides incompressibles où ρ est constante. La conservation de la masse se simplifie alors en : div(u) = 0. (1.2) — Conservation de la quantité de mouvement : ρ (∂tu + (u · ∇)u)) + ∇p = divτ + f, sur Ω (1.3) avec p la pression, τ le déviateur du tenseur des contraintes qui s’exerce sur le fluide et f les forces volumiques extérieures. Par ailleurs, la conservation des moments nous donne que τi,j = τj,i. Remarque 1 Le tenseur total des contraintes est le tenseur qui réunit la pression et τ : σ = −pId + τ avec p = −T r(σ) et τ = σ − 1 d T r(σ)Id p représente la moyenne des contraintes qui sont des contraintes normales tandis que τ est appelé déviateur des contraintes. — Conservation de l’énergie : redondante avec les deux premières dans le cas isotherme, cadre dans lequel nous allons travailler. Remarque 2 Dans tout le manuscrit, on se place dans le cas incompressible. Il est à noter que la pression joue un rôle bien différent dans ce cas par rapport au cas compressible où elle est juste une fonction de la densité. Dans le cadre incompressible, la pression est liée à div(u) = 0 qui est ici vu comme une contrainte. Dans la formulation faible et les méthodes de dualité, la pression sera en fait vue comme un multiplicateur de Lagrange associé à cette contrainte. Cependant, il est possible de passer de la formulation compressible à celle incompressible par une limite à faible nombre de Mach [75]. La limite à faible nombre de Mach consiste à prendre les équations adimensionnées et à prendre une limite à petite vitesse et petit temps (t˜ = t/ε, u˜ = u/ε et ε → 0). Dans ce cas, les solutions de Stokes et Navier-Stokes compressibles tendent vers leur pendant incompressible. Une justification de ce passage à la limite a été également montrée dans le cas Stokes notamment par Danchin [41] et dans le cas Navier-Stokes notamment par Feireisl et Novotný [48], pour ne citer que ces articles. Les équations de conservation du mouvement sous leur forme adimensionnée se réécrivent : Re (∂tu + (u · ∇)u)) + ∇p = div(τ ) + f, sur Ω (1.4) où les contraintes et les forces extérieures ont été adimensionnées par les contraintes visqueuses ηU0 L et où le nombre de Reynolds Re est le rapport ρU0L η (le rapport des forces inertielles ρU2 0 sur la contrainte visqueuse). Ici, ρ, η, U0 et L représentent respectivement la densité, la viscosité, la vitesse caractéristique du fluide et la longueur caractéristique du domaine. À petites vitesses, les termes convectifs deviennent négligeables et le système (1.4) devient le système linéaire de Stokes (avec l’ajout de la contrainte d’incompressibilité) : ( div(u) = 0, sur Ω ∇p = div(τ ) + f, sur Ω (1.5) Ces équations ne sont pas suffisantes pour décrire totalement le mouvement d’un fluide. En effet, il y a d+1 équations (en dimension d) pour d + 1 + d(d+1) 2 inconnues. Il y a donc besoin de d(d+1) 2 équations supplémentaires pour fermer le système. Par ailleurs, si les lois de conservation suffisaient à rendre compte du mouvement d’un matériau quelconque, il serait alors impossible de rendre compte de la richesse des comportements et phénomènes des différents fluides et matériaux. Pour rendre compte de certains phénomènes, les équations manquantes vont faire le lien entre les contraintes et les autres composantes du fluide. Notamment, comme le tenseur des contraintes est symétrique, il ne va dépendre que du tenseur symétrique du gradient de vitesse. On pose D(u) = ∇u+(∇u) T 2 . Remarque 3 D(u) est appelé le taux de déformation du fluide. Une raison à cela est donnée par la propriété suivante : Propriété 1 Soit une fonction f : Ω → R d , C 2 telle que pour tout x ∈ Ω, ∇f(x) est antisymétrique. Alors ∇f est une fonction constante. Preuve : Soit i, j, k 6 d. Alors : ∂ 2 i,jfk = ∂ 2 j,ifk = −∂ 2 j,kfi = −∂ 2 k,jfi = ∂ 2 k,ifj = ∂ 2 i,kfj = −∂ 2 i,jfk Et donc, ∂ 2 i,jf = 0 et f est une fonction affine.  On appelle dans la suite R l’ensemble de ces mouvements (les mouvements dits rigides) : R = {f : x 7→ Ax + b, avec A ∈ An(R), b ∈ R d }, la partie Ax représente les mouvements de rotation rigide et b les mouvements de translation. 1.1.2 Les lois constitutives Ces lois ont vocation à décrire des comportements observables macroscopiquement et dûs à la structure interne microscopique du fluide. Dès que les éléments d’un matériau sont soumis à une force extérieure, ils se déplacent et la force de contrainte interne évolue en fonction de ce déplacement. Étudier précisément l’évolution de chaque particule au cours du temps serait beaucoup trop long, surtout dans le cas de fluides complexes très inhomogènes composés de particules de différentes formes et tailles. De fait, une loi constitutive pour passer du discret au continu permet de rendre compte de cette spécificité microscopique tout en gardant un point de vue macroscopique (cf. le chapitre 1 de [39]). Les fluides newtoniens La première loi constitutive a été proposée historiquement par Newton dans son « Principes mathématiques de la philosophie naturelle ». L’hypothèse principale de sa neuvième section est ainsi traduite : « La résistance qui vient du défaut de lubricité des parties d’un fluide doit être, toutes choses égales, proportionnelle à la vitesse avec laquelle les parties de ce fluide peuvent être séparées les unes des autres. » Le « défaut de lubricité » étant justement ce que l’on entend aujourd’hui par viscosité et la « résistance » comme cette force surfacique provenant des contraintes. La viscosité est synonyme de friction interne, de résistance que le matériau oppose à son propre écoulement. Autrement dit, si l’on suppose que l’on se place dans un écoulement unidimensionnel, la relation constitutive serait de la forme : τxy = η∂yu. 4 CHAPITRE 1. INTRODUCTION – ÉTAT DE L’ART La formulation tensorielle a été proposée par Stokes dans les années 1850 : τ = 2ηD(u), où η est la viscosité du fluide, ici supposée constante. Dans ce cadre, on peut réinjecter la loi constitutive dans l’équation (1.3) et retomber sur l’équation traditionnelle de Navier-Stokes. Généralisation : de nouvelles lois constitutives Par la suite, de nombreuses couches de complexité ont été rajoutées pour rendre compte de certains phénomènes. Cependant, la formulation tensorielle d’une loi constitutive est généralement difficile à déduire d’un écoulement unidimensionnel. Elle doit respecter un certain principe d’invariance et plus précisément, la forme générale doit être indépendante du repère de référence (voir [95, 67, 6]). Le tenseur σ, étant symétrique, est donc diagonalisable et peut être représenté par ses trois valeurs propres (σi)i63 dans R 3 , les σi étant appelées les contraintes principales et les vecteurs propres associées à ces valeurs propres étant les directions principales du tenseur. Cette condition d’invariance signifie qu’une loi constitutive doit dépendre de ces σi et ce, de manière symétrique. En particulier, dans la suite, pour un tenseur σ, on peut définir ses trois invariants principaux. — Le premier invariant : I1 = tr(σ) n’est autre que la pression hydrostatique dans le cas qui nous intéresse. — Le second invariant : I2 = −(σ1σ2 + σ2σ3 + σ3σ1). — Le troisième invariant : I3 = σ1σ2σ3. Les contraintes principales sont en fait les racines du polynôme X3 − I1X2 − I2X − I3. Exprimé en termes des contraintes déviatoriques principales, le premier invariant est identiquement nul et son second invariant deviennent J2 = −(τ1τ2 + τ2τ3 + τ3τ1) = 1 2 (τ 2 1 + τ 2 2 + τ 2 3 ), qui n’est autre que la norme du tenseur des contraintes. Les fluides quasi-newtoniens et rhéo-fluidifiants/épaississants La première généralisation a été celle des fluides quasi-newtoniens. Il y a encore en tout point une proportionnalité entre les directions et valeurs propres des deux tenseurs mais la viscosité n’est plus supposée constante. Les modèles les plus connus sont ceux où η suit une loi puissance : τ = 2α ∗ |D(u)| r−1D(u), avec r > 0. Dans le cas où r < 1, on parle de fluide rhéo-fluidifiant (comme des émulsions ou des suspensions) : plus le fluide se déforme, plus sa viscosité apparente diminue. Dans le cas où r > 1, on parle de fluide rhéo-épaississant (comme le sable humide) : la viscosité va croître au fur et à mesure que le cisaillement croît. Viscoélasticité et thixotropie Jusqu’à présent, on n’a considéré que des situations où les équations constitutives sont indépendantes du temps. Cependant, il peut arriver également que l’épaississement ou la fluidification du fluide ne soit pas un effet instantané mais arrive progressivement avec un certain retard. Un matériau présente un comportement thixotropique si ce temps de réponse du matériau est non nul (voir par exemple [39]). Par ailleurs, la loi de Hooke traditionnelle qui donne la relation constitutive pour un solide peut être adaptée dans le cas de fluides. Elle mélange le caractère visqueux du matériau avec le fait que les particules microscopiques forment un système élastique. Ces comportements peuvent s’observer macroscopiquement par un phénomène d’irréversibilité ou d’hystéresis au niveau de la viscosité. Les matériaux gardent en mémoire dans leur structure microscopique ce qu’ils ont vécu. Dans les matériaux que nous considérerons, on négligera ces effets. 

Visco-plasticité et fluides à seuil Phénoménologie

Un fluide est dit à seuil s’il ne se déforme que lorsque les contraintes qui s’exercent sur lui dépassent un seuil donné. Si à un endroit ces contraintes sont en dessous de cette limite, le matériau ne se déforme pas et se comporte comme un solide rigide. De nombreux exemples sont présents dans la nature, dans les laboratoires et les applications industrielles sont également nombreuses. En fait, une bonne partie des fluides que l’on rencontre au quotidien présentent cet aspect à seuil : — Produits alimentaires et domestiques : pâtes (ketchup, dentifrice,…), gel à cheveu… — Matériaux géophysiques : argiles, boues, neiges denses… — Matières biologiques : Mucus, sang… — Produits industriels : Huiles, pétrole, béton… Ces matériaux présentent plusieurs caractéristiques phénoménologiques communes : — ils ne se déforment qu’aux endroits où une contrainte suffisante s’exerce sur eux, — ils se décomposent en une zone où ils se déforment (comme un fluide) et une zone où ils se déplacent de manière rigide (comme un solide) — ils sont immobiles si aucune force majeure ne s’exerce sur eux. Ce dernier point peut se voir par exemple dans l’expérience suivante : si on fait tomber un de ces matériaux sur un plan incliné, il va s’arrêter au bout d’un temps assez long. Équation constitutive L’équation constitutive qui nous intéresse a été d’abord dérivée de manière empirique par Eugène Bingham [20] et Théodore Schwedoff [117] pour des écoulements unidimensionnels. τ = 2ηγ˙ + τy γ˙ |γ˙ | , ssi γ˙ 6= 0, |τ | 6 τy, ssi γ˙ = 0, où γ˙ est le cisaillement du matériau. À l’inverse, cette loi représente une viscosité infinie pour des petites contraintes (inférieures à la contrainte seuil), tandis qu’elle est finie au delà de ce seuil. Il faudra attendre Prager [68, 99, 100] et Il’iushin [71] pour proposer une généralisation de cette loi en une formulation tensorielle pour un écoulement multi-dimensionnel : τ = 2ηD(u) + τy D(u) |D(u)| , ssi D(u) 6= 0, |τ | 6 τy, ssi D(u) = 0, (1.6) où l’on utilise la convention suivante pour la norme d’un tenseur de dimension d ainsi que pour le produit scalaire qui en découle : |T| 2 = 1 d X i,j T 2 i,j , T : T 0 = 1 d X i,j Ti,jT 0 i,j 6 CHAPITRE 1. INTRODUCTION – ÉTAT DE L’ART La version adimensionnée du modèle de Bingham est : τ = 2D(u) + B D(u) |D(u)| , ssi D(u) 6= 0, |τ | 6 B, ssi D(u) = 0, (1.7) où B = τyR ηU est le nombre de Bingham (avec U et R respectivement la vitesse et la taille caractéristiques du système). Il s’agit du rapport entre τy la contrainte seuil et ηU R la pression caractéristique du système. Remarque 4 On peut généraliser cette condition de seuil et l’exprimer en fonction (entre autres) du deuxième invariant du tenseur (voir par exemple [67, 6]). Géométriquement, ces relations peuvent se voir en représentant les contraintes comme des vecteurs de R 3 : à un tenseur σ, on associe le point S : (σ1, σ2, σ3) et à son déviateur, le point P : (τ1, τ2, τ3). Dans cet espace R 3 , on peut définir une « zone fluide » et une « zone rigide », la frontière entre les deux étant appelée surface limite (de plasticité). Cette surface limite est une ligne de niveau de la forme f(I1, I2, I3) = 0, les zones f > 0 et f < 0 représentant respectivement les parties fluide et rigide. De par la relation τ1 + τ2 + τ3 = 0, on peut se restreindre à considérer la surface limite dans le plan V ect(1, 1, 1). Dans la suite, on travaille avec le critère de von Mises : f(J2) = p J2 − τy mais d’autres critères comme celui de Tresca ou de Mohr-Coulomb (qui dépendent du troisième invariant) nous donneraient un critère de seuil différent et seraient représentés différemment (voir figure 1.1). Après avoir défini le critère de seuil, il faut définir ce qui se passe dans la zone Figure 1.1 – Représentation géométrique des contraintes principales (tiré des figures 1 et 3 de [67]). À droite, interprétation géométrique des conditions de Tresca et von Mises. fluide (la loi constitutive proprement dite). Deux règles (principalement empiriques) régissent le comportement du fluide (voir [67]) : — le principe de coaxialité (i.e. la colinéarité des directions principales de taux de déformation D(u) et du déviateur des contraintes), — l’associate normal flow rule (i.e. la proportionnalité entre le taux de déformation et l’excès de contraintes au delà de τy). Ces deux règles amènent à la loi D(u) = λ( √ J2 − τy)∇f (où λ est un réel) ce qui donne la loi D(u) =    λ 2 ( √ J2 − τy) J2 τ ssi f > 0 0 sinon. = λ 2 ( √ J2 − τy)+ J2   pour von Mises, i.e. exactement la loi de Bingham lorsque λ est constante. En pratique, très peu de fluides suivent la loi de Bingham, quand bien même ils présentent un effet seuil. La loi de Bingham peut être généralisée (voir par exemple [39]) en choisissant une viscosité dépendant de certains paramètres (la déformation, la pression…). En particulier, la loi de Herschel-Bulkley [21] propose une viscosité en loi puissance de la déformation : τ = 2K|D(u)| n−1D(u) + τy D(u) |D(u)| , ssi D(u) 6= 0, |τ | 6 τy, ssi D(u) = 0. On retrouve les caractères respectivement rhéo-fluidifiants et rhéo-épaississants suivant la valeur de n ainsi que le cas Bingham lorsque n = 1. Le modèle de Casson est un autre exemple de fluide à seuil moins répandu et plus utilisé pour modéliser des particules en suspension. On ne fera que l’évoquer : τ = √ η + s τy |D(u)| !2 D(u), ssi D(u) 6= 0, |τ | 6 τy, ssi D(u) = 0. Remarque 5 De manière équivalente, on peut exprimer la déformation en fonction des contraintes en inversant la relation précédente, ce qui donne dans le cas de la loi de Bingham : 2η(D(u))D(u) =  1 − τy |τ |  τ (1.8) Remarque 6 Une autre manière d’écrire cette loi constitutive serait de dire que τ est dans le sous différentiel (voir chapitre 2 pour la définition) d’une fonctionnelle (cf. [46]) : τ ∈ ∂W(D(u)), où W(E) = η|E| 2 + τy|E| (1.9) 

Détermination du seuil et controverse

L’un des problèmes principaux posés par ce genre de modèle est donc la forme de ces zones rigides ainsi qu’une manière de les déterminer, que ce soit expérimentalement, analytiquement ou numériquement. Pour ce faire, il faut déterminer expérimentalement le seuil τy ce qui est une tâche délicate encore maintenant (le lecteur intéressé pourra lire l’article de synthèse de Bonn et co-auteurs [23] sur ce sujet). Plus encore, un débat remontant à plusieurs décennies porte sur l’existence même de « vraies » régions rigides. Barnes et Walters [11] ont, avec un titre quelque peu provocateur, exprimé les premiers doutes sur cette existence, arguant que les zones rigides observées n’étaient en réalité que des zones où la déformation était petite, trop petite pour la précision des rhéomètres de l’époque. Cette controverse fut également alimentée par le paradoxe de lubrification que l’on définira section 

La controverse a été forte les décennies suivantes  pour ne citer que quelques contributions)

Le but ici n’est pas d’entrer dans le débat mais juste de le citer pour son importance historique. Nous considérerons ici, que quand bien même il y a écoulement visqueux pour de petites contraintes, la viscosité est suffisamment élevée pour que le fluide soit considéré comme un fluide à seuil idéal (et donc qu’il est pertinent de considérer des modèles de Bingham et Herschel-Bulkley