Running inversion multiple times give different result

[makeabhishek]

Problem description

1. I'm running traveltime inversion. While running the model 'velInv ' multiple times its generating different results. I'm not sure why? Does the model inverting previous mesh and not initialise the model from zero slowness?
2. How do we check, what is the initial model for the inversion and plot it?
3. Also, calculating misfit is giving error.

Your environment

       OS : Windows
        CPU(s) : 64
       Machine : AMD64
  Architecture : 64bit
           RAM : 511.9 GiB
   Environment : Jupyter

Python 3.9.18 | packaged by conda-forge | (main, Dec 23 2023, 16:29:04) [MSC
v.1929 64 bit (AMD64)]

       pygimli : 1.4.3
        pgcore : 1.4.0
         numpy : 1.25.0
    matplotlib : 3.7.2
         scipy : 1.11.2
       IPython : 8.12.0
       pyvista : 0.43.2

Steps to reproduce

mgr = tt.TravelTimeManager(data=data, verbose=True)
velInv = mgr.invert(data, mesh=mesh, lam=1, useGradient=1, zWeight=15.0, secNodes=5, absoluteError=error,
                    blockyModel = True, robustData = True,
                    lambdaFactor=1,
                    # cType=1,
                    #startModel=1/2200, 
                    isReference=False,
                    showProgress= False,
                    #verbose=True,
                   )

 # Misfit calculation
misfit = mgr.inv.response / mgr.data['t'] * 100 - 100
pg.show(mgr.data, misfit, cMap="bwr", cMin=-10, cMax=10, label="misfit (%)")

ERROR - None.show(C:\Users\376189\Documents\pyGimli\pygimli\gimli\pygimli\viewer\showmesh.py:159)
Can't interprete obj: Data: Sensors: 16 data: 120, nonzero entries: ['err', 'g', 's', 't', 'valid'] to show.
(None, None)

Expected behavior

Running velInv should generate same inverted result every time I run it. But its giving me different results. Seems like inversion is updating the earlier model without clearing the memory.

[halbmy]

Thank you for reporting. This behaviour should, of course, not happen. Please make sure it also occurs with the newest version, v1.5.0

Can you, please, specify the full code and all required files so that we can reproduce the problem?

[makeabhishek]

1. • Thanks for your reply @halbmy. Updating to latest version (1.4.6) still have the same issue.
2. • In addition to that, plotting figure is also an issue. When I'm plotting results in cell 2, its not generating mesh. However running cell 2 second time, it generate figure.
3. • How can I assign an initial model during inversion , say for lambda=1. Assume if I run the inversion once and the inverted model to another inversion with different parameters, say for lambda=5.

gimli.zip

[halbmy]

1. The latest version is v1.5.0.
2. estimating the error after initializing the manager and calling invert is useless. Better do it before only once. Apparently the error in the file is taken.
3. Your traveltimes are in the order of 6-50µs. So an absolute error of 1e-3 (1ms) does not make much sense as a homogeneous model is already fitting the data.
4. zWeight=0 is not recommended as there will be no smoothness in vertical direction.
5. A gradient model does not make much sense for a circular geometry, here a constant model is better.
6. You should use the highest lambda that is able to fit your data. robustData only makes sense in case of outliers.
7. When using meaningful parameters, I am getting the same result. No idea where your difference comes from.
8. You can specify a starting model by startModel= either as a float or a vector.

I don't understand your second point. All works well.

[makeabhishek]

I used conda update -c gimli -c conda-forge pygimli and it updated my pygimli version from v1.4.3 to v1.4.6. Not sure why it didn't update to v1.5.0.

If I'm runnig cell 2. It's generating figure like that. It only generate the second subplot after running the cell twice.

What is the meaning of nBounds in constraint matrix constraint matrix of size(nBounds x nModel) 3833 x 2614
What is the abbreviation means CGLSCDWWtrans

[halbmy]

A core update (1.4 to 1.5) often means some update to other libraries that is sometimes not done unless explicitly written (pygimli=1.5), but then it is better just to create a new environment.

[halbmy]

For me, both plots are shown in both jupyter lab and VSCode

[halbmy]

nBounds is the number of internal boundaries (edged between two triangles) used for defining the smoothness.

[halbmy]

CSLGCDWWtrans is the inverse solver used, see Günther et al. (2006) or Günther(2004)

[makeabhishek]

It is strange that I installed the pygimli new version in a new environment using conda create -n pg -c gimli -c conda-forge "pygimli>=1.5.0" but it is still showing me v1.4.6. However, pgcore : 1.5.0. Outputting print(pygimli.Report())

Date: Wed May 08 13:55:34 2024 Mountain Daylight Time

                OS : Windows (10 10.0.17763 SP0 Multiprocessor Free)
            CPU(s) : 64
           Machine : AMD64
      Architecture : 64bit
               RAM : 511.9 GiB
       Environment : Jupyter

  Python 3.11.9 | packaged by conda-forge | (main, Apr 19 2024, 18:27:10) [MSC
  v.1938 64 bit (AMD64)]

           pygimli : 1.4.3
            pgcore : 1.5.0
             numpy : 1.26.4
        matplotlib : 3.8.4
             scipy : 1.13.0
           IPython : 8.22.2
           pyvista : 0.43.7

whereas outout of print(pygimli.__version__) v1.4.6

• Also, installing new version ging error in plotting mgr.showData(ax=ax1, cMap="Spectral_r") AttributeError: 'Colorbar' object has no attribute 'draw_all'

---

As suggested, changing the parameters and error part in the inversion is still giving in consistent result. As shown below,

# Parameter 1
for run in range(2):
    velInv = mgr.invert(data, mesh=mesh, lam=10, useGradient=0, zWeight=10.0, secNodes=5,
                        blockyModel = True, robustData = True,lambdaFactor=1,
                        isReference=False,showProgress= False,verbose=True)
    mgr.showResult(cMap="Spectral_r",cMin=2800,cMax=3600); plt.title('Run ' + str(run));plt.show();

# Parameter 2
for run in range(2):
    error = 1e-3 
    pg.info("Traveltime Inversion")   
    velInv = mgr.invert(data, mesh=mesh, lam=1, useGradient=1, zWeight=0.0, secNodes=5, absoluteError=error,
                        blockyModel = True, robustData = True,
                        lambdaFactor=1,
                        # cType=1,
                        #startModel=1/2200, 
                        isReference=False,
                        showProgress= False,
                        verbose=True,
                       )
    mgr.showResult(cMap="Spectral_r",cMin=2800,cMax=3600); plt.title('Run ' + str(run));plt.show();

---

I followed the reference, where I understood CGLSCD, CGLSI, CGLAN, but didnot find the details about CGLSCDWWtrans
CG/Lanczos least squares solver (CGLAN)
Conjugate Gradients, weighted least squares (CGLSCD)
Conjugate Gradients, applied to damped least squares (CGLSI)
CGLSCDWWtrans = ?

[halbmy]

The inverse solver just solves the system of equations, but does not determine the solution. Some of the solvers in my PhD thesis do just not apply and others just incorporate additional weightings. If is all explained there. CGLSCDWWtrans is the one that is used.

[halbmy]

As discussed in the #689, both values of 0 and 10 make much sense.

Nevertheless, the results of two runs with the same parameters should be identical, at least within numerical accuracy and we don't know where the difference comes from but keep digging. Looking at the ratio between two runs, the lower the impact of regularization is (low lamda or lamdaFactor, blockyModel), i.e. the more we are inverting noise.