# Inversion

## Operating principle

The goal of the inversion process is to determine a permeability field that, combined with a flow simulation, returns the measured piezometric data.

The piezometric data can be:

- values provided at ponctual location, called pilot points (points pilotes),

- values provided on the whole study area, through a piezometric layer, called reference potential (potentiel de référence). This layer usually is obtained by geostatistical methods from ponctual measures.

The inversion process provided by the plugin is limited to permanent regime (constant infiltration) and is based on the pilot points method [RJ08], [RT16].

The operating principle is to ajust permeability on several “pilots zones”, using an iterative computation to minimize the mean error between the piezometric simulated values and the measured values :

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

where :

> - \(n\) is the total pilot points number,

> - \(h_{k}^{i}\) is the piezometric value calculated at the location \(k\) for the iteration \(i\),

> - \(h_{k}^{m}\) is the piezometric value measured at the location \(k\).

For each iteration, each pilot zone permeability is ajusted as follows :

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

with :

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

where :

> - \(\alpha\) is strictly between 0 et 1,

> - \(K_{j}^{i}\) is the permeability of the zone \(j\) for the iteration \(i\),

> - \(n_{j}\) is the number of pilot points in the zone \(j\).

> - \(\delta_{max}^{i}\) is the maximal value of \(|\delta_{j}^{i}|\) for the iteration \(i\).

\(\alpha\) parameter is initiated to 0.5. The value seems to be a good compromise between precision and speed convergence. When the error increasing in a to important way, \(\alpha\) parameter is divided by two and a new iteration starts from the before last iteration.

For each permeability zone, permeability is constrained by min and max value :

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

Iterative process is stopped in several cases :

- when the maximal number of iteration is reached,

- when the calculated error is less than a limit value \(e_{max}\) : \(e_{i}\leq e_{max}\),

- when the \(\alpha\) parameter is less than a minimal value : \(\alpha\leq \alpha_{min}\),

- when the number of attempts to decrease the \(\alpha\) parameter has been reached.

## Implementation

Pilot points must first be defined in the points_pilote layer, and the zone where they belong and the measured piezometric value (altitude_piezo_mesuree) must be set. Zones must be numbered from 1 without empty number. Several pilot points can belong to the same zone.

It is possible to import delimited text files in QGIS (QGIS menu Layer>Add a layer>Add a delimited text layer) and do copy/paste actions between the imported layer and the points_pilote layer. To ease the copy/paste action, it is recommended to name the columns in the imported layer according to the points_pilote table columns (nom, zone, altitude_piezo_mesuree) (`inversion_pilotPointsTable`).

Figure 76 Pilot points table

Then the parameter values of the hynverse_parametres table must be set (`inversion_parametersTable`) :

- infiltration : yearly mean infiltration in m/s. For heterogeneous infiltration, please enter 0 value and provide a [u_infiltration.node.sat or u_infiltration.elem.sat] file.

- icalc : calcul type (default value : 3),

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

- errlim : error limit : process stops when the error is less than this value (default value : 0.1 m),

- alfmin : minimal value for \(\alpha\) : process stops when \(\alpha\) is less than this value (default value : \(10^{-5}\)),

- alfa : initial \(\alpha\) value (default value : 0.5),

- terr : initial value of the maximum permissible growth rate of error (default value : 1.2) : if the growth rate of error is larger, the last iteration is not taken into account and a new one is done by dividing \(\alpha\) by 2. This can be done nbessaismax times (see below),

- nessaismax : number of possible iterations by dividing \(\alpha\) by 2 (default value : 10),

- permini : initial permeability value (default value : \(2.10^{-5}\)),

- permin : minimal permeability value (default value : \(2.10^{-7}\)),

- permax : maximal permeability value (default value : \(2.10^{-4}\)),

- nv_npp : number of closest nodes from the pilot points, used for pilot points potential estimation (default value : 4),

- nv_ppm : number of closest pilot points from mesh centers, used for meshes permeability estimation (default value : 10),

- dv_pp : radius (m) around pilot points to smooth permeability field by moving average (default value : 300 m),

- d_mesh : min length from mesh center to node with imposed potential or to existing pilot point, used to add the mesh center as pilot point (in case of an external reference potential) (default value : 100 m).

Figure 77 Inversion parameters table

Inversion is triggered by clicking on button, which opens the dialog of Figure 163.

Select the first item to launch pilot points inversion.

Inversion process will display a new panel with two graphs (`inversion_canvas`) :

> - the first one displays error evolution at each iteration,

> - the second displays correlation evolution between simulated and measured piezometric values at each iteration, grouped by zone.

> The process updates the permeability map at each iteration and the layer is displayed if it has been activated before the calculation.

Figure 78 QGIS windows during an inversion process

Note that in the (`inversion_canvas`), nodes with imposed potential are highlighted. This can be obtained by using the DBManager (QGIS menu Database>DBManager…), selecting the current calculation database and executing the following SQL request :

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

When the result is loaded, check the Load as a new layer checkbox, choose the GEOMETRY item for the geometry column, and name the new layer :