empymod icon indicating copy to clipboard operation
empymod copied to clipboard

ENH: Implementing current density as output of empymod

Open prisae opened this issue 1 year ago • 8 comments

Message from @efinden

To understand more of the free and polarization currents in the ground it would be useful to adjust the output of empymod simply by multiplying the electric field by the eta H-function, in order to obtain the total current density in the frequency domain and then use model.tem to transform it to the time-domain.

In other words, just to use Ohms law in the frequency domain: J(w) = E(w)*etaH(w).

prisae avatar Sep 11 '24 17:09 prisae

Thanks @efinden for this suggestion. I made a very rough proof-of-concept implementation in #220 - implemented, but NOT TESTED yet. Please have a go and try it out, I'd appreciate feedback!

You can install it via

python -m pip install git+https://github.com/emsig/empymod@ecurrent 

Note: No idea if this is anywhere close to correct. I will have to read into it more, and also compare it to other results. It was just a quick try!

prisae avatar Sep 11 '24 17:09 prisae

Thank you very much for the new implementation!

I did some testing today and have a question on which locations in space that etaH are computed.

Is the etaH working on the output of the electric field for every receiver-location?

I tested one case with resisivity=200 Ohm-m, for a halfspace over air of resistivity 1e40. The air layer has a nearly zero conductivity, and thus etaH for the air layer should be very low (I used res=1e-40).

The configuration is a rectangular bipole source in the air layer at z=-0.1 m. Electric source. Electric receiver. Receivers are spaced out between x,y,z = [-5,5]. The current density of the air layer is very low, as expected, but the current density is also very low for the 200 Ohm-m halfspace.

I would perhaps expect the y-component of the electric field multiplied with the etaH-function in the ground of 200 Ohm-m to be much larger. Maybe it is something in my missing understanding of the implementation that is causing this results.

efinden avatar Sep 12 '24 12:09 efinden

Good point. I multiplied it with the etaH from the source layer (as that is what is used in loop with zetaH), but it does not make sense here.

Then I have to put it into the loop where it goes over receiver depths. I will push an update.

prisae avatar Sep 12 '24 12:09 prisae

I pushed a change, could you test again? Now it should multiply it with the etaH from the receiver layer

prisae avatar Sep 12 '24 13:09 prisae

Thanks! It is fine, then I understand it as well!

A loop over the receiver depths will make it!

I will run a couple of tests, it takes some time with the filters I currently use and a lot of rec. positions, so perhaps I will not report here before tomorrow!

Looking forward to try!

efinden avatar Sep 12 '24 13:09 efinden

A test of a pure resistive 200 Ohm-m halfspace, in the diffusional model (which means setting the permittivity to zero) is done. When asking for receivers to show the y-component of the electric field, and with the ecurrent = True, the modified code seems to reproduce exactly the current density as obtained when the conductivity is multiplied with the electric field outside empymod (on the output of empymod).

Thank you very much!

I look forward to investigate the polarization currents now for an IP case. Also, I will try to model a larger domain to see if the "smoke ring" of currents caused by the TEM system behaves as expected with depth.

efinden avatar Sep 13 '24 07:09 efinden

Great! As I said, I only implemented according to your description - I haven't done any background research/read-up yet or thought deeper about it, so it is good to keep being curious about the feasibility of the outcome.

prisae avatar Sep 13 '24 09:09 prisae

Thank you! I did a bit of background check on the physics before I asked and it seems all right, but other considerations are always welcome.

I post below an image of the current density due to a TEM switch off from a 2x2 m loop lying 10 cm above the surface for a 200 Ohm-m halfspace. Only some points in space are plotted, and the surrounding points are given the same intensity in the pixels, therefore it looks like the surface layer has zero current, this is only due to the sparse plotting.

In this case I use the quasi-static approximation and force a constant permittivity to zero, therefore there is only free current shown, despite the title says total current density.

1-200-Ohm-m_Jy_xz_plane1 000000e-07s_ecurrent 2-200-Ohm-m_Jy_xz_plane4 328761e-07s_ecurrent Uploading 3-200-Ohm-m_Jy_xz_plane1.873817e-06s_ecurrent.png…

efinden avatar Sep 16 '24 12:09 efinden

Hi @efinden - just to follow up, can we close this issue, did it work alright?

prisae avatar Jan 09 '26 23:01 prisae

Hi @prisae ! Yes, it do work very well! So it can be closed. I have not had the time to do more numerical expermiments on this phenomena, but I plan to continue after submission of my PhD in fall 2026!

efinden avatar Jan 12 '26 10:01 efinden

Great, thanks for the feedback. If you ever get around to submit an example of it to the gallery that would be great of course.

Success finishing your PhD!

prisae avatar Jan 12 '26 13:01 prisae