Inversion

Principe

Le principe de l’inversion est de déterminer un champ de perméabilité tel qu’un calcul d’écoulement restitue les données piézométriques fournies en entrée de la procédure.

Ces données piézométriques peuvent être :

  • des valeurs fournies en des points donnés, dénommés “points pilotes”,

  • des valeurs fournies sur l’ensemble du domaine, au travers d’une carte piézométrique, dénommée “potentiel de référence”. Cette carte est généralement obtenue par des méthodes géostatistiques à partir de mesures ponctuelles.

La procédure d’inversion intégrée dans le plugin est limitée au régime permanent (i.e. infiltration constante) et fondée sur la méthode des points pilotes [RJ08], [RT16].

Le principe est d’ajuster les perméabilités de plusieurs “zones pilotes”, par un processus itératif visant à minimiser un critère d’erreur mesurant l’écart des valeurs piézométriques calculées à celles mesurées :

\[ e_{i}=\frac{1}{n}\sum_{k=1}^{n}\left|h_{k}^{i}-h_{k}^{m}\right| \]

où :

  • \(n\) désigne le nombre total de points pilotes,

  • \(h_{k}^{i}\) désigne le niveau piézométrique calculé au point \(k\) à l’itération \(i\),

  • \(h_{k}^{m}\) désigne le niveau piézométrique mesuré au point \(k\).

À chaque itération, les perméabilités de chaque zone pilote sont ajustées de la manière suivante :

\[ K_{j}^{i+1}={\displaystyle (1+\alpha \frac{\delta_{j}^{i}}{\delta_{max}^{i}}) K_{j}^{i}} \]

avec :

\[ \delta_{j}^{i}=\frac{1}{n_{j}}\sum_{k=1}^{n_{j}}\left(h_{k}^{i}-h_{k}^{m}\right) \]

où :

  • \(\alpha\) est compris strictement entre 0 et 1,

  • \(K_{j}^{i}\) désigne la perméabilité de la zone \(j\) à l’itération \(i\),

  • \(n_{j}\) désigne le nombre de points pilotes de la zone \(j\).

  • \(\delta_{max}^{i}\) désigne la valeur maximale des \(|\delta_{j}^{i}|\) à l’itération \(i\).

Le paramètre \(\alpha\) est initialisé à 0.5, valeur qui semble un bon compromis entre la précision et la vitesse de convergence. Lorsque l’erreur augmente dans une proportion trop importante, ce paramètre \(\alpha\) est diminué d’un facteur deux et une nouvelle itération reprend à partir de l’avant dernière itération.

Pour chaque zone de perméabilité, la perméabilité est contrainte par des valeurs minimale et maximale :

\[ \forall i,j\;\;\;K_{j}^{min}\leq K_{j}^{i}\leq K_{j}^{max} \]

Le processus itératif est stoppé dans plusieurs cas :

  • lorsque le nombre maximal d’itérations est atteint,

  • lorsque l’erreur devient inférieure à une valeur limite \(e_{max}\) : \(e_{i}\leq e_{max}\),

  • lorsque le paramètre \(\alpha\) devient inférieur à une valeur minimale : \(\alpha\leq \alpha_{min}\),

  • lorsque le nombre de tentatives de diminution du paramètre \(\alpha\) atteint une valeur maximale.

Implementation

Commencer par spécifier les points pilotes dans la couche points_pilote, en indiquant obligatoirement leur zone d’appartenance et le niveau piézométrique mesuré (altitude_piezo_mesuree). Les zones doivent être numérotées à partir de 1 et sans trou. Plusieurs points pilotes peuvent appartenir à une même zone.

Il est possible d’importer des fichiers texte (menu QGIS Couche>Ajouter une couche>Ajouter une couche de texte délimité) et de procéder par copier-coller entre la couche importée et la couche points_pilote. Afin de faciliter le copier-coller, il est conseillé de nommer les colonnes du fichier importé selon les noms utilisés dans la table points_pilote (nom, zone, altitude_piezo_mesuree).

../../_images/inversion_pilotPointsTable.png

Figure 164 Spécification des points pilotes

Puis spécifier les valeurs des différents paramètres de la table hynverse_parametres (inversion_parametersTable) :

  • infiltration : valeur de l’infiltration moyenne annuelle en m/s. Si l’infiltration est hétérogène, indiquer une valeur nulle et fournir un fichier [u_infiltration.node.sat or u_infiltration.elem.sat].

  • icalc : type de calcul (3 par défaut),

  • niter : maximum iteration number (default value : 20, but higher values are useful to obtain a good convergence),

  • errlim : valeur objectif du critère d’erreur : le calcul s’interrompt lorsque l’erreur est inférieure à errlim (0.1 m par défaut),

  • alfmin : valeur minimale du paramètre \(\alpha\) : le calcul s’interrompt lorsque \(\alpha\) devient inférieur à alfmin (\(10^{-5}\) par défaut),

  • alfa : valeur initiale du paramètre \(\alpha\) (0.5 par défaut),

  • terr : valeur initiale du taux d’accroissement maximal admissible de l’erreur (1.2 par défaut) : si le taux d’accroissement de l’erreur est supérieur à ce taux, la dernière itération n’est pas prise en compte et une nouvelle itération est réalisée, en diminuant le paramètre \(\alpha\) d’un facteur 2. Cette procédure peut être renouvelée un nombre de fois inférieur ou égal à nbessaismax (voir ci-dessous),

  • nessaismax : nombre maximal d’itérations possibles en réduisant le paramètre \(\alpha\) d’un facteur 2(10 par défaut),

  • permini : valeur initiale de la perméabilité (\(2.10^{-5}\) par défaut),

  • permin : valeur minimale de la perméabilité (\(2.10^{-7}\) par défaut),

  • permax : valeur maximale de la perméabilité (\(2.10^{-4}\) par défaut),

  • nv_npp : nombre de noeuds les plus proches des points pilotes, pour le calcul du potentiel aux points pilotes (4 par défaut),

  • nv_ppm : nombre de points pilotes les plus proches des centres de mailles, pour le calcul du champ de perméabilité aux mailles (10 par défaut),

  • dv_pp : rayon (en mètres) autour des points pilotes pour lisser la perméabilité par moyenne mobile (300 m par défaut),

  • d_mesh : distance minimale d’un centre de maille à un noeud à potentiel imposé, pour ajouter le centre de maille comme point pilote (dans le cas d’un potentiel de référence externe) (100 m par défaut).

../../_images/inversion_parametersTable.png

Figure 165 Table des paramètres d’inversion

L’inversion est déclenchée en cliquant sur le bouton ../../_images/inverse.png, qui ouvre le dialogue de la Figure 75.

Sélectionner le premier item pour lancer l’inversion.

L’inversion fait apparaître une fenêtre ancrable continuant deux courbes (inversion_canvas) :

  • la première montre l’évolution de l’erreur à chaque itération,

  • la seconde montre l’évolution de la corrélation entre les valeurs piézométriques mesurées et calculées aux différents points pilotes, organisés par groupes.

La carte de perméabilités est également mise à jour à chaque itération, et s’affiche sous réserve qu’elle a été sélectionnée avant le début du calcul.

../../_images/inversion_canvas.png

Figure 166 Fenêtres QGIS au cours de l’inversion

À noter que, dans la (inversion_canvas), les noeuds à potentiel imposé sont mis en évidence. Pour cela, il suffit de faire appel au Gestionnaire de base de données (menu QGIS Database>DBManager…), de sélectionner la base de calcul courante et d’effectuer la requête SQL :

select * from noeuds, potentiel_impose where noeuds.OGC_FID=potentiel_impose.nid

Il faut également sélectionner Charger en tant que nouvelle couche, en spécifiant GEOMETRY pour la colonne de géométrie, et attribuer un nom à la couche :

../../_images/noeuds_potentiel_impose.png