Biot-Savart volume integral

[ddudt]

Computes B_plasma in a volume as:

𝐁ᵥ(𝐫) = μ₀/4π ∫ 𝐉(𝐫') × (𝐫 − 𝐫')/|𝐫 − 𝐫'|³ d³𝐫'

This is useful for when you need B_total = B_plasma + B_external in the full computational domain, not just on the boundary surface.

[codecov[bot]]

Codecov Report

Attention: Patch coverage is 61.29032% with 12 lines in your changes missing coverage. Please review.

Project coverage is 95.21%. Comparing base (87592bc) to head (814942d).

Additional details and impacted files

@@            Coverage Diff             @@
##           master     #996      +/-   ##
==========================================
- Coverage   95.26%   95.21%   -0.05%     
==========================================
  Files          87       87              
  Lines       21719    21743      +24     
==========================================
+ Hits        20691    20703      +12     
- Misses       1028     1040      +12     

Files	Coverage Δ
desc/magnetic_fields/_core.py	96.40% <100.00%> (ø)
desc/objectives/_coils.py	99.14% <100.00%> (ø)
desc/compute/_core.py	99.83% <90.90%> (-0.17%)	:arrow_down:
desc/singularities.py	88.00% <31.25%> (-3.29%)	:arrow_down:

[github-actions[bot]]

|             benchmark_name             |         dt(%)          |         dt(s)          |        t_new(s)        |        t_old(s)        | 
| -------------------------------------- | ---------------------- | ---------------------- | ---------------------- | ---------------------- |
 test_build_transform_fft_lowres         |     +0.65 +/- 3.17     | +3.28e-03 +/- 1.60e-02 |  5.08e-01 +/- 1.0e-02  |  5.05e-01 +/- 1.3e-02  |
 test_build_transform_fft_midres         |     +0.48 +/- 1.36     | +2.80e-03 +/- 8.02e-03 |  5.91e-01 +/- 4.7e-03  |  5.89e-01 +/- 6.5e-03  |
 test_build_transform_fft_highres        |     +0.19 +/- 2.26     | +1.84e-03 +/- 2.21e-02 |  9.79e-01 +/- 1.0e-02  |  9.77e-01 +/- 2.0e-02  |
 test_equilibrium_init_lowres            |     +0.51 +/- 0.61     | +1.86e-02 +/- 2.21e-02 |  3.63e+00 +/- 1.5e-02  |  3.61e+00 +/- 1.6e-02  |
 test_equilibrium_init_medres            |     +0.10 +/- 1.63     | +4.19e-03 +/- 6.64e-02 |  4.09e+00 +/- 5.4e-02  |  4.09e+00 +/- 3.9e-02  |
 test_equilibrium_init_highres           |     +0.32 +/- 0.52     | +1.76e-02 +/- 2.87e-02 |  5.50e+00 +/- 2.1e-02  |  5.48e+00 +/- 1.9e-02  |
 test_objective_compile_dshape_current   |     +0.18 +/- 0.92     | +6.74e-03 +/- 3.45e-02 |  3.76e+00 +/- 2.5e-02  |  3.75e+00 +/- 2.4e-02  |
 test_objective_compile_atf              |     +0.08 +/- 1.25     | +6.46e-03 +/- 1.02e-01 |  8.12e+00 +/- 7.9e-02  |  8.11e+00 +/- 6.4e-02  |
 test_objective_compute_dshape_current   |     +1.45 +/- 5.20     | +1.82e-05 +/- 6.52e-05 |  1.27e-03 +/- 5.9e-05  |  1.26e-03 +/- 2.7e-05  |
 test_objective_compute_atf              |     -0.32 +/- 4.53     | -1.38e-05 +/- 1.92e-04 |  4.23e-03 +/- 1.5e-04  |  4.25e-03 +/- 1.2e-04  |
 test_objective_jac_dshape_current       |     -5.67 +/- 11.28    | -2.18e-03 +/- 4.33e-03 |  3.62e-02 +/- 3.1e-03  |  3.84e-02 +/- 3.0e-03  |
 test_objective_jac_atf                  |     +0.14 +/- 2.41     | +2.62e-03 +/- 4.53e-02 |  1.89e+00 +/- 3.8e-02  |  1.88e+00 +/- 2.4e-02  |
 test_perturb_1                          |     +0.50 +/- 0.52     | +6.47e-02 +/- 6.73e-02 |  1.30e+01 +/- 4.8e-02  |  1.30e+01 +/- 4.7e-02  |
 test_perturb_2                          |     +0.51 +/- 0.55     | +9.10e-02 +/- 9.81e-02 |  1.79e+01 +/- 8.3e-02  |  1.78e+01 +/- 5.3e-02  |
 test_proximal_jac_atf                   |     +0.16 +/- 0.96     | +1.18e-02 +/- 7.01e-02 |  7.31e+00 +/- 5.7e-02  |  7.30e+00 +/- 4.1e-02  |
 test_proximal_freeb_compute             |     +1.08 +/- 0.70     | +1.90e-03 +/- 1.23e-03 |  1.78e-01 +/- 7.3e-04  |  1.76e-01 +/- 9.8e-04  |
 test_proximal_freeb_jac                 |     +0.39 +/- 0.99     | +2.83e-02 +/- 7.19e-02 |  7.31e+00 +/- 5.2e-02  |  7.28e+00 +/- 4.9e-02  |
 test_solve_fixed_iter                   |     +0.80 +/- 8.77     | +1.17e-01 +/- 1.29e+00 |  1.48e+01 +/- 9.3e-01  |  1.47e+01 +/- 8.9e-01  |

[dpanici]

use this fxn https://github.com/PlasmaControl/DESC/blob/b10f8c8be9d3e751684e5d34abecaf5cd40b2205/desc/magnetic_fields/_core.py#L31

[dpanici]

scipy.integrate.tpl_quad to check this

[dpanici]

IAEA2020_JGeiger_Manuscript_8p_finalversion (1).pdf Contributions to Plasma Physics - 2010 - Geiger - Effects of Net Currents on the Magnetic Configuration of W7‐X (2).pdf

[f0uriest]

Several people at the scidac meeting today are using the BMW code: https://github.com/ORNL-Fusion/BMW/tree/master

refs: https://scipub.euro-fusion.org/wp-content/uploads/eurofusion/WPS1CPR17_17542_submitted-4.pdf https://bpb-us-e2.wpmucdn.com/wordpress.auburn.edu/dist/5/118/files/2022/05/ISHW2017.pdf

basically it does a biot-savart volume integral to compute A (not B), then takes the finite difference curl of A to get B. As far as I can tell from the code the source grid is hardcoded to ntheta=101, ns, nzeta taken from the vmec/mgrid file. Allegedly it works well and gives accurate results?

Could try something similar (and maybe use spectral derivatives for going A->B)

[ddudt]

Try implementing Eq. 17-18 of Peraza-Rodriguez et al. Phys. Plasmas 2017 It should be a simple linear system solve to get the Fourier-Zernike coefficients of the vector potential as: $\nabla \times (\nabla \times \mathbf{A}) = \mu_0 \mathbf{J}$ Eq. 16 might be needed for boundary conditions on $\mathbf{A}$ at $\rho=0$ and $\rho=1$. Then the magnetic field is: $\mathbf{B} = \nabla \times \mathbf{A}$ And all of this can be done in the native $(\rho,\theta,\zeta)$ coordinate system.