-- Simulation numérique en mécanique des fluides --

TP n°3

 Représenter un champ de données fluide

simon.marie@lecnam.net
In [3]:
import numpy as np
import matplotlib.pyplot as plt
# la librairie time nous sera utile pour calculer le temps d'execution des scripts
import time 
#Option pour afficher les figures dans le notebook et eviter le plt.show():
%matplotlib inline   

On se propose dans ce TP, de représenter et d'analyser des données issues d'une simulation numérique.

Pour cela on dispose de 3 fichiers:

  • CM_rho.dat qui contient les données de pression
  • CM_ux.dat qui contient les données de vitesse longitudinal
  • CM_uy.dat qui contient les données de vitesse transversal

Ces données sont issues d'une simulation numérique de l'écoulement autour d'un cylindre carré de coté $D$ réalisée sur un maillage régulier de $600\times600$ mailles. Les carrée est au centre du domaine et son coté est représenté à l'aide de 5 maille.($D=5$).

Les objectifs de ce TP sont:

  • Représenter des données en 2D
  • Extraires des informations de ces données
  • Calculer d'autres grandeurs issues de ces données.

1 - Chargement des données

D'abord, afin de pouvoir utiliser les données et les représenter, il faut les charger en mémoire c'est à dire lire le fichier de données et les stocker dans un tableau. Pour cela il existe de nombreux moyen en Python. Pour ce TP nous utiliserons la fonction fromfile de numpy qui permet de lire des données binaires.

In [5]:
nx,ny=600,600
rho=np.fromfile('../CM_Rho.dat').reshape(nx,ny)
ux=np.fromfile('../CM_UX.dat').reshape(nx,ny)
uy=np.fromfile('../CM_UY.dat').reshape(nx,ny)

## Ici la fonction reshape (reformeze) permet de redimensionner le vecteur ux en une matrice 600x600
## En effet la fonction fromfile renvoi un vecteur et non une matrice

Les données issues de simulations sont souvent normalisées par des grandeurs de références. Ici la masse volumique et la vitesse ont été normalisé par les grandeurs de références suivantes:

  • $\rho_0=1.22$ $kg.m^{-3}$
  • $U_0=588.89$ $m/s$

A partir des données chargées, retouver les valeurs min, max et moyenne pour la vitesse longitudinale en unité physiques.

In [6]:
rho0=1.22
U0=588.89
rho=rho0*rho
ux=U0*ux
uy=U0*uy
In [7]:
print(ux.min(),ux.max(),ux.mean())
-16.2694990623 128.459282087 101.084098258

Ainsi, la vitesse max est de 128 m/s et la vitesse min de -16 m/s c'est à dire 16 m/s dans le sens des $x$ négatifs.

2 - Visualisation des données

Maintenant que les données sont en mémoire, nous pouvons les "visualiser" en utilisant la fonction pcolormesh de matplotlib.

In [8]:
plt.figure(figsize=(12,6),dpi=92)
plt.pcolormesh(rho.T,cmap="jet")
plt.colorbar()
plt.axis('equal')
plt.xlim([250,450])
plt.ylim([200,400]) # On peut "zoomer " sur la partie qui nous intéresse
Out[8]:
(200, 400)
In [9]:
plt.figure(figsize=(12,6),dpi=92)
plt.pcolormesh(ux.T,cmap="jet")
plt.colorbar()
plt.axis('equal')
plt.xlim([250,450])
plt.ylim([200,400]) # On peut "zoomer" sur la partie qui nous intéresse
Out[9]:
(200, 400)

On peut également utiliser la fonction contour pour visualiser les lignes de même vitesse (lignes de potentiel) en précisant le nombre de contour que l'on souhaite (ici 25).

In [14]:
plt.figure(figsize=(12,6),dpi=92)
plt.contour(ux.T,25)
plt.axis('equal')
plt.xlim([250,400])
plt.ylim([250,350]) # On peut "zoomer" sur la partie qui nous intéresse
Out[14]:
(250, 350)

La fonction streamplot trace les lignes de courant permettant de visualiser la trajectoire d'une particule de fluide.

In [11]:
x=np.arange(nx)
y=np.arange(ny)
plt.figure(figsize=(12,6),dpi=92)
plt.streamplot(x,y,ux.T,uy.T,7,color=rho.T,cmap="jet")
plt.axis('equal')
plt.xlim([250,400])
plt.ylim([250,350]) # On peut "zoomer " sur la partie qui nous intéresse
Out[11]:
(250, 350)

Commentaires:

Il peut également être utile de représenter une certaine partie des données.

A partir des données chargées, tracer le profil de pression à $x=320$ et le profil de vitesse à $y=310$:

In [15]:
plt.plot(rho[320,:])
Out[15]:
[<matplotlib.lines.Line2D at 0x7ff4d8cd9630>]
In [16]:
plt.plot(ux[:,310])
Out[16]:
[<matplotlib.lines.Line2D at 0x7ff4d807a978>]

Commentaires: On observe la dépression derrière le carré et l'accélération sur le haut du carré.

3 - Extraction d'informations

Les données fournies ($\rho$, $Ux$, $Uy$) sont les solutions de la simulation. Il peut être utile d'en extraire des informations.

En vous aidant des résultats précédents ainsi que de vos connaissances en mécanique des fluides, retrouver les infrmations suivante en partant des données fournies uniquement:

  • Pression dynamique
  • Pression total
  • Coefficient de trainée $C_x$ (On doit trouver pour un carré $C_x\sim 1$
In [54]:
Pdyn=0.5*rho.T*(ux**2+uy**2) # Pression Dynamique
plt.figure(figsize=(12,6),dpi=92)
plt.pcolormesh(Pdyn.T,cmap="jet")
plt.colorbar()
plt.axis('equal')
plt.xlim([250,450])
plt.ylim([200,400])
Out[54]:
(200, 400)
In [17]:
plt.figure(figsize=(12,6),dpi=92)
plt.pcolormesh(340**2*rho.T+0.5*rho.T*(ux**2+uy**2).T,cmap="jet") # Pression totale
plt.colorbar()
plt.axis('equal')
plt.xlim([250,450])
plt.ylim([200,400])
Out[17]:
(200, 400)

La pression totale se conserve le long d'une ligne de courant dans le cas non-visqueux. Ici, nous avons visiblement un écoulement visqueux engendrant un écoulement décollé derrière le carré. Nous reviendront largement sur ces points physique dans la suite du cours.

Pour le coefficient de trainée, la force étant définie en fonction de la pression dynamique:

$$ F_D=\dfrac{1}{2}\rho U_{\infty}^2 SC_D $$

et étant reliée au diférentiel de pression entre les 2 face du carrée:

$$ F_D=S\Delta p $$

On peut écrire:

$$ C_D=\dfrac{\Delta p}{P_{dyn}} $$

Si on trace les deux profils de pression de part et d'autres du carré ($x=296$ et $x=302$) on observe:

In [56]:
x1=297
x2=301
plt.plot(rho[x1,:],'r',label='x='+str(x1))
plt.plot(rho[x2,:],'b',label='x='+str(x2))
plt.legend()
plt.xlim([290,310])
Out[56]:
(290, 310)

Ainsi on peut calculer le différentiel de pression et en déduire la coefficient de trainée du Carré.

In [57]:
P1=np.mean(rho[x1,x1:x2+1])
P2=np.mean(rho[x2,x1:x2+1])
Delta=340**2*(P1-P2)
Cd=Delta/Pdyn[100,100]
print(Cd)
1.07654067276

Conclusion

Dans ce TP nous nous sommes entrainé à charger puis à visualiser des données issues d'un calcul et à analyser les résultat en traçant d'autres grandeurs dérivant des grandeurs enregistrées. (Pdyn, Ptot, Cp, Cx...)

In [1]:
from IPython.core.display import HTML
style=open('../notebooks.css', "r").read()
HTML(style)
Out[1]: