S'abonner à un flux RSS
 

Utilisateur:Jean-Michel Tanguy/SujetENTPE2020/BRETON-NGUYEN-SALLES : Différence entre versions

De Wikibardig
(Solution semi-analytique (homotopie))
(Solution semi-analytique (homotopie))
Ligne 280 : Ligne 280 :
  
 
[[File:cas5_biais.jpg]]
 
[[File:cas5_biais.jpg]]
 
+
[[Fichier:Cas5 biais.jpg|800px|vignette|gauche]]
  
 
Et de faire un changement de repère de façon à obtenir :
 
Et de faire un changement de repère de façon à obtenir :

Version du 22 juin 2020 à 14:54


Vidéo d'introductionpar Wikhydro Modélisation de la houle par homotopie : une brève introduction

Notre travail a pour but de modéliser l'impact du changement climatique sur les côtes, en quantifiant l'impact des houles sur le littoral.

Pour modéliser les houles, nous utilisons le modèle de Berkhoff aussi appelé équation de pente douce, obtenu en 1972. Ce modèle est une équation aux dérivées partielles, que nous nous proposons de résoudre par une méthode analytique pour les cas simples et semi-analytique pour les cas plus complexes.

Nous nous placerons dans l'hypothèse d'ondes longues et nous étudierons différents cas concrets de géométries et de conditions aux limites.

Sommaire

Cas n°1 : canal uniforme unidimensionnel avec sortie libre en amont

Cette situation traite un canal uniforme unidimensionnel plat de longueur L avec l'entrée d'une onde de fréquence unitaire $ ϕ=1 $ par l'aval et une sortie libre en amont $ ϕ_x=ikϕ $.

Dans ce cas précis, l'équation de Berkhoff se simplifie comme suit: $ \frac{\partial^2 ϕ}{\partial x^2} +k^2ϕ = 0 $


Solution analytique

Nous résolvons l'équation caractéristique associée : $ X^2 + k^2 = 0 $
$ Δ = -4k^2 $ d'où $ λ_1 = -ik $ et $ λ_2 = ik $
donc la solution est de la forme $ Ae^{-ikx}+ Be^{ikx} $, avec A et B des constantes.

Détermination des constantes A et B grâce aux conditions limites :

  • $ ϕ(x=0)=1 $$ A+B=1 $
  • $ \dfrac{\partial ϕ}{\partial x}(x=0)=-ik $$ A=1 $ et $ B=0 $


donc $ ϕ(x)=e^{ikx} $
$ ϕ(x,t)=e^{i(kx-\omega t)} $
or $ H=\Re(ϕ) $ donc $ \fbox{h(x,t)=cos(kx-ωt)} $

N.B. : sur les graphiques, la valeur $ h(x)=0 $ représente le niveau moyen de la surface.

Houle cas1.gif

Solution semi-analytique (homotopie)

Le principe de la méthode d'homotopie est brièvement expliqué dans cette vidéo. Il s'agit de "passer" d'une fonction à une autre au travers d'un paramètre p variant de 0 à 1. Lorsque p vaut 0 et 1, $ h(x) $ est confondue avec respectivement $ f(x) $ et $ g(x) $.

Homotopie-p2-1-compressor.gif

Dans ce cas, la relation d'homotopie (ayant pour paramètre p entre 0 et 1) s'exprime comme ceci (on choisit la dérivée seconde comme opérateur linéaire et on part d'une solution initiale nulle): $ (1−p)(ϕ_{xx}-u_{0,xx})+p(ϕ_{xx}+k^2ϕ)=0 $
Or:

  • $ ϕ(x,p) = \sum_{n=0}^∞ p^nϕ_n = ϕ_0(x)+pϕ_1(x)+p^2ϕ_2(x)+p^3ϕ_3(x)+... $
  • $ ϕ_{xx}(x,p) = ϕ_{0,xx}(x)+pϕ_{1,xx}(x)+p^2ϕ_{2,xx}(x)+p^3ϕ_{3,xx}(x)+... $

D'où : $ (1−p)((ϕ_{0,xx}(x)+pϕ_{1,xx}(x)+p^2ϕ_{2,xx}(x)+p^3ϕ_{3,xx}(x)+...)-u_{0,xx})+p[ϕ_{0,xx}(x)+pϕ_{1,xx}(x)+p^2ϕ_{2,xx}(x)+p^3ϕ_{3,xx}(x)+...+k^2(ϕ_0(x)+pϕ_1(x)+p^2ϕ_2(x)+p^3ϕ_3(x)+...)]=0 $

On développe ensuite en fonction des ordres de p.
Ordre 0 :
$ ϕ_{0,xx}-u_{0,xx}=0 $$ \dfrac{\partial^2 ϕ_0}{\partial x^2} = \dfrac{\partial^2 u_0}{\partial x^2} $$ \dfrac{\partial ϕ_0}{\partial x} = \dfrac{\partial u_0}{\partial x} + A $$ ϕ_0(x) = Ax + u_0 + B $


Détermination des constantes A et B :

  • $ ϕ(0) = 1 $$ 1 = u_0+B $ or $ u_0=1 $ donc $ B=0 $
  • $ \dfrac{\partial ϕ}{\partial x} = ikϕ $ or $ ϕ = Ax + B + u_0 $ et $ \dfrac{\partial ϕ}{\partial x} = A $
    d'où pour $ x=L $, $ A = ik(AL+B+u_0) $$ A = \dfrac{ik}{1-ikL} $


Donc $ ϕ_0(x) = u_0+\dfrac{iku_0x}{1-ikL} = 1 + \dfrac{ikx}{1-ikL} $
or $ kL=1 $ ie $ L=\dfrac{1}{k} $donc $ ϕ_0(x) = 1 - \dfrac{kx}{2} + \dfrac{ikx}{2} $

On a : $ H = \Re(ϕ_0) = 1+\dfrac{kx}{2} $
donc $ \boxed{h(x,t)= 1+\dfrac {kx}{2} cos(\omega t)} $


Ordre 1 :
$ \dfrac{\partial^2 ϕ_1}{\partial x^2} = -\dfrac{\partial^2 u_0}{\partial x^2} - k^2 ϕ_0 $ $ \Rightarrow $ $ \dfrac{\partial ϕ_1}{\partial x} = -\dfrac{\partial u_0}{\partial x} - k^2 \int ϕ_0 + A_1 $ donc $ ϕ_1(x)=-u_0-k^2 \iint ϕ_0 +A_1x + B_1 $
or $ \iint ϕ_0 = -\dfrac{u_0x^2}{2}+\dfrac{iku_0x^3}{6(1-ikL)} $
donc pour $ u_0=1 $: $ ϕ_1(x)=-1-\dfrac{k^2x^2}{2}+\dfrac{ik^3x^3}{6(1-ikL)} + A_1x + B_1 $

  • $ x=0 \rightarrow \phi_1(x)=1 \rightarrow B_1=0 $
  • $ x=L : \dfrac{\partial ϕ}{\partial x}=ik\phi=-\dfrac{17}{12}ik+\dfrac{k}{12}+A_1i~~(1) $

Or on a également $ \dfrac{\partial ϕ}{\partial x}=-k^2x+\dfrac{k^3}{4}x^2-\dfrac{ik^3}{4}x^2 + A_1 $ donc pour $ x=L : \dfrac{\partial ϕ}{\partial x}= -k+\dfrac{k}{4}-\dfrac{ik}{4}+A_1~~(2) $
Donc $ (1)=(2) \Rightarrow -\dfrac{17}{12}ik+\dfrac{k}{12}+A_1L=-k+\dfrac{k}{4}-\dfrac{ik}{4}+A_1~ie~A_1=k-\dfrac{1}{6}iK $


D'où $ \phi_1(x)= -1-\dfrac{k^2x^2}{2}+\dfrac{(k^3-ik^3)}{12}x^3+(k-\dfrac{ik}{6})x $
On continue pendant un nombre n d'itérations. On pourra ensuite calculer $ \boxed{h(x,t)=\Re((\phi_0+...+\phi_n)e^{-i\omega t})} $
Nous faisons ce calcul sur Maxima et nous obtenons les courbes suivantes (mises en regard de la solution analytique pour $ k=1 $ et $ \omega=1 $, à $ t=0 $) :


Cas1 homotopie ordres0à10 sans module.gif



Etude de sensibilité des paramètres

Il y a deux paramètres qui interviennent dans notre équation : la fréquence $ \omega $ et la profondeur $ H $. On peut voir sur le graphique ci-dessous l'impact de la variation de la fréquence $ \omega $ sur la forme de la houle.

Houle cas1 TriplePériode.gif

En effet, on remarque que plus la période est faible, plus les vagues sont rapprochées. Une brève discussion sur la période d'une houle est faite en fin d'article.

Sur le prochain graphique, on observe l'impact de la profondeur $ H $ sur la houle.

Houle cas1 TripleHauteur.gif

On remarque que plus la profondeur est faible plus les vagues se rapprochent et plus elles ralentissent, à la différence d'un simple changement de période. L'amplitude est supposée constante ici, mais ce n'est pas le cas dans la réalité : à faible profondeur non seulement les vagues se rapprochent et ralentissent mais aussi leur hauteur augmente.

Cas n°2 : canal uniforme unidimensionnel avec réflexion totale en amont

Nous étudions maintenant un canal uniforme unidimensionnel, avec une réflexion totale en amont $ ϕ_x=0 $ et une condition de flux aval $ ϕ_x=ik(2−ϕ) $

Comme dans le cas n°1, l'équation de Berkhoff se simplifie comme suivant: $ \frac{\partial^2 ϕ}{\partial x^2} +k^2ϕ = 0 $ Seules les conditions limites changent.

Solution analytique

L'équation caractéristique est la même que dans le cas n°1 : $ X^2 + k^2 = 0 $
La solution est de la forme $ Ae^{-ikx}+ Be^{ikx} $, avec A et B des constantes.

Détermination des constantes $ A $ et $ B $ grâce aux conditions limites :

  • $ A=e^{2ikL} $
  • $ B=1 $

Donc $ ϕ(x)=e^{ik(2L-x)}+e^{ikx} $

La solution est simplifiée à une amplitude et une phase nulles. Cependant on constate qu'ajouter une amplitude $ U $ et une phase $ \Phi $ ne change pas cette solution.

Ainsi $ ϕ(x,t)=e^{i[k(2L-x)-\omega t]}+e^{i(kx-\omega t)} $

D'où la hauteur d'eau est : $ \boxed{h(x,t)=\Re(ϕ(x,t))=cos(k(2L-x)-\omega t)+cos(kx-\omega t)} $


Houle3.gif

Solution semi-analytique (homotopie)

Comme dans le cas n°1, la relation d'homotopie s'exprime comme suivant : $ (1−p)(ϕ_{xx}-u_{0,xx})+p(ϕ_{xx}+k^2ϕ)=0 $


Grâce au même raisonnement que ci-dessus, on obtient donc : $ (1−p)((ϕ_{0,xx}(x)+pϕ_{1,xx}(x)+p^2ϕ_{2,xx}(x)+p^3ϕ_{3,xx}(x)+...)-u_{0,xx})+p[ϕ_{0,xx}(x)+pϕ_{1,xx}(x)+p^2ϕ_{2,xx}(x)+p^3ϕ_{3,xx}(x)+...+k^2(ϕ_0(x)+pϕ_1(x)+p^2ϕ_2(x)+p^3ϕ_3(x)+...)]=0 $


On développe de même en fonction des ordres de p (seules les conditions limites changent).
Ordre 0 :
On a $ ϕ_{0}=Ax+B+u_0 $


Détermination des constantes $ A $ et $ B $ grâce aux conditions limites :

  • en aval : $ x=L $ $ \Rightarrow $ $ \dfrac{\partial ϕ_0}{\partial x} = 0 \rightarrow A=0 $
  • en amont : $ x=0 $ $ \Rightarrow $ $ \dfrac{\partial ϕ_0}{\partial x} = ik(2-\phi_0) = ik(2-Ax-B-u_0)=A~ie~ik(2-B-u_0)=0~ie~B=1 $

donc $ \boxed{ϕ_{0}=1+u_0=2} $


Ordre 1 :
On a $ ϕ_{1}=u_0-k^2\iint ϕ_0+A_1x+B_1 $
or $ \iint ϕ_0= x^2 $ donc $ ϕ_1(x) = u_0-k^2x^2+A_1x+B_1 $

or

  • en aval : $ x=L $$ \Rightarrow $ $ \dfrac{\partial ϕ_1}{\partial x} = 0 \Rightarrow A_1-2k^2L=0 \Rightarrow A_1=2k^2L \Rightarrow A_1=2k $ avec $ kL=1 $
  • en amont : $ x=0 $ $ \Rightarrow $ $ \dfrac{\partial ϕ_1}{\partial x} = ik(2-\phi_1) = ik(2-u_0-k^2x^2+2kx+B_1) = 2k-2k^2x $

Or $ x=0 $ donc $ ik(2-u_0+B_1)=2k = 2ik-iku_0+B_1k=2k $
$ \Rightarrow B=\dfrac{2k-2ik+ik}{ik}={2ik^2+k^2}{-k^2}~ie~B=2i+1 $

D'où $ \boxed{ϕ_1(x)=-k^2x^2+2kx+2i+2} $ avec $ u_0=1 $

On continue pendant un nombre n d'itérations. On pourra ensuite calculer $ \boxed{\Re((\phi_0+...+\phi_n)e^{-i\omega t})} $
Nous faisons ce calcul sur Maxima et nous obtenons les résultats suivants :

Cas2 gif homotopie.gif

N.B. La phase et l'amplitude de la solution analytique sont ajustées respectivement à -3,45 et 2,42 afin d'obtenir une bonne superposition des données. Les raisons de cet ajustement nous sont inconnues.

Cas n°3 : canal uniforme unidimensionnel avec pente de fond constante et sortie libre en amont

Nous cherchons à modéliser un canal uniforme unidimensionnel avec une pente de fond constante avec l'entrée d'une onde de fréquence unitaire par l'aval $ ϕ=1 $ et une sortie libre en amont $ ϕ_x=ikϕ $.

Dans ce cas précis, l'équation de Berkhoff se simplifie comme suit: $ \nabla (CC_g\nabla\phi)+k^2CC_g\phi=0 $
or $ CC_g=gH $ et $ H $ n'est pas constante car il s'agit d'une pente.

On prend donc $ H(x)=-ax+h $ avec $ h $ la hauteur initiale et $ a $ la pente.

D'où $ H(x)\phi_{xx}(x)+H'(x)\phi_x(x)+k^2H(x)\phi(x)=0 $.


Solution analytique

On cherche à mettre cette équation sous forme d'équation de Bessel de type : $ x^2 \displaystyle\frac{\partial^2 ϕ}{\partial x^2} +x\frac{\partial ϕ}{\partial x} + x^2ϕ = 0 $
On prend comme changement de variable $ X=\dfrac{K}{P}(h-px) $
On a alors : $ pkX\dfrac{\partial^2\phi}{\partial X^2}+pK\dfrac{\partial\phi}{\partial X}+pkX\phi=0 $
En multipliant par X et en divisant par pK, on a $ X^2\phi_{XX}+X\phi_X+X^2\phi=0 $

On a donc la solution générale $ \phi(x)=A.J_0(x)+B.Y_0(x) $ avec $ J_0 $ une fonction de Bessel de 1ère espèce et $ Y_0 $ une fonction de Bessel de 2ème espèce. D'après les conditions limites : $ \phi(0)=1 $ et $ \phi(x)=ik\phi $

On trouve :$ \phi(x) = \frac{J_0(x)*(Y_1(x_L)-iY_0(x_L)) + Y_0(x)*(J_1(x_L)+iJ_0(x_L))}{J_0(x_0)*(Y_1(x_L)-iY_0(x_L))+Y_0(x_0)*(J_1(x_L)+iJ_0(x_L))} $ , soit, avec $ \widetilde{Y} = Y_1(x_L)-iY_0(x_L) $ et $ \widetilde{J}=J_1(x_L)+iJ_0(x_L) $ : $ \phi(x) = \frac{J_0(x)*\widetilde{Y} + Y_0(x)*\widetilde{J}}{J_0(x_0)*\widetilde{Y}+Y_0(x_0)*\widetilde{J}} $

Cas3 analytique4.png

N.B. La pente du fond n'est pas représentée ici.

Solution semi-analytique (homotopie)

On a donc $ H(x)\phi_{xx}(x)+H'(x)\phi_x(x)+k^2H(x)\phi(x)=0 $ avec $ H=-ax+h $
D'où $ (-ax+h)\phi_{xx}(x)-a\phi_x(x)+k^2(-ax+h)\phi(x)=0 $


On a donc la relation d'homotopie : $ (1-p)(\phi_{xx}-u_{0,xx})+p((h-ax)\phi_{xx}-a\phi_x(x)+k^2(h-ax)\phi(x)=0 $


On prend $ u_0=0 $ pour simplifier les calculs. On obtient donc :
$ (1-p)(\phi_{xx}+p((h-ax)\phi{xx}-a\phi_x(x)+k^2(h-ax)\phi(x)=0 $


Ordre 0 :

  • $ \phi_{0,xx}=0 \Rightarrow \dfrac{\partial ϕ_0}{\partial x} = A \Rightarrow \phi_0=Ax+B $
  • $ \phi_0(0)=B=1 \Rightarrow B=1 \Rightarrow \phi_0=Ax+B $
  • $ \dfrac{\partial ϕ_0}{\partial x}(L) = ik(AL+B)=A \Rightarrow ik(AL+1)=A \Rightarrow A=\dfrac{ik}{1-ikL} $

D'où : $ \boxed{\phi_0(x)=\dfrac{ikx}{1-ikL}x+1} $


Ordre 1 :

  • $ \phi_{1,xx}+(h-ax-1)\phi{0,xx}-a\phi_{0,x}+k^2\phi_0=0 $
  • $ \dfrac{\partial ϕ_1^2}{\partial x^2}= -(h-ax-1)\dfrac{\partial ϕ_0^2}{\partial x^2}+a\dfrac{\partial ϕ_0}{\partial x}-k^2\phi_0 $
  • $ \dfrac{\partial ϕ_1}{\partial x}= -(h-ax-1)\dfrac{\partial ϕ_0}{\partial x}+a ϕ_0-\int{k^2}\phi_0+A $
  • $ ϕ_1= -(h-ax-1)ϕ_0+a \int{ϕ_0}-k^2\iint{\phi_0}+Ax+B $

or $ ϕ_0=\dfrac{ikx}{1-ikL}x+1 ~ie~ \int{ϕ_0}=\dfrac{ikx^2}{2(1-ikL)}x^2+x ~ie~ \iint{ϕ_0}=\dfrac{ikx^3}{6(1-ikL)}x^3+\dfrac{x^2}{2} $

D'où : $ \boxed{\phi_1(x)=(ax+1-h)\phi_0+a(\dfrac{ikx^2}{2(1-ikL)}x^2+x)-k^2(\dfrac{ikx^3}{6(1-ikL)}x^3+\dfrac{x^2}{2})+Ax+B} $

Les graphes donnés par Maxima sont déconcertants (pour l'instant):

Graphe avec prise en compte du module de $ \phi $ La solution homotopique n'est même pas calculée.

Cas3 loupé.png


Graphe sans prise en compte du module de $ \phi $

Cas3 sans module.png

Cas n°5 : domaine bidimensionnel plat avec direction préférentielle des houles

Le dernier cas concret que nous traitons est un domaine bidimensionnel plat avec une direction préférentielle des houles.


L'équation de Berkhoff se simplifie comme suivant pour un canal bidimensionnel : $ \nabla (gH\nabla\phi)+k^2gH\phi=0 $
Dans ce cas précis, le fond est plat, d'où : $ \nabla \phi+k^2\phi=0 $
Ici, on ne peut pas trouver la solution analytique. On peut seulement approcher la solution grâce à une méthode de différences finies ou grâce à l'homotopie. Nous allons ici résoudre cette équation par homotopie (solution semi-analytique).

Solution semi-analytique (homotopie)

On a $ \nabla \phi(x,y)+k^2\phi(x,y)=0 ~ie~ \dfrac{\partial^2\phi}{\partial x^2}+\dfrac{\partial^2\phi}{\partial y^2}+k^2\phi(x,y)=0 $
Donc on obtient la relation homotopique suivante : $ (1-p)(\phi_{xx}+\phi_{yy})+p(\phi_{xx}+\phi_{yy}+k^2\phi)=0 $


Une idée de développement peut être de considérer l'équation de Berkhoff, simplifiée en Helmholtz :


Cas5 biais.jpg

Cas5 biais.jpg

Et de faire un changement de repère de façon à obtenir :


Cas5 droit.jpg

Les limites du modèle et applications

Nous rencontrons 2 grands types de limitations : les limitations physiques ainsi que les limitations mathématiques

Limites physiques

  • Le fond du canal influence la houle et réciproquement. Or, notre modèle ne prend notamment pas en compte des phénomènes tels que le transport des sédiments et du sable qui pourraient modifier le fond du canal et donc avoir une influence sur la houle. Ces phénomènes sont cependant limités car se déroulant sur un temps long.
  • Le fond est en général très irrégulier, et la présence de hauts-fonds ou d'autres irrégularités est de nature à modifier la houle ; c'est d'ailleurs ce qui est en général utilisé pour limiter l'érosion, en installant des épaves de bateau en guise de digue par exemple.

Epaves Australie Moreton Island.jpeg (Épaves protectrices - Moreton Island, Australie - Photo Ludovic BRETON - travail personnel)

  • Le trait de côte est très irrégulier lui aussi, et les phénomènes de réflexion et diffraction sont nombreux. Modéliser la houle en simplifiant le trait de côte peut ainsi induire des erreurs.
  • L'onde est assimilée à une onde sinusoïdale pure. Bien que la houle a tendance à devenir une onde pure au bout d'une certaine distance (en ne conservant que les ondes les plus longues), celle-ci n'est pas forcément totalement pure, et notre modèle ne le prend pas en compte.
  • Nous pouvons questionner le réalisme de certaines données (kL = 1 par exemple).
  • Les vagues deviennent de la houle à partir du moment où elles quittent leur aire génératrice (en général loin au large), c'est-à-dire où elles sortent de la zone où souffle le vent qui les a fait naître. Il y a donc très souvent une différence entre "vagues de mer" et houle au niveau de leurs direction, célérité et amplitude. Ce modèle n'est pas capable de prendre en compte des zones où la mer est fréquemment sujette à des "croisements" entre vagues et houle.

Ile de ré Mer croisée Michel Griffon.JPG (Mer croisée à l'île de Ré - Photo Michel Griffon - Source wikipedia)


Enfin, considérer la hauteur d'eau comme un proxy pour la quantification de l'érosion est réducteur, puisque seule l'énergie potentielle des vagues est considérée. Leur énergie cinétique n'est pas prise en compte ; or on peut raisonnablement penser que celle-ci joue fortement sur l'érosion du trait de côte, ainsi que sur les submersions marines.

Limites mathématiques

  • Beaucoup de variables sont considérées constantes ($ k,\omega,L $etc), afin de limiter la complexité du modèle.
  • on suppose une uniformité spatiale et temporelle du phénomène : La zone génératrice de la houle est supposée être à un seul endroit et constante. On peut ainsi prédire l'effet d'une houle moyenne (en ondes longues, donc à faible profondeur), mais le pouvoir érosif de grandes crêtes n'est pas calculé.
  • l'ordre d'homotopie "n" influe beaucoup sur la qualité du modèle : La divergence est rapide si "n" est petit.

Propositions d'applications

Cette méthode permet de modéliser à la fois la réflexion des ondes, ainsi que leur réfraction ou diffraction sur les digues. Elle permet ainsi de prédire non seulement les phénomènes de submersion marine, mais aussi et surtout l’impact de l’érosion sur le trait de côte.
De plus, la méthode de l'homotopie méthode s’adapte à une grande variété de géométries et donc une grande variété de cas. Cela peut être très utile pour modéliser des formes de littoral complexes avec des estuaires, des rades ou des baies par exemple afin de piloter leur aménagement par la suite.

Outils personnels