Skip to content

Commit 9987e61

Browse files
committed
DOC: extend inversion guide a bit
1 parent 818c62b commit 9987e61

File tree

2 files changed

+83
-47
lines changed

2 files changed

+83
-47
lines changed

doc/BUILD_LINUX.md

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -29,11 +29,9 @@ The workflow have several targets which you might try:
2929

3030
If something goes wrong you can try the manual compilation.
3131

32-
3332
## Manual compiling with a virtual environment
3433

35-
If you don't want to use the conda environment we encourage the use of a
36-
virtual python environment.
34+
If you don't want to use the conda environment we encourage the use of a virtual environment.
3735
Assuming you have a proper build toolchain and the required libraries (see Installation on Ubuntu below) installed.
3836

3937
First we need to create a root directory for our installation,
@@ -107,6 +105,7 @@ Last try on Ubuntu 22.04.03 (23-11-14)
107105
libomp-dev doxygen \
108106
libcppunit-dev clang
109107
```
108+
110109
Rest see above.
111110

112111
## Useful cmake settings
@@ -119,23 +118,27 @@ CLEAN=1 cmake ../gimli
119118
```
120119

121120
Use an alternative C++ compiler:
121+
122122
```bash
123123
CC=clang CXX=clang++ cmake ../gimli
124124
```
125125

126126
Build the library with debug and profiling flags:
127+
127128
```bash
128129
cmake ../gimli -DCMAKE_BUILD_TYPE=Debug
129130
```
130131

131132
Build the library with GCC address sanitizer check:
133+
132134
```bash
133135
cmake ../gimli -DCMAKE_BUILD_TYPE=Debug -DASAN=1
134136
```
135137

136138
## Useful make commands
137139

138140
Show more verbose build output:
141+
139142
```bash
140143
make VERBOSE=1
141144
```

doc/user-guide/inversion.md

Lines changed: 77 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -135,6 +135,7 @@ from pygimli.physics import ves
135135
ab2 = np.logspace(0, 2.5, 21)
136136
synth = [10, 10, 100, 300, 30]
137137
data = ves.VESModelling(ab2=ab2).response(synth)
138+
data *= (np.random.randn(len(data))*0.03 + 1)
138139
thk = np.logspace(0, 1.8, 23)
139140
fop = ves.VESRhoModelling(ab2=ab2, thk=thk)
140141
```
@@ -198,21 +199,80 @@ print(inv.chi2(newresponse))
198199
which fits a parabola through the old point ($\tau$=0), the full step ($\tau$=1) and a test ($\tau$=0.3) and optain an optimum line search parameter of about 0.6-0.65.
199200
As method, we can also use `'exact'` (forward calculations) or `'inter'` (interpolation), yielding almost the same results.
200201
The latter is the simplest one and the former takes the most effort.
201-
In total, the chi-square misfit decreases, but slowly.
202+
In total, the chi-square misfit, computed by $\Phi_d/N$, decreases slowly.
202203

203-
### Controlling inversion
204+
As gradient-based minimization converges much slower, we switch to a Gauss-Newton framework.
205+
After initialization, we set the model transformations as strings (we can also create instances)
206+
We can choose `lin`, `log`, `logL` (L being lower bound), `logL-U` (two bounds), `cotL-U` or `symlogT` (T being the linear threshold).
207+
The model roughness vector (including model transformation and weighting) can be
208+
accessed by `inv.roughness()`.
204209

205-
The data errors $\epsilon_i$ play a crucial role in the inversion process as they
206-
define the data weighting matrix $\mathbf{W_d}$. Sometimes they are well known from
207-
statistical considerations, but often they need to be estimated. We usually consider
208-
a so-called error model consisting of a relative and absolute errors.
210+
```{code-cell}
211+
from pygimli.frameworks.inversion import GaussNewtonInversion
212+
213+
inv = GaussNewtonInversion(fop=fop)
214+
inv.modelTrans = 'log' # already default
215+
inv.dataTrans = 'log' # default linear
216+
```
217+
218+
Like the transformations, there are a lot of options that can be set directly to the inversion instance:
219+
220+
- `fop` - the forward operator
221+
- `robustData`, `blockyModel` - use L1 norm for data misfit and model roughness
222+
- `verbose` - to see some output
223+
- `model` - the current model
224+
- `response` - the model response
225+
- `dataVals`, `errorVals` - data and error vectors
209226

227+
Most of them can also be passed to the inversion run and should better
228+
229+
- `maxIter` - maximum iteration number
230+
- `lam` - the overall regularization strength
231+
- `zWeight` - the vertical-to-horizontal regularization ratio (2D/3D problems)
232+
- `startModel` - the starting model as float or array
233+
- `relativeError` and `absoluteError` to define the error model
234+
- `limits` - list of lower and upper parameter limits (overriding `inv.modelTrans`)
235+
236+
After running the inversionq
237+
238+
```{code-cell}
239+
model = inv.run(data, relativeError=0.03, verbose=True)
240+
```
241+
242+
we observe that the data are fitted within noise in very few iterations.
243+
The chi-square value can be accessed by `inv.chi2()`, its convergence is stored in
244+
`inv.chi2History`. The data, model and total objective function values can be retrieved
245+
by `inv.phiData()`, `inv.phiModel()` and `inv.phi()`. By default, the current model and
246+
its response are used, alternatively you can pass `model=` to `phiModel()` or `phi()`
247+
and `response=` to `phiData()` and `phi()`.
248+
The important measure of data fit is the chi-square value
249+
250+
$$ \chi^2 = \frac{\Phi_\text{d}}{N} = \frac1N\sum_{i=1}^{N} \left( \frac{f_i(\mathbf{m}) - d_i}{\epsilon_i} \right)^2 $$
251+
252+
as it includes the error model (and the data transformation).
253+
In many cases, one has a better feeling by computing the (untransformed)
254+
root-mean-square (RMS), either absolute
255+
256+
$$ \text{ARMS} = \sqrt{\frac1N \sum_{i=1}^{N} (f_i(\mathbf{m}) - d_i)^2} $$
257+
258+
by using `inv.absrms()`, or relative
259+
260+
$$ \text{RRMS} = \sqrt{\frac1N \sum_{i=1}^{N} \left( \frac{f_i(\mathbf{m}) - d_i}{d_i}\right)^2} $$
261+
262+
by using `inv.relrms()`.
263+
Traveltime tomography is a good example for looking at ARMS (e.g. in ms), while in ERT
264+
one usually has a good feeling for RRMS (same for voltage, resistance or apparent
265+
resistivity) in %.
266+
267+
### Data errors
268+
269+
The data errors $\epsilon_i$ play a crucial role in the inversion process as they
210270
In case the errors are known, the regularization strength $\lambda$ needs to be chosen
211271
so that the mean squared data misfit $\chi^2$ reaches a value of about 1, indicating
212272
that the data are explained within their errors. To this end, the regularization
213273
strength should be adjusted so that, according to Occams razor principle, the simplest
214274
model explaining the data within their errors is found. In many cases, simple means
215-
smooth, so that the smoothest model reaching a $\chi^2$ of about 1 is sought.
275+
smooth, so that the smoothest model reaching $\chi^2\approx 1$ is sought.
216276

217277
However, one can adapt the regularization in a way that the meaning of simple reflects
218278
the prior knowledge or assumption about the subsurface. This could be different weights
@@ -233,43 +293,16 @@ Examples that are already implemented in pyGIMLi are for example:
233293
- **Structurally coupled cooperative inversion** of disparate data based on structural similarity (e.g., {cite}`RonczkaHelGueWisDah2017NSG`
234294
- **Structure-based inversion** using layered 2D models {cite}`AttwaAkcBasGue2014JAG`
235295

236-
+++
237-
238-
## Input data
239-
240-
### Data weights / errors
241-
242-
### Data transforms
243-
244-
+++
245-
246-
## Model parametrization
247-
248-
### Mesh-free inversion (0-D)
249-
250-
### Mesh inversion
251-
252-
#### 1-D
253-
254-
#### 2-D
255-
256-
#### 3-D
257-
258-
+++
259-
260-
## Regularization - Including prior information
261-
262-
### Starting model
263-
264-
### Reference model
265-
266-
### Parameter limits
267-
268-
### Damping
269-
270-
### Smoothing
271-
272-
### Advanced regularization
296+
<!-- Model parametrization
297+
Mesh-free inversion (0-D)
298+
Mesh inversion (1D, 2D, 3D)
299+
Regularization - Including prior information
300+
Starting model
301+
Reference model
302+
Parameter limits
303+
Damping
304+
Smoothing
305+
Advanced regularization -->
273306

274307
+++
275308

0 commit comments

Comments
 (0)