gimli-org / gimli

Geophysical Inversion and Modeling Library :earth_africa:
https://www.pygimli.org
Other
344 stars 131 forks source link

why output file 'i' and 'u'is zero by code of 2D ERT modelling and inversion #706

Open subenyu opened 2 weeks ago

subenyu commented 2 weeks ago

Problem description

Please describe your issue here.

Your environment

*Operating system Windows, Linux or Mac? Python version: 3.10 import matplotlib.pyplot as plt import numpy as np import pygimli as pg import pygimli.meshtools as mt from pygimli.physics import ert world = mt.createWorld(start=[-50, 0], end=[50, -50], layers=[-1, -5], worldMarker=True) block = mt.createCircle(pos=[-5, -3.], radius=[4, 1], marker=4, boundaryMarker=10, area=0.1) poly = mt.createPolygon([(1,-4), (2,-1.5), (4,-2), (5,-2), (8,-3), (5,-3.5), (3,-4.5)], isClosed=True, addNodes=3, interpolate='spline', marker=5) geom = world + block + poly scheme = ert.createData(elecs=np.linspace(start=-15, stop=15, num=21), schemeName='dd') for p in scheme.sensors(): geom.createNode(p) geom.createNode(p - [0, 0.1])

Create a mesh for the finite element modelling with appropriate mesh quality.

mesh = mt.createMesh(geom, quality=34)

Create a map to set resistivity values in the appropriate regions

[[regionNumber, resistivity], [regionNumber, resistivity], [...]

rhomap = [[1, 100.], [2, 75.], [3, 50.], [4, 150.], [5, 25]]

data = ert.simulate(mesh, scheme=scheme, res=rhomap, noiseLevel=1, noiseAbs=1e-6, seed=1337)

u=np.array(data('u')) # I save u values in a array i=np.array(data('i')) # I save u values in a array

uu=u/uk

np.savetxt('u.csv',u)

print(u) print(i) ########################################################

2/05/24 - 14:46:00 - pyGIMLi - INFO - Data error estimate (min:max) 0.010000294838286121 : 0.01056761917552525 02/05/24 - 14:46:00 - pyGIMLi - INFO - 0.13172095457607066 956.0725038881354 02/05/24 - 14:46:00 - pyGIMLi - INFO - Simulated data Data: Sensors: 21 data: 171, nonzero entries: ['a', 'b', 'err', 'k', 'm', 'n', 'rhoa', 'valid'] 02/05/24 - 14:46:00 - pyGIMLi - INFO - The data contains: ['a', 'b', 'err', 'i', 'ip', 'iperr', 'k', 'm', 'n', 'r', 'rhoa', 'u', 'valid'] 02/05/24 - 14:46:00 - pyGIMLi - INFO - Simulated rhoa (min/max) 42.792496044211 104.19891011722967 02/05/24 - 14:46:00 - pyGIMLi - INFO - Selected data noise %(min/max) 1.0000294838286121 1.056761917552525 relativeError set to a value > 0.5 .. assuming this is a percentage Error level dividing them by 100 [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.

                                              1. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
    1. 0.] [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
                                              1. 0.
    1. 0.]
halbmy commented 2 weeks ago

Dear @subenyu, we had this discussion already quite often: When creating an issue, there is a template where you should fill in some information (see below). This helps other users to understand the problem and eventually to help you. I (this is my personal opinion) consider just throwing your code (not even inserting it into the reserved field where it will be displayed correctly) as impolite and not motivating others to help you. In your own interest, you should change your behaviour.

Problem description

Please describe your issue here.

Your environment

Please provide the output of print(pygimli.Report()) here. If that does not work, please give provide some additional information on your:

Operating system: e.g. Windows, Linux or Mac? Python version: e.g. 3.9, 3.10, etc.? pyGIMLi version: Output of print(pygimli.__version__) Way of installation: e.g. Conda package, manual compilation from source, etc.

Steps to reproduce

Tell us how to reproduce this issue. Ideally, you could paste the code that produces the error:

import pygimli as pg
...

Expected behavior

Tell us what should happen or what you want to achieve.

Actual behavior

Tell us what happens instead and/or provide the output of your script.

Paste your script output here.

If possible, please add one or more labels to your issue, e.g. if you expect that your issue is rather a question than a problem with the code, please add the label "question".

halbmy commented 2 weeks ago

You are modelling ERT and wonder why there is no current AND voltage. As both cannot be independently known, please guess where the modelling results are. Indeed, there could be a keyword ert.simulate(..., current=0.1) leading to fill both I and U fields.

subenyu commented 2 weeks ago

Thank you very much for your reply. But I add current=0.1, 'i' and 'u'are also 0 . data = ert.simulate(mesh, scheme=scheme, res=rhomap, noiseLevel=1,current=0.1, noiseAbs=1e-6, seed=1337)

halbmy commented 2 weeks ago

Please, try at least to read first

there could be a keyword

which means there is none but we could implement it easily.

subenyu commented 2 weeks ago

Solved it by data = ert.simulate(mesh, scheme=scheme, res=rhomap, Current = 10, sr=True,calcOnly=True, verbose=True). One anthor question that changing the Current =1 to Current =100, why the response do not change. Thank you very much.

halbmy commented 2 weeks ago

Again: There is no keyword current yet as I already explained above but there will be one in the future.

calcOnly fills the voltage and current (assuming a unit current of 1A) so you will have to scale it anyway (like you will have to do without calcOnly).