Simulation de fluide, utilisation de la bibliothèque

Date de publication : 9/11/2006 , Date de mise à jour : 9/11/2006

Dans ce tutoriel, je présente la classe de simulation de fluide et son utilisation à travers des exemples et des dessins explicatifs.

I. Introduction

La bibliothèque permettant de simuler les écoulements de fluide a été réalisé en C++. Afin de visualiser le résultat, j'ai utilisé la bibliothèque graphique GLUT, mais vous pouvez très bien utiliser la classe Fluide avec une autre bibliothèque graphique. Elle peut avoir des applications graphiques intéressantes, notamment dans le domaine du jeu-vidéo.

Cette classe permet de réaliser plusieurs choses :

Simuler les vecteurs vitesse du fluide dans un domaine rectangulaire
L'ajout d'obstacle rectangulaire
L'ajout d'un fluide de couleur différent (problème de diffusion)
Ajout d'une force pour la simulation de fumée
Calcul de la pression
Détection des tourbillons
Ajout ou suppression des 4 murs autour du fluide (haut, gauche, droit, bas)

Nous allons voir chacune de ces fonctionnalités à travers des exemples.

II. Création d'un fluide

Dans un premier temps, nous allons créer un fluide rectangulaire. Cela se fait avec l'appel au constructeur de Fluide.

Fluide fluide(int taillex, int tailley,  float discretisation);

Le constructeur dispose de trois paramètres. La taille horizontale (taillex), la taille verticale (tailley) et le pas de discrétisation temporelle (discretisation). Typiquement, il faut prendre un pas de discrétisation temporelle petit (cela joue un peu sur la vitesse du fluide), par exemple, on peut prendre : 0.1f.

J'ai essayé au maximum d'optimisé les calculs d'évolution de fluide, il semble que l'ordinateur commence à peiner pour un fluide de taille supérieur à 200 * 200 sur un 1.5ghz. Mais cela peut souvent suffir à faire de belles animations.

Initiallement, lorsque le fluide est crée, il n'y a aucun contour, aucun obstacle et la vitesse du fluide est nulle en tout point.

La taille du fluide après création peut être obtenu par deux méthodes.

 
/*prototype de la méthode*/
int tailleX() const; /*retourne la taille horizontale*/
int tailleY() const; /*retourne la taille verticale*/

int n = fluide.tailleX();

III. Lecture et écriture directe de la vitesse

Nous allons voir comment lire et modifier la vitesse en un point du fluide.

La lecture s'effectue simplement par l'appel à une méthode.

 
/*prototype de la méthode*/
std::pair<float, float> vitesseLire(int i, int j) const;

/*exemple*/
std::pair<float float> vectvitesse;

vectvitesse = fluide.vitesseLire(1,2);
std::cout<<"Vitesse X au point (1,2) = "<<vectvitesse.first<<std::endl;
std::cout<<"Vitesse Y au point (1,2) = "<<vectvitesse.second<<std::endl;

for(int j = 1; j < fluide.tailleY(); j++)
  for(int i = 1; i < fluide.tailleX(); i++)
  {
    vectvitesse = fluide.vitesseLire(i,j);
    /*traitement*/
  }

il faut noter que i et j doivent être dans les bornes du fluide.

Nous allons à présent voir comment forcer la vitesse en un point.

 
/*prototype de la méthode*/
void vitesseEcrire(int i, int j, float ix, float jy);

/*exemple*/
fluide.vitesseEcrire(1,2, 0.1f, 0.2f);

Ici, nous avons forcer au point (1,2) le vecteur vitesse (0.1, 0.2).

IV. Evolution du champs de vecteur vitesse

Le champs de vitesse se laisse évoluer par l'appel à une simple méthode.

 
/*prototype de la méthode*/
void evoluerVitesse();

/*exemple*/
for(int i = 0; i<10; i++)
{
  fluide.evoluerVitesse();
}

Cette exemple fait évoluer pendant 10 tours de boucle la vitesse. En combinant les écritures de vitesse en plusieurs points et l'évolution du fluide. Nous pouvons lire l'ensemble du champs de vitesse et obtenir ce type de résultat :

V. Ajout d'obstacles et bords

La classe Fluide permet de gérer facilement les obstacles rectangulaires et le fait que le fluide soit entouré de mur ou non.

 
/*prototype des méthodes*/
void ajouterObstacle(int x, int y, int w, int h);

/*Exemple*/
fluide.ajouterObstacle(10,20, 10, 30);

/*prototype*/
    void setMurGauche();
    void setMurDroit();
    void setMurHaut();
    void setMurBas();
    void unsetMurGauche();
    void unsetMurDroit();
    void unsetMurHaut();
    void unsetMurBas();

/*exemple*/
fluide.setMurHaut();
fluide.setMurDroit();

L'ajout d'obstacle ajoute un obstacle dont la coordonnée inférieur gauche est en (10,20) et de taille horizontale 10 et de taille verticale 30. À noter qu'une fois un obstacle placé, il n'est pas possible de le retirer seul. Tous peuvent être retiré d'un coup à l'aide de la méthode obstacleDetruire() (voir à la fin)

Les méthodes suivantes permettent d'ajouter ou de supprimer un mur (haut, gauche, droite et haut). Initiallement, tous les murs sont absents.

VI. Détermination et lecture de la pression

La pression est un champs de scalaires (c'est à dire un tableau à deux dimensions d'élément de type float). La pression peut se calculer à un moment donné par l'appel à une méthode.

 
/*prototype des méthodes*/
void calculPression();
float pressionLire(int i, int j) const;

/*Exemple*/
fluide.calculPression();
float p = fluide.pressionLire(2,3);

fluide.evoluerVitesse();
fluide.calculPression();

/*maintenant : p != fluide.pressionLire(2,3)

Il faut bien noter qu'il n'est pas nécessaire de faire appel à evoluerVitesse si par exemple on souhaite connaître la pression que tous les 10 tours de boucle. Attention, le calcul de la pression est gourmand au niveau de la rapidité. La vitesse peut presque chuter par deux.

Ici, la pression négative est bleue et la pression positive est en rouge. En générale, la pression est très basse, j'ai multiplié par un facteur important pour obtenir une couleur plus nette.

VI. Détection des tourbillons

Connaissant les vecteurs vitesse d'un fluide, il est assez facile à détecter les tourbillons dans celui-ci. Afin de visualiser cela, j'ai ajouté deux méthodes supplémentaires à la classe pour permettre ceci. Ces méthodes fonctionnent de la même manière que pour la pression.

 
/*prototype des méthodes*/
float vortexLire(int i, int j) const;
void calculVortex();


/*Exemple*/
fluide.calculVortex();
float p = fluide.vortexLire(2,3);

De la même manière que pour la pression, il n'est pas utile d'appeler calculVortex à chaque fois si l'on ne souhaite pas connaître la valeur du champs à tous les tours de boucles.

Le signe du float retourné par vortexLire indique le sens de rotation du fluide (ainsi que la couleur dans ce cas là).

VIII. Diffusion

Une application intéressante et la résolution du problème de diffusion. Imaginez que vous placiez une goutte d'encre noir dans une bassine (il faut imaginer qu'elle soit en deux dimensions). On souhaiterait faire évoluer cette goutte au cours du temps et connaître sa concentration. J'ai mis à disposition 3 méthodes : une qui permet de déposer du fluide à un endroit, une qui permet de faire évoluer à travers l'autre et une qui permet de connaître la concentration de ce fluide en tout point.

 
/*prototype des méthodes*/

/*evolution du fluide*/
void evoluerDiffusion();

/*lecture de la concentration*/
float diffusionLire(int i, int j) const;

/*depot du fluide en un point*/
void diffusionEcrire(int i, int j, float valeur);


/*Exemple*/
fluide.diffusionEcrire(10,12, 200);

for(int i= 0; i<20; i++)
{
  fluide.evoluerVitesse();
  fluide.evoluerDiffusion();
}

float d = fluide.diffusionLire(10,20);

Ici, contrairement à la pression et au vortex, il est obligatoire d'appeler evoluerDiffusion à chaque tour de boucle si l'on souhaite faire évoluer l'encre noir. Je ne l'ai pas intégré à evoluerVitesse pour permettre des optimisations de vitesse si l'on ne souhaite par exemple uniquement travailler sur la vitesse.

IX. Point automatique de diffusion

Afin de simplifier, j'ai ajouté une méthode permettant de déposer un fluide à un endroit automatiquement à chaque tour de boucle (cela permettra notamment de gérer la fumée automatiquement).

 
/*prototype des méthodes*/
void diffusionAutomatique(int x, int y, float valeur);

fluide.diffusionAutomatique(10,20, 200);

Cela ajoute donc une source de diffusion constant de valeur 200 à cette position. Il est possible d'en ajouter plusieurs. Ainsi, il ne sera pas nécessaire d'ajouter l'appel à la fonction diffusionEcrire à chaque tour de boucle, elle sera automatiquement appelée pour toutes les sources que vous avez ajoutées.

X. Simulation de fumée

Une application intéressante des fluides et de permettre la simulation de fumée. Pour cela, j'ai ajouté une méthode qui permet de diriger la fumée (en général vers le haut) avec une certaine vitesse de montée. Pour connaître la concentration de fumée, il faudra modifier le champs de diffusion, via diffusionEcrire ou via diffusionAutomatique.

 
/*prototype des méthodes*/
void ecrireForceFumee(float x, float y);

fluide.ecrireForceFumee(0.001f, 0.0f);

Ici, nous avons crée une force de fumée dirigée vers le haut. Il faut en général des nombres assez faible pour que la fumée ne monte pas trop vite. Ainsi, on peut facilement simuler de la fumée en utilisant ce code.

 
/*prototype des méthodes*/
fluide.ecrireForceFumee(0.0f, 0.0001f);
fluide.diffusionAutomatique(30,10);

for(int i=0; i<20; i++)
{
  fluide.evoluerVitesse();
  fluide.evoluerDiffusion();
}

/*lecture de la concentration de la fumée*/
float d = fluide.diffusionLire(20,20);

Voici un exemple avec :

 
  Fluide fluide(114, 84,  dt);

 fluide.ecrireForceFumee(0.001f,0.001f);
 fluide.diffusionAutomatique(30,10, 200);

XI. Classe Fluide

J'ai écrit ici l'ensemble des méthodes publiques de la classe Fluide.

 
class Fluide
{

   public:

    Fluide(int n, int m,  float dtin);
    ~Fluide();

     int tailleX() const;
     int tailleY() const;

    void ajouterObstacle(int x, int y, int w, int h);

    void evoluerVitesse();
    void calculPression();
    void evoluerDiffusion();
    void calculVortex();

    void vitesseEcrire(int i, int j, float ix, float jy);
    std::pair<float, float> vitesseLire(int i, int j) const;

    float diffusionLire(int i, int j) const;
    float pressionLire(int i, int j) const;
    float vortexLire(int i, int j) const;
    void diffusionEcrire(int i, int j, float valeur);

    void diffusionReinitialiser(void);
    void vitesseReinitialiser();


    void setMurGauche();
    void setMurDroit();
    void setMurHaut();
    void setMurBas();
    void unsetMurGauche();
    void unsetMurDroit();
    void unsetMurHaut();
    void unsetMurBas();


   /*fumée*/
   void ecrireForceFumee(float x, float y);
   void diffusionAutomatique(int x, int y, float valeur);

   void diffusionAutomatiqueDetruire();

   void obstacleDetruire();
};

Je précise l'existence des deux méthodes diffusionReinitialiser et vitesseReinitialiser qui permettent de rendre nulle les champs de vitesses et de diffusion.

Les méthodes obstacleDetruire permet de détruire tous les obstacles et la méthode diffusionAutomatiqueDetruire permet de détruire la liste des diffusions automatiques.

Les sources présentées sur cette page sont libres de droits et vous pouvez les utiliser à votre convenance. Par contre, la page de présentation constitue une œuvre intellectuelle protégée par les droits d'auteur. Copyright © 2006 Florent HUMBERT. Aucune reproduction, même partielle, ne peut être faite de ce site ni de l'ensemble de son contenu : textes, documents, images, etc. sans l'autorisation expresse de l'auteur. Sinon vous encourez selon la loi jusqu'à trois ans de prison et jusqu'à 300 000 € de dommages et intérêts.