Closed soehag closed 2 months ago
I don't think the C matrix is the problem, but rather the line:
manager.setData(np.random.random)
The argument of setData
should be a DataContainerERT
and not a function, not even a vector as a result of the function. It should be
manager.setData(syn_data)
(and we should check the type).
Hi, sorry for my silly mistake - that might be something left behind from me experimenting. This is obviously not supposed to be there. I edited the code - replacing np.random.random with syn_data.
Hi, sorry for my silly mistake - that might be something left behind from me experimenting. This is obviously not supposed to be there. I edited the code - replacing np.random.random with syn_data.
My intention was to say - the error persists. This seemed to be a copying mistake from my side.
Dear Hagen @soehag,
interesting issue. The example runs fine if you use the pyGIMLi data type for identity matrices (C=pg.matrix.IdentityMatrix(len(final_mesh.cells()))
) or a SparseMapMatrix
. Other matrix types seem to segfault.
@carsten-forty2 ?
C_npy = np.eye(len(final_mesh.cells()))
from scipy import sparse
C1 = pg.matrix.IdentityMatrix(len(final_mesh.cells())) # works
C3 = pg.utils.numpy2gmat(C_npy) # segfaults, similiar to original issue
C4 = pg.utils.toSparseMatrix(C_npy) # segfaults
C5 = pg.utils.toSparseMatrix(sparse.csr_matrix(C_npy)) # segfaults
C6 = pg.matrix.SparseMapMatrix(C5) # works
The problem is the lazy evaluation of the mesh with its associated regularization settings that overwrites your custom C.
The final mesh for the forward operator FOP is set inside the manager.invert
call. After that, all regularization parameter are forwarded to the regionmanager. This clears all prior settings like your C. The log info shows when a new mesh is set to the FOP.
The FIX added a custom C into these lazy evaluation process. Until the FIX (dev branch) comes with the next release, you can force the lazy evaluation by calling the mesh once before you set C with manager.fop.mesh()
. Also you should remove the mesh setting inside the invert call because this restarts the C cleaning again.
Thank you very much for your contributions and the quick fix - that seams to work well! Thank you also for your hint of not passing the mesh along again - that was also a remainder of trying different settings.
On a side note - if you do not set the regularisation with a custom constraints matrix by inv.setRegularization(), but pass it as manager.invert(..., C=C), pyGIMLi throws the error:
"Core - ERROR - no cWeights defined. You should create constraints matrix first."
but procedes with the inversion. Afterwards the constraints matrix can be found in the forward-operator, so it seems like it was used in the inversion. Maybe, if the matrix is properly used, it is reasonable to just put an info and set the cWeights automatically to 1 or throw an error and not conduct the inversion? Right now it is not entirely clear, what pyGIMLi does in this situation especially since manager.invert() doesnt accept arguments for cWeights.
Cheers and thank you again for you help!
Problem description
Hi everyone,
I am running ERT-inversions with custom constraint matrices to account for the specific geometry in my setup. Running these inversions on the my local HPC, I run into segmentation fault errors when the constraint matrix is being build. I have attached a small minimum working example, where an identity matrix as constraint matrix crashes the jupyter kernel before it can properly throw an exception. Running as a script I encounter the same error in a mwe on my local machine.
Your environment
Date: Wed Feb 28 11:09:28 2024 CET
Python 3.10.13 | packaged by conda-forge | (main, Oct 26 2023, 18:07:37) [GCC 12.3.0]
Steps to reproduce
Tell us how to reproduce this issue. Ideally, you could paste the code that produces the error:
Expected behavior
Run properly.
Actual behavior
... Building constraints matrix Segmentation fault (core dumped)