Vizualisation of infiltration into a sand column
Water is infiltrated into a sand column and ERT datasets are collected throughout the infiltration. In this example we will invert them with ResIPy using a custom mesh (build withing gmsh UI).
import sys
sys.path.append('../src')
from resipy import Project
import numpy as np
import pyvista as pv
testdir = '../src/examples/'
API path = /media/jkl/data/phd/tmp/resipy/src/resipy
ResIPy version = 3.2.2
cR2.exe found and up to date.
R3t.exe found and up to date.
cR3t.exe found and up to date.
Let’s import the surveys, electrodes (with correct label
<string> <elecNumber>) and vizualise the mesh. The mesh actually has
two regions: the top region is the sand where the ERT measurements are
done. The bottom region is gravels which are just used for drainage.
Note that the last two rings of electrodes which are in the gravels are
not actually used in the sequence. However, for correct ERT modelling
the entire column needs to be modelled.
k = Project(typ='R3t')
k.createTimeLapseSurvey(testdir + 'dc-3d-timelapse-column/data', ftype='ProtocolDC')
k.importElec(testdir + 'dc-3d-timelapse-column/elec.csv')
k.importMesh(testdir + 'dc-3d-timelapse-column/mesh.msh')
k.showMesh()
Working directory is: /media/jkl/data/phd/tmp/resipy/src/resipy
clearing dirname
4/4 imported
Reading ../src/examples/dc-3d-timelapse-column/mesh.msh
Gmsh version == 3.x
reading node coordinates...
Determining element type...Tetrahedra
Reading connection matrix...
ignoring 4866 elements in the mesh file, as they are not required for R2/R3t
Finished reading .msh file
Converting np.character to a dtype is deprecated. The current result is np.dtype(np.str_) which is not strictly correct. Note that np.character is generally deprecated and 'S1' should be used. Converting np.character to a dtype is deprecated. The current result is np.dtype(np.str_) which is not strictly correct. Note that np.character is generally deprecated and 'S1' should be used.
Now invert the data. Note that to avoid automatic cropping of the mesh
to the region of interest (defined by the maximum extend of the
electrodes locations), we need to setup num_xy_poly to 0 and
zmin and zmax to \(-\infty\) and \(+\infty\)
respectively.
k.param['num_xy_poly'] = 0
k.param['zmin'] = -5
k.param['zmax'] = 1
k.invert(parallel=True)
Writing .in file and protocol.dat... Matching quadrupoles between surveys for difference inversion...100 in common...done in 0.012795s
done!
------------ INVERTING REFERENCE SURVEY ---------------
>> R 3 t E R T M o d e l v 2.30 <<
>> Date: 26-01-2021
>> My beautiful 3D survey
>> I n v e r s e S o l u t i o n S e l e c t e d <<
>> A d v a n c e d M e s h I n p u t <<
>> T e t r a h e d r a l E l e m e n t M e s h <<
>> Reading mesh file
>> Determining storage needed for finite element conductance matrix
>> Generating index array for finite element conductance matrix
>> Reading resistivity model from res0.dat
>> L o g - D a t a I n v e r s i o n <<
>> N o r m a l R e g u l a r i s a t i o n <<
>> Memory estimates:
For 1000 measurements the memory needed is: 0.083 Gb
For 2000 measurements the memory needed is: 0.162 Gb
For 5000 measurements the memory needed is: 0.399 Gb
For 10000 measurements the memory needed is: 0.794 Gb
>> Forming roughness matrix
>> Number of measurements read: 70
>> Total Memory required is: 0.009 Gb
Processing frame 1 - output to file f001.dat
Iteration 1
Initial RMS Misfit: 76.32 Number of data ignored: 0
Alpha: 5.743 RMS Misfit: 0.74 Roughness: 4.644
Step length set to 1.000
Final RMS Misfit: 0.74
Final RMS Misfit: 1.03
Solution converged - Outputing results to file
Calculating sensitivity map
End of data: Terminating
>> Program ended normally
----------------- Computing d-d0+f(m0) ---------------
Writing .in file and mesh.dat... done!
Writing protocol.dat... done!
Running forward model...
>> R 3 t E R T M o d e l v 2.30 <<
>> Date: 26-01-2021
>> My beautiful 3D survey
>> F o r w a r d S o l u t i o n S e l e c t e d <<
>> A d v a n c e d M e s h I n p u t <<
>> T e t r a h e d r a l E l e m e n t M e s h <<
>> Reading mesh file
>> Determining storage needed for finite element conductance matrix
>> Generating index array for finite element conductance matrix
>> Reading resistivity model from resistivity.dat
>> Memory estimates:
For 1000 measurements the memory needed is: 0.003 Gb
For 2000 measurements the memory needed is: 0.003 Gb
For 5000 measurements the memory needed is: 0.003 Gb
For 10000 measurements the memory needed is: 0.003 Gb
>> Number of measurements read: 100
>> Total Memory required is: 0.003 Gb
>> Program ended normally
60/100 reciprocal measurements found.
0 measurements error > 20 %
60/100 reciprocal measurements found.
0 measurements error > 20 %
Forward modelling done.Matching quadrupoles between surveys for difference inversion...100 in common...done in 0.017472s
--------------------- MAIN INVERSION ------------------
________________System-Check__________________
Kernel type: Linux
Processor info:
4 Threads at <= 2900.0 Mhz
Total memory = 7.7 Gb (usage = 74.5)
Wine version = 4.0
3/3 inversions completed
----------- END OF INVERSION IN // ----------
4/4 results parsed (4 ok; 0 failed)
k.showResults(index=0, pvslices=[[0],[0], np.linspace(0.175, 0.455, 6)]) # background model
Converting np.character to a dtype is deprecated. The current result is np.dtype(np.str_) which is not strictly correct. Note that np.character is generally deprecated and 'S1' should be used. Converting np.character to a dtype is deprecated. The current result is np.dtype(np.str_) which is not strictly correct. Note that np.character is generally deprecated and 'S1' should be used.
empty mesh
k.showResults(index=1, attr='difference(percent)', vmin=-50, vmax=0, color_map='Blues_r')
k.showResults(index=2, attr='difference(percent)', vmin=-50, vmax=0, color_map='Blues_r')
Download python script: nb_3d-column-time-lapse.py
Download Jupyter notebook: nb_3d-column-time-lapse.ipynb
View the notebook in the Jupyter nbviewer
Run this example interactively: 