12.5. TP problème inverse avec des PINNs#

Marc Buffat Dpt mécanique, UCB Lyon 1

problème inverse

# initialisation: (ne pas re-executer sauf en cas de restart)
import numpy as np
import os, time
from IPython.display import display,Markdown
def printmd(text):
    display(Markdown(text))
from validation.validation import info_etudiant, bib_validation, test_function
from validation.valide_markdown import test_code
try: NUMERO_ETUDIANT
except NameError: NUMERO_ETUDIANT = None 
if type(NUMERO_ETUDIANT) is not int :
    printmd("**ERREUR:** numéro d'étudiant non spécifié!!!")
    NOM, PRENOM, NUMERO_ETUDIANT = info_etudiant()
    #raise AssertionError("NUMERO_ETUDIANT non défini")
# parametres spécifiques
_uid_    = NUMERO_ETUDIANT
np.random.seed(_uid_)
printmd("**Login étudiant {} {} uid={}**".format(NOM,PRENOM,_uid_))
bib_validation('cours','IntroIA')
from Brinkman import Brinkman
BM = Brinkman(_uid_)
# initialisation GPU (ne pas re-executer sauf en cas de restart)
from validation.libIA_GPU import Init_torchGPU
try: cuda_dev
except NameError: cuda_dev = Init_torchGPU(4, _uid_%2)
if cuda_dev is None: cuda_dev = "cpu" 
print("Cuda GPU device: ",cuda_dev)
_notebook_="TP4_PINNinvpble.ipynb"

ERREUR: numéro d’étudiant non spécifié!!!

Login étudiant Marc BUFFAT uid=137764122

Max threads : 4 / used threads 4
Attention: no GPU available! using CPU
Cuda GPU device:  cpu
/home/buffat/venvs/jupyter/lib/python3.10/site-packages/torch/cuda/__init__.py:107: UserWarning: CUDA initialization: The NVIDIA driver on your system is too old (found version 9010). Please update your GPU driver by downloading and installing a new version from the URL: http://www.nvidia.com/Download/index.aspx Alternatively, go to: https://pytorch.org to install a PyTorch version that has been compiled with your version of the CUDA driver. (Triggered internally at ../c10/cuda/CUDAFunctions.cpp:109.)
  return torch._C._cuda_getDeviceCount() > 0

12.5.1. Objectif#

L’objectif du TP est l’utilisation d’un modèle de PINN (Physical Induced Neural Network) pour résoudre un problème inverse en mécanique des fluides.

Pour des écoulements en milieu poreux, à l’aide de quelques points de mesure et d’un modèle, on veut déterminer la valeur d’un paramètre physique difficilement mesurable : la viscosité efficace \(\nu_e\).

12.5.1.1. modèle de Brinkman: Écoulement en milieu poreux#

milieu poreux

Dans les milieux poreux, on utilise classiquement la loi de Darcy (1856):

\[ \frac{\nu}{K} u = f \]

\(\nu\) est la viscosité du fluide, \(K\) sa perméabilité et \(f\) le terme de force extérieure (gradient de pression et force de gravité): $\( f = -\frac{1}{\rho} \frac{\partial p}{\partial x} + \vec{g}.\vec{e_x}\)$

Le modèle de Brinkman est une extension de cette loi de Darcy pour les milieux poreux qui tient compte du cisaillement visqueux lorsque la perméabilité est élevée et impose la condition de non-glissement sur les parois solides. Dans le cas d’un écoulement 2D plan stationnaire dans un milieu poreux entre 2 plaques distantes de H, l’équation de Brinkman pour la vitesse \(u(x)\) s’écrit:

\[ \frac{\nu}{K} u = f + \frac{\nu_e}{\epsilon} \frac{d^2 u}{d x^2} \hspace{1cm} x \in [0,H] \]

associée aux 2 conditions aux limites sur les parois: \(u(0) =u(H) =0\)

Dans le terme supplémentaire de Brinkman, \(\epsilon\) est la porosité du milieu poreux et \(\nu_e\) la viscosité effective. Cette viscosité \(\nu_e\) est reliée à la structure poreuse du milieu et est malheureusement très difficile à mesurer.

A partir de 5 points de mesure de vitesse, on va déterminer la « meilleur approximation » de la viscosité effective \(\nu_e\) en utilisant le modèle de Brinkman.

12.5.1.2. Principe: modèle inverse#

  1. On se donne les valeurs expérimentales aux points de mesure: \(xd_i\), \(ud_i\) pour \(i=1,Nd=5\)

  2. Pour une valeur de \(\nu_e\) fixé, on résout numériquement le problème de Brinkman pour obtenir une solution \(\hat{u}(x)\)

  3. On calcul l’erreur aux points de mesure: $\(Err = \sum_{i=1}^{Nd} \left | \hat{u}(xd_i) - ud_i \right |^2\)$

  4. On corrige la valeur de \(\nu_e\) pour minimiser cette erreur

  5. On itère en 2 jusqu’à une précision fixée.

La résolution numérique et le problème de minimisation inverse utilise une approche PINN

# bibliothèque 
import numpy as np
import matplotlib.pyplot as plt
import torch
import torch.nn as nn

print(torch.cuda.is_available())
device = torch.device(cuda_dev)
print("CUDA device: ",cuda_dev)
print("Torch CUDA device: ",device," threads:",torch.get_num_threads(),torch.get_num_interop_threads())
False
CUDA device:  cpu
Torch CUDA device:  cpu  threads: 4 4

12.5.2. Paramètres du problème#

Les paramètres du problème étudié sont les suivants:

printmd(f"**paramètres: H={BM.H}m  nu={BM.nu} m^2/s  eps={BM.eps}  K={BM.K} m^2  F={BM.F} m/s^2**")

paramètres: H=1.0m nu=0.001 m^2/s eps=0.4 K=0.001 m^2 F=1.0 m/s^2

12.5.3. Résolution de l’équation de Brinkman à nue fixé#

En choisissant \(\nu_e=\nu\) , on va résoudre le modèle de Brinkman avec une approche PINN en utilisant les valeurs des paramètres précédents.

Définir les paramètres du problème dans la cellule suivante en utilisant les noms de variables ci-dessous:

  • H nu eps K F nu_param

# paramètres
H = nu = eps = K = F = None
nue_param = None
### BEGIN SOLUTION
### END SOLUTION

12.5.3.1. création du réseau de neurones PINN#

Définir un réseau de neurones à 3 couches linéaires 1x32 32x32 32x1 avec une fonction d’activation en tanh. On créera une nouvelle classe PINN dérivée de la classe nn.Module.

Créer le modèle sur le GPU en utilisant cette nouvelle classe et le mettre dans la variable model

print("GPU ",device)
torch.set_default_device(device)
if torch.cuda.is_available(): torch.cuda.init()
NRN = 32
model = None
## BEGIN SOLUTION
## END SOLUTION
GPU  cpu

12.5.3.2. variables d’entrée#

on choisit comme entrée Nc = 100 points de collocation equi-répartis entre 0 et H à l’exclusion des bornes.

Créer le tenseur pytorch dans x _phys en indiquant que l’on veut calculer les gradients par rapport à ce tenseur (requires_grad) et sélectionner la bonne méthode d’optimisation dans la variable optimiser.

# points de collocation
Nc = 100
x_phys = None
optimiser = None
## BEGIN SOLUTION
## END SOLUTION

12.5.3.3. calcul du résidu#

Définir la fonction résidu(x,y,nue) qui calcule le résidu de l’équation pour une approximation y de la solution en un point x (qui sont tous les deux des tenseurs pytorch) pour une valeur de la viscosité effective nue

Vérifier en calculant le résidu à partir d’une prédiction du modèle

# résidu de l'équation
## BEGIN SOLUTION
## END SOLUTION

12.5.3.4. boucle de minimisation#

Ecrire la boucle de minimisation avec une variable itérative epoch pour minimiser la fonction perte et obtenir la solution approchée par PINN.

La fonction loss doit prendre en compte les conditions aux limites imposées;

# boucle de minimisation
loss = None
nue_param = nu
## BEGIN SOLUTION
## END SOLUTION

12.5.3.5. Comparaison avec la solution de référence#

Comparer la solution obtenue avec la solution de référence fournit par la fonction BM.Uex(x). On evaluera ces solutions sur 300 points equi-répartis, que l’on mettra dans la variable x_test . Pour ces points, on calcule la solution de référence dans y_exact et la solution prédite par PINN dans y_pred.

Tracer ces 2 solutions sur un même graphe avec des légendes pour comparer.

# comparaison solution exacte
x_test = None
y_exact = None
y_pred = None
## BEGIN SOLUTION
## END SOLUTION

12.5.4. Problème inverse#

Lire les points de mesure dans le fichier data.txt qui contient sur 2 colonnes x_data et y_data.

Tracer les mesures pour vérifier, puis les convertir en tenseur pytorch.

# lecture des données
x_data = y_data = None
## BEGIN SOLUTION
## END SOLUTION

12.5.4.1. Construction du model PINN#

De la même façon que précédemment, construire le modèle PINN mais en ajoutant un paramètre supplémentaire : nue_param

model = None
nue_param = None
optimizer = None
## BEGIN SOLUTION
## END SOLUTION

12.5.4.2. Boucle d’optimisation#

Ecrire la boucle d’optimisation en ajoutant dans la fonction perte l’écart avec les données mesurées.

# boucle de minimisation
loss = None
## BEGIN SOLUTION
## END SOLUTION

12.5.4.3. Analyse de la solution#

Analyser le résultat en traçant sur un même figure, la solution PINN et les données mesurées.

# comparaison donnee
x_test = y_pred = None
## BEGIN SOLUTION
## END SOLUTION 

12.5.5. Analyse et conclusion#

  • dans le fichier CompteRendu.md écrire votre analyse et vos conclusions en markdown

  • puis générer la version html avec la commande suivante

# génération de la version html du CR
!genereTPhtml CompteRendu

12.5.6. FIN#