Penalisation active

Equation de la chaleur 2D

Le champ de temperature initial est le dipole comme dans la these de Keetels, mais les parois sont plus proches (si la taille du domaine est 1, les parois sont a 0.2 et 0.8). De cette facon on a la place de bien lisser le masque (methode basee sur une fonction de Bessel) a une resolution faible N=128, et les calculs sont rapides.
Le schema temporel est RK3 explicite. Le terme de controle est formule comme un simple produit alpha(y,t)*chi(x), ou chi est le masque et alpha est le parametre de controle, discretise suivant y avec 64 points, avec une interpolation lineaire pour les points manquants.
La fonctionelle a minimiser est estimee par une interpolation en spline bicubiques (mais si j'ai bien compris il faudrait changer en une simple interpolation polynomiale), et la derivee temporelle de alpha est estimee par une difference finie d'ordre 1. Normalement dans le schema RK3 il faudrait utiliser les valeurs de alpha aux temps intermediaires mais on fait l'hypothese (fausse) que alpha ne varie pas au milieu d'un pas de temps.
Le coefficient de diffusion est fixe a 1.
Apres essais, on trouve que les pas de temps inferieurs a 0.05 sont instables, ce qui est consistent avec le temps de propagation de l'information du masque lisse a la paroi.
Voici un film avec dt=0.06 qui montre que le controle a l'air de bien se passer mais que la solution n'est plus reguliere en temps dans la paroi.
Pour obtenir un meilleur resultat, on prend dt=0.03 pour avancer en temps mais on garde dt=0.06 lorsqu'on prevoit la valeur a la paroi (sinon c'est instable).
On obtient effectivement quelque chose de regulier en temps. Ce n'est pas clair pour le moment quelle solution est la meilleure, il faudra comparer avec la solution exacte.

Equation de Stokes 2D

On peut faire exactement la meme etude que ci-dessus avec l'equation de Stokes. (En fait on prend Navier-Stokes avec Re ~ 0.01) Elle n'est pas equivalente a l'equation de la chaleur car les conditions aux limites de Dirichlet portent sur u et non sur omega. On applique rajoute donc les termes de penalisation dans les equations pour ux et uy et on prend le rotationnel pour avoir celle de omega. (On pourrait peut-etre penaliser directement l'equation de omega pour simplifier ?)

Les films obtenus sont la et la.

Equation de Navier-Stokes 2D (Decembre 2009)

Pour pouvoir traiter le cas Navier-Stokes, on utilise deux couches de controle dans la paroi, l'une tres proche de la paroi, et l'autre nettement plus a l'interieur. De facon intuitive, on peut dire que cela permet de controler la vitesse a la paroi par l'intermediaire du champ de pression qui sont non locales, puisque on ne peut plus la controler par les contraintes visqueuses qui sont a trop courte portee.
On utilise ici une resolution de 256x256, et une resolution verticale de 32 le long des parois (donc 128 points de controle en tout en comptant les deux couches de chaque cote). Le temps processeur pour chacun de ces calculs a ete de l'ordre de 23h (eh oui c'est beaucoup).
On voit sur les films de vitesse horizontale et verticale que la vitesse reste bien nulle a la paroi, par contre elle atteint des vitesse assez importantes dans la paroi.
Sur la simulation Re~130, on voit aussi qu'il y a un debut d'instabilite dans la paroi. Dans les deux simulations, la symetrie verticale est brisee aux temps longs, donc il y a quand meme eu une erreur quelque part que je n'arrive pas bien a comprendre.
Cas Re ~ 100
Vorticite
Vitesse horizontale
Vitesse verticale (attention le signe est inverse)
Pression
Cas Re ~ 130
Vorticite
Vitesse horizontale
Vitesse verticale (attention le signe est inverse)
Pression