[WIP] TDEM loop source

[lheagy]

Summary

Implement the B field computation from a closed-loop wire source in the TDEM code. We first compute the vector potential using the Biot-Savart Law (implemented in geoana) and then take the curl to get the initial B-field

PR Checklist

• [x] If this is a work in progress PR, set as a Draft PR
• [x] Linted my code according to the style guides.
• [ ] Added tests to verify changes to the code.
• [ ] Added necessary documentation to any new functions/classes following the expect style.
• [ ] Marked as ready for review (if this is was a draft PR), and converted to a Pull Request
• [ ] Tagged @simpeg/simpeg-developers when ready for review.

Reference issue

#1072

What does this implement/fix?

Adds the computation of B with a step-off waveform using the B or H formulations.

Additional information

Note that this current implementation would not work for a permeable background.

• [ ] todo: add a check + not implemented error if mu != mu_0 (this can be a separate pr)
• [ ] include tests: check with a circular loop source that is small enough and the fields should match. Also add derivative tests for the inversion

Input appreciated!

@dccowan and @aguspesce: this PR connects with work that you are both doing, and if you have the time, your thoughts would be appreciated! I have done a quick test, and the results pass the "eyeball" test, but I know that you both have notebooks where you have computed dbdt from a triangular waveform, which should be the same thing as b from a step-off waveform. If you have the time to compare, that would be awesome!

This notebook is adapted from our consortium example and runs a 3D forward problem on a tensor mesh, now using B-field sensors with a step-off waveform: 3-forward-Tensor.ipynb.zip

[codecov[bot]]

Codecov Report

:x: Patch coverage is 88.28829% with 13 lines in your changes missing coverage. Please review. :white_check_mark: Project coverage is 81.47%. Comparing base (4378511) to head (545bce6). :warning: Report is 1 commits behind head on main.

Files with missing lines	Patch %	Lines
simpeg/electromagnetics/time_domain/sources.py	73.33%	5 Missing and 7 partials :warning:
tests/em/tdem/test_large_loop.py	98.48%	0 Missing and 1 partial :warning:

Additional details and impacted files

@@           Coverage Diff           @@
##             main    #1651   +/-   ##
=======================================
  Coverage   81.47%   81.47%           
=======================================
  Files         419      420    +1     
  Lines       54969    55056   +87     
  Branches     5227     5240   +13     
=======================================
+ Hits        44784    44856   +72     
- Misses       8786     8791    +5     
- Partials     1399     1409   +10     

:umbrella: View full report in Codecov by Sentry.
:loudspeaker: Have feedback on the report? Share it here.

:rocket: New features to boost your workflow:

• :snowflake: Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[aguspesce]

@lheagy ,

Thanks for taking the time to code this.

I adapted one of my notebooks where I compare the step-off and UTEM solutions for square and circular loops.

Simulation 1: UTEM Survey

• Receiver: db/dt
• Two source loops: one circular and one square, both with a triangular waveform.

# Define a list with all the receivers
receivers_dbdt = [
    tdem.receivers.PointMagneticFluxTimeDerivative(
        receiver_locations, time_channels, s
    ) for s in components_title
]
# If mamake_check is activate, Add otherr source as a circular loop
source_utem_line = tdem.sources.LineCurrent(
            receivers_dbdt, 
            location=tx_points, 
            waveform=general_waveform, 
            srcType="inductive", 
            current=max_current
)
source_utem_circular = tdem.sources.CircularLoop(
        receivers_dbdt, location=[0,0,0], radius=tx_halfwidth,
        waveform=general_waveform, srcType="inductive", current=max_current
    )
# Create survey
survey_utem = tdem.Survey([source_utem_line, source_utem_circular])

Simulation 2: Step-Off Survey

• Receiver: B
• Two source loops: one circular and one square, both with a step-off waveform.

tep_off = tdem.sources.StepOffWaveform()
receivers_b = [
    tdem.receivers.PointMagneticFluxDensity(
        receiver_locations, time_channels, s
    ) for s in components_title
]
source_stepoff_circular = tdem.sources.CircularLoop(
    receiver_list=receivers_b,
    location=[0, 0, 0],
    radius=tx_halfwidth,
    waveform=step_off, 
    srcType="inductive", 
    current=2*max_current
)
source_stepoff_line = tdem.sources.LineCurrent(
    receivers_b, 
    location=tx_points, 
    waveform=step_off, 
    srcType="inductive", 
    current=2*max_current
)
# Create survey
survey_stepoff = tdem.Survey([source_stepoff_line, source_stepoff_circular])

I compared the bz obtained from the UTEM simulation, the step-off simulation, and the analytical solution calculated using geoana:

bz_geoana = geoana.em.tdem.vertical_magnetic_flux_horizontal_loop(
        time_channels, radius=tx_halfwidth, sigma=background_conductivity, 
        current=2*max_current, turns=1
)

In the case of the sources with a circular shape, the solutions are approximately the same. However, for the square loop, the result differs — by roughly a factor of 2. I think the current is not being taken into account in the calculation.

I attached my notebook: fwd_3d_check_utem_stepoff-square_loop.zip

[lheagy]

Thanks @aguspesce! Yes, you are totally right, it was missing the current. I have added that now. When running your notebook, I now get much better agreement

[santisoler]

@lheagy, I see some lines not being covered by tests. Let me know if you need any help here.

[aguspesce]

Hello @lheagy, I updated the notebook to plot the results for a line of receivers in all three components. I also added a circular loop generated by 30 segments using a LineCurrent to compare the prediction with the analytical solution.

When I plot the difference between the different methods of calculating the B-field for a circular loop (specifically the Z-component at the center of the loop), I noticed a discrepancy of about 1 pT between the step-off and UTEM methods. This is probably because, to generate the secondary field, I’m subtracting the result from the last time channel — but I’m not entirely sure.

Sorry if the notebook is a bit messy: fwd_3d_check_utem_stepoff-square_loop_v2.zip

[lheagy]

Let's address #1575 before this one. I would like to make a few updates here:

• [x] use the mesh.project_face/edge_vector rather than assuming x, y, z coordinates (thanks @jcapriot for the suggestion!)
• [x] double check the boundary conditions when using the vector potential
• [x] improve testing