Modèle à dissipation
Reprenons le modèle à dissipation donnée par l’équation (25) :
$ \nabla \big( cc_g \nabla \Psi \big)+ k^2c c_g \Psi = -j \omega \big(W_{def}+W_{fon}W_{poreux}\big)\Psi $ (30)
où Wdef, Wfond et Wporeux sont la dissipation par respectivement déferlement, frottement sur le fond et percolation dans les fonds poreux.
Dissipation par déferlement
Le déferlement est un phénomène dissipatif de l’énergie de la houle. Il est conditionné par la cambrure, rapport de la hauteur de houle sur sa longueur d’onde. La cambrure des vagues augmente au fur et à mesure que la houle s’approche du rivage jusqu’à une valeur maximale appelée cambrure limite. Il y a apparition du déferlement lorsque la cambrure atteint cette valeur limite.
Le déferlement engendre un écoulement complexe turbulent ; aujourd’hui encore il fait l’objet de nombreux travaux de recherche. On distingue cinq types de déferlement :
- Le moutonnement : au large, l’action du vent augmente la cambrure et fait déferler la houle ;
- Le déferlement glissant (spilling breaking) : en faible profondeur sur une pente faible, la vague s’écroule progressivement sur sa face avant ;
- Le déferlement plongeant (plunging breaking) : sur des pentes un peu plus fortes, la vague forme des rouleaux plongeants, bien connus des surfeurs ;
- Le déferlement frontal (surging breaking) : sur des pentes très prononcées, il se manifeste par l’écoulement d’un front d’onde ;
- Le mascaret : il concerne les ondes très longues.
Il existe dans la littérature de nombreuses formulations de la hauteur maximale de la houle avant déferlement Hd. Toutes obtenues de façon expérimentale, elles sont plus ou moins complexes et ont généralement un domaine de validité. Elles sont toutes fonction de la profondeur h, souvent de la longueur d’onde $ \lambda $ et parfois de la pente des fonds m.
Pour les introduire dans un modèle numérique, une première méthode consiste à égaliser les hauteurs calculées par le code de calcul à celle obtenue par la formule de déferlement. C’est la méthode dite du simple écrêtage. Une seconde méthode consiste à calculer un terme de dissipation d’énergie que l’on ajoute dans l’équation de pente douce (dite de Berkhoff).
Déferlement par écrêtage
La prise en compte du déferlement est conditionnée par le choix de la formule d’écrêtage. Les formules d’écrêtage les plus courantes sont les suivantes :
- La formule de Miche (1944) égale à $ H_d = 0,142 \lambda \tanh \Big(2\Pi\frac{h}{\lambda}\Big) $
- La formule de Munk (1949) en eaux plus profondes égale à $ H_d =0,78 h $
- La formule de Goda (1975) sans effet de pente égale à $ H_d = 0,17\Big(1-\exp \Big[-\frac{1,5\Pi h}{\lambda}\Big]\Big) $
- La formule de Goda (1975) avec effet de pente égale à$ H_d = 0,17\Big(1-\exp \Big[-\frac{1,5\Pi h}{\lambda}\big(1+15m^(\frac{4}{3})\big)\Big]\Big) $
Déferlement par dissipation
- La formulation du déferlement est décrite par Battjes et Janssen (1978). Les coefficient Wdef conditionne l’intensité de la dissipation. Il est fonction de la période de houle T, du rapport b de la hauteur de houle calculée sur la hauteur maximale Hd. Cette-dernière est donnée par les formules précédentes.
$ W_{def}=\alpha \frac{Q_b}{T}b^2 $ (31)
avec $ \alpha\, $ un coefficient de calage proche de l’unité et Qb la fraction des vagues déferlantes. Cette fraction vérifie l’équation suivante :
$ 1-Q_b=-b^2\ln Q_b $ (32)
Dissipation par frottement sur les fonds
Dingemans (1997) donne la formule suivante pour la dissipation de la houle par frottement sur les fonds :
$ W_{fond}=\frac{4}{3\Pi}\frac{f_w}{h}\frac{\hat u_s}{\cosh ^2(kh)} $ (33)
avec $ f_w $ le coefficient de frottement. Pour des fonds lisses, certaines formules donnent $ f_w $ en fonction du nombre de Reynolds. Pour des fonds ridés la dissipation par frottement sur le fond est beaucoup plus grande. Le nombre d’onde et la profondeur sont notés k et h respectivement. L’amplitude effective $ \hat u_s $de la vitesse horizontale à la surface u (égale à est obtenue après application de la procédure de linéarisation de Lorentz.
$ \hat u_e = \frac{3\Pi}{8}\frac{\int_{0}^{T}\left|u\right|^3 dt }{\int_{0}^{T}\left|u\right|^2 dt} = \sqrt{a+b}\Big[E(m)-\frac{1}{4}\big(1-\frac{b}{a}\big)K(m)\Big] $ (34)
et $ a= \frac{1}{T}\int_{0}^{T}{\left|u\right|}^3 dt $
$ b= \max_{t\in [T]}\Big(\left|u\right|^2-a\Big) $
$ m = \frac{2b}{2+b} $
K(m) et E(m) représentent respectivement les intégrales elliptiques de première et seconde espèce.
Dissipation par percolation dans les milieux perméables
Les structures maritimes de protection des ports contre la houle sont des digues à talus en enrochements, des digues verticales ou des digues composites. Elles arrêtent la houle par dissipation locale de leur énergie ou par réflexion vers le large. Les digues à talus sont formées en particulier d’un empilement d’enrochements qui monte jusqu’à la surface libre. Les structures de protection du trait de côte comme les brise-lames détachés ou les épis constituent également des ouvrages en enrochements émergés ou submersibles. De nombreux modèles ont été développés pour décrire la transformation des vagues sur un fond poreux de profondeur variable. Les premiers modèles théoriques ont été présentés par Sollit et Cross (1972) qui ont évalué la transmission et la réflexion des vagues en présence d’un brise-lames perméable et ont validé les résultats avec des données expérimentales. L’écoulement est modélisé par l’équation de Forchheimer qui est linéarisée à l’aide de l’hypothèse de Lorentz de travail équivalent. Madsen (1974) et Dalrymple et al. (1991) ont mis au point des modèles voisins appliqués à des brises-lames émergés. Des modèles numériques tridimensionnels ou bidimensionnels vertical prenant en compte des couches poreuses ont aussi été présentés. Gu et Wang (1992) résolvent un écoulement potentiel à l’aide de la méthode intégrale d’éléments de frontière et étudient l’interaction des vagues avec un brise-lames composite dont la base est constituée d’enrochements multicouches. Ropert (1999) étudie également l’écoulement à l’intérieur d’un caisson Jarlan posé sur une base en enrochements à l’aide d’un modèle bidimensionnel vertical d’écoulement potentiel aux éléments finis. Parmi les modèles intégrés sur la verticale, le modèle de Boussinesq de Cruz et al (1997) étendu à la transformation des vagues sur un lit poreux a été validé sur des données expérimentales de transformation de la houle sur une barre poreuse triangulaire et autour de l’ouverture d’un brise-lames poreux submersible. Les modèles horizontaux les plus communs sont les modèles elliptiques dérivés de l’équation de pentes douces de Berkhoff. Rojanakamthorn et al. (1990) ont les premiers adapté l’équation elliptique pour des brise-lames submersibles poreux en incluant le phénomène de déferlement.
Fig. 13 : Croquis ancien de digue composite Le milieu poreux saturé est composé de deux phases : un squelette rigide et un fluide remplissant l’espace interstitiel. Le milieu poreux est inerte : il n‘y a pas de création de matière pour les deux phases. Le squelette est supposé immobile et le fluide incompressible. Les variables dépendent de la position mais leur dépendance est omise pour simplifier les équations. Les équations de conservation de la masse et de la quantité de mouvement sont respectivement :
$ div (\epsilon \vec V_s) = 0 $ (35)
avec $ \epsilon \, $ la porosité, $ \vec V_s $ la vitesse interstitielle, $ \epsilon \vec V_s $ la vitesse de filtration.
$ \frac{\partial \vec V_s}{\partial t}=-\frac{1}{\rho}\vec \nabla (p+ \rho g z)-\vec F_r - \vec F_i $ (36)
où $ \rho\, $ représente la densité du fluide, p la pression, z la coordonnée verticale mesurée positivement vers le haut et $ \vec \nabla $ l’opérateur Gradient. Le terme de résistance $ \vec F_r $s’exprime en fonction de la vitesse interstitielle :
$ \vec F_r = \frac{v\epsilon}{K}\vec V_s + \frac{C_f \epsilon ^2}{\sqrt{k}}\left|\vec V_s\right| V_s $ (37)
où $ \nu $ est la viscosité, $ C_f $ le coefficient de résistance turbulente et K la perméabilité du milieu poreux. Le terme de résistance inertielle $ \vec F_i $dépend du coefficient de masse virtuelle $ C_m $ et de l’accélération dans le sens de l’écoulement (Sollit et Cross (1972)).
$ F_i= \frac{1-\epsilon }{\epsilon} C_m \frac{\partial \vec V_s}{\partial t} $
Le coefficient $ C_m $ est estimé pour un grain individuel de forme régulière mais reste en général inconnu pour un ensemble de grains composant un corps poreux. Remplaçant les expressions $ \vec F_r $et $ \vec F_i $dans l’équation de conservation de la quantité de mouvement (36), nous obtenons l’équation de Forchheimer qui est utilisée en général par les ingénieurs pour décrire l’écoulement dans un milieu poreux soumis à l’action des vagues.
$ S \frac{\partial \vec V_s}{\partial t}= - \frac{1}{\rho}\vec \nabla (p + \rho gz)- \frac{\nu \epsilon }{K}\vec V_s - \frac{C_f \epsilon ^2 }{\sqrt{K}}\left|\vec V_s\right| \vec V_s $ (39)
avec $ S = 1+ \frac{1-\epsilon}{\epsilon}C_m $. Le terme linéaire, le terme de Darcy, gouverne l’équation pour les petits nombres de Reynolds dans un milieu poreux fin. Le terme non-linéaire, le terme de Forchheimer, est une extension du terme de Darcy pour des milieux poreux plus grossiers avec de grands nombres de Reynolds pour lesquels l’effet d’inertie dépasse l’effet de viscosité.
En considérant une excitation monochromatique, les termes dissipatifs, linéaires ou non-linéaires de l’équation de Forchheimer sont remplacés par un terme équivalent linéaire :
$ \frac{\nu \epsilon}{K}\vec V_s +\frac{C_f \epsilon ^2}{\sqrt{K}}\left|\vec V_s\right|\vec V_s \to f\omega \vec V_s $ (40)
Dans cette expression $ f \, $ représente un coefficient de frottement ou d’atténuation, $ \omega \, $ est la pulsation de l’excitation monochromatique de période T. L’évaluation du coefficient $ f \, $ repose sur l’hypothèse de Lorentz de travail équivalent qui consiste à égaler, sur une période de houle, le travail moyen du terme non-linéaire avec le travail moyen d’un terme linéaire. Le terme $ f \, $ est :
$ f_{\omega} = \frac{v\epsilon}{K} \left(1+\frac{C_f \epsilon \sqrt {k}}{v}\bar f\right) $ (41)
avec $ \bar f = \frac{\int_{0}^{T} \left|\vec V_s\right|^2 \,dt }{\int_{0}^{T} \left|\vec V_s\right|^3 \,dt} = \frac{8}{3\pi} \sqrt {a+b}\Big[E(m)- \frac{1}{4}\left(1-\frac{b}{a}\right)K(m) \Big] $ et $ \qquad \qquad \qquad $ $ a= \frac{1}{T}\int_{0}^{T} \left|\vec V_s\right|^2 \,dt $
$ b=\max _{t\in [0,T]}\Big(\left|\vec V_s\right|^2 -a \Big) $,$ \qquad \qquad \qquad $ $ m=\frac{2a}{a+b} $ $ \qquad \qquad \qquad $ $ K(m) $ $ \qquad \qquad \qquad $et $ \qquad \qquad \qquad $$ E(m) $ représentent respectivement les intégrales elliptiques de première et seconde espèce.
Les solutions monochromatiques s’écrivent avec $ p^* = p+\rho gz\, $
$ \vec V_s = Re\big(\vec V_s e^{-j\omega t}\big) $
$ P= Re\big(P e^{-j\omega t}\big) $ (42)
L’expression linéaire de la conservation de la quantité de mouvement en fonction des variables $ \big(\vec V ,P\big) $ et du coefficient de frottement $ f\, $ s’écrit :
$ -j\omega S \vec V = -\frac{1}{\rho}\vec {\nabla P}-f\omega \vec V $ (43)
L’équation peut aussi s’écrire sous la forme donnant la vitesse de filtration :
$ \epsilon \vec V = -\frac{1}{p j \omega}\beta \vec {\nabla P} $ (44)
avec $ \beta = \frac{\epsilon}{S+jf} $. L’équation de conservation de la masse permet d’obtenir l’équation finale :
$ div \big(\beta \vec \nabla P\big) $ (45)
Cette équation, avec les conditions aux limites à la surface libre et sur le fond rigide, définit un modèle avec la pression comme inconnue principale. Ce modèle représente la propagation de vagues linéaires monochromatiques au-dessus d’un fond poreux. La résolution de cette équation par un schéma itératif donne les champs de vitesse et de pression dans le milieu poreux. Cette équation est commune aux milieux fluide et poreux en considérant un milieu hétérogène avec une fonction $ \beta \, $ variable dans l’espace. Cette fonction est constante et égale à 1 dans le milieu fluide et est variable dans le milieu poreux. De façon à simplifier les équations, nous choisissons un écoulement suivant l’axe x. Pour une excitation monochromatique nous écrivons des nouvelles expressions condensées pour les coefficients a et b :
$ a = \frac{1}{2}\Big(\left|\vec V_x\right|^2 + \left|\vec V_z\right|^2 \Big) $
$ b = \frac{1}{2}\Big(V_x^2 + V_z^2 \Big) $ (46)
avec $ V_x\, $ et $ V_z\, $ les composantes horizontale et verticale respectivement de la vitesse complexe $ \vec V\, $
Fig. 14 : Schéma d’un bicouche.
Un bicouche, représenté sur la figure 13, est composé d’un milieu fluide submergeant un milieu poreux qui repose sur un fond rigide. Les axes orthogonaux Ox et Oy définissentun plan horizontal qui coïncide avec la surface libre. La profondeur constante du fluide est notée $ h_p\, $ et la profondeur constante totale est notée $ h\, $. Le déplacement de la surface libre $ \mu (x,t)\, $ est écrit pour des solutions monochromatiques :
$ \mu = Re(Ae^{-j\omega t}) $ (47)
A la surface libre, les conditions linéaires dynamiques et cinématiques avec l’équation de Bernouilli linéarisée donnent en z = 0 :
$ \frac{\partial P}{\partial z}- \frac{\omega ^2}{ g}P = 0 $ et $ P= \rho g A_, $ (48)
Sur le fond plat rigide, la composante verticale de la vitesse $ \vec V $est nulle de telle sorte que $ \frac{\partial P}{\partial z} = 0 $. Le système d’équations satisfait par la pression P est finalement :
$ \begin{cases} div \big(\beta \vec {\nabla P}\big)= 0 \qquad\qquad -h\le z \le 0 \\ \frac{\partial P}{\partial z}- \frac{\omega ^2}{g}P = 0 \qquad\qquad z = 0 \\ \frac{\partial P}{\partial z} = 0 \qquad\qquad\qquad z = -h \end{cases} $ (49)
Avec un bicouche, la quantité $ \beta \, $varie avec le coefficient de frottement $ f\, $. Le coefficient de frottement $ f\, $ est fonction des coefficients $ a\, $ et $ b\, $ qui sont eux-mêmes des fonctions de la vitesse interstitielle $ \vec V \, $et de la pression P. Comme la pression P varie spatialement, c’est aussi le cas de $ \beta \, $ qui varie en particulier le long de l’axe vertical.
Nous cherchons une solution du système d’équations (49) sous la forme de variables séparées :
$ P(x,z) = \phi (x) Z(z) $ (50)
Dans la couche poreuse et la couche fluide, la première équation du système (49) donne, avec le nombre d’onde k identique dans les deux couches :
$ \frac{1}{\phi}\frac{d^2\phi}{d x^2}= \frac{\big(\beta z'\big)'}{\beta Z}= k^2 $ (51)
Le coefficient $ \beta $ est moyenné dans chaque couche horizontale de telle sorte qu’il est supposé dépendre seulement de la coordonnée verticale z. Nous réécrivons l’équation (51) sous la forme suivante :
$ \frac{d^2\phi}{d x^2}+ k^2\phi = 0 $
$ \Big(\beta Z'\Big)'-k^2\beta Z = 0 $ (52)
La seconde équation demande la continuité de $ \beta Z'\, $ et de Z. Nous choisissons Z(0) égal à 1 ce qui donne $ \phi (x) = \rho gA(x) = 0\, $ avec l’équation (50). L’amplitude de houle vérifie donc également l’équation :
$ \frac{d^2A}{dx^2}+k^2 A = 0 $ (53)
Nous choisissons une houle plane incidente se propageant selon l’axe des x avec une amplitude $ A_0\, $ en x=0.
$ A(x) = A_0 e^{jkx}\, $ (54)
La pression P prend donc la forme finale :
$ P(x,y,z) = \rho g A_0 e^{jkx} Z(z)\, $ (55)
Avec la célérité c égale à $ \frac{\omega}{k}\, $, nous en déduisons que :
$ V_x = -\frac{g}{c}\frac{\beta}{\epsilon}A_0 Z(z)e^{jkx} $ (56)
$ V_z = -\frac{g}{j \omega}\frac{\beta}{\epsilon}A_0 Z'(z)e^{jkx} $ (57)
$ a = \frac{1}{2}\left(\frac{g}{\omega}\frac{\beta}{\epsilon}\right)^2 A_0^2 \Big(\left|K(z) Z\right|^2 + \left|Z'(s)\right|^2 \Big) \left|e^{jkx}\right|^2 $ (58)
$ a = \frac{1}{2}\left(\frac{g}{\omega}\frac{\beta}{\epsilon}\right)^2 A_0^2 \left|K^2(z) Z^2 + Z'^2(s)\right| \left|e^{jkx}\right|^2 $ (59)
Le coefficient $ beta $ est moyenné dans chaque couche horizontale sur un segment de longueur $ [0,L_x]\, $ avec la longueur $ L_x\, $ supposée petite de telle sorte que $ \left|e^{jkL_x}\right|\approx 1 $.
Figure 15 : Procédure de résolution.
La solution de Darcy $ (\bar f = 0)\, $ est d’abord calculée. La procédure consiste à trouver le vecteur propre $ Z(z)\, $ associé avec la première valeur propre et le nombre d’onde complexe $ k\, $. En donnant l’amplitude de houle incidente $ A_0\, $, les coefficients a et b sont calculés et utilisés pour actualiser la valeur de $ \bar f\, $. Cette opération est répétée aussi longtemps que la différence relative $ \bar f\, $ de entre deux itérations successives excède $ 10^{-6}\, $ . Le schéma de la procédure est synthétisé sur la figure 14.
Le système (45) est réécrit en fonction de la fonction Z :
$ \begin{cases} \Big(\beta Z'\Big)' - k^2 \beta Z = 0 \qquad -h \le z \le 0 \\ Z' - \frac{\omega ^2}{g}Z \qquad\qquad\qquad z = 0 \\ Z'= 0 \qquad\qquad\qquad z = -h \end{cases} $ (60)
La pulsation $ \omega \, $ et la quantité $ \beta \, $, fonction de $ z\, $, sont données. Nous cherchons les vecteurs propres $ Z\, $ et valeurs propres $ K, $ du système (54) appelé aussi problème de Sturm-Liouville. Dans ce système, les fonctions $ Z'\, $ et $ K $ ne sont pas continues. Il y a en particulier un saut de $ Z' $ et de $ \beta \, $ à l’interface entre le milieu fluide et le milieu poreux. Cependant $ \beta Z'\, $ et $ Z\, $ sont des fonctions continues. Cela implique que $ \epsilon V_z \, $et $ P\, $ sont des fonctions continues. La solution du système (60) dans une couche $ i\, $ avec $ \beta \, $ constant par couche est :
$ Z(z) = A_i \cosh (kz) + B_i \sinh (kz)\, $ (61)
Losada et al. (1996) utilisent cette méthode avec une couche unique du milieu poreux. Ils développent l’équation de pente douce pour l’évolution des vagues sur un brise-lames poreux submersible. Ils ont enfin mené une analyse complète des résultats intégrant les effets d’oblicité, de géométrie du brise-lames et des caractéristiques du matériau poreux. Silva et al. (2002) ont ajouté de nouveaux termes à l’équation de pente douce étendue au milieu poreux de façon à prendre en compte des profondeurs rapidement variables. Houari (2002) montre que cette équation de pente douce modifiée revient à l’équation de pente douce avec dissipation (30) lorsque la dissipation par percolation est modérée. Dans tous ces travaux, le coefficient de frottement f est obtenu numériquement par itérations et est choisi moyenné sur le volume. Houari et al. (2003)utilisent un coefficient de frottement variable sur un plan horizontal mais moyenné sur l’axe vertical.
Pour un grand nombre de couches, la méthode basée sur l’équation (57) devient complexe et lourde en temps de calcul. C’est la raison pour laquelle Sergent et Duhamel (2006) préfèrent écrire la forme variationnelle de (56) avec une fonction test $ \beta Z\, $. Ils évaluent ainsi les variations du coefficient d’atténuation f sur l’axe vertical, la nouvelle relation de dispersion et les variations du profil de pression pour un coefficient d’atténuation variable.