Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Dp/current potential field #592

Merged
merged 65 commits into from
Oct 14, 2023
Merged

Dp/current potential field #592

merged 65 commits into from
Oct 14, 2023

Conversation

dpanici
Copy link
Collaborator

@dpanici dpanici commented Jul 21, 2023

Implements CurrentPotentialField and FourierCurrentPotentialField classes, which allow for computation of the magnetic field from a surface current density given by
$$\mathbf{K} = \mathbf{n} \times \nabla \Phi(\theta,\zeta) $$

where for the CurrentPotentialField, the current potential $\Phi$ is a user-specified, arbitrary function of $\theta$ and $\zeta$, with its derivatives also specified by the user.
For FourierCurrentPotentialField, the current potential is assumed to have the form
$$\Phi(\theta,\zeta) = \Phi_{sv}(\theta,\zeta) + \frac{G\zeta}{2\pi} + \frac{I\theta}{2\pi} $$

where $(\theta,\zeta)$ are poloidal and toroidal angles on the winding surface, $G$ is current linking the coil surface poloidally (which for modular or helical coils is determined entirely by the target plasma equilibrium, where one finds that $G = \frac{1}{\mu_0}\int_0^{2\pi} B_{\zeta}d\zeta$), and $I$ is the current linking the coil surface toroidally. $\Phi_{sv}$ is the single-valued part of $\Phi$ and is given by a DoubleFourierSeries:
$$\Phi_{sv}(\theta,\zeta) = \sum_m \sum_n \Phi_{mn} \mathcal{F}_m(\theta) \mathcal{F}_n(\zeta)$$

So for the FourierCurrentPotentialField, the free parameters of the field are the Fourier coefficients $\Phi_{mn}$, and the net currents $I$, and $G$.


TODO

  • subclass with FourierRZToroidalSurface so that can have its own Rb_lmn etc, for winding surface optimization (look at Coils classes)
  • Modify the compute magnetic fieild to take advantage of fact that it has discrete symmetry (based off Rory's FB https://github.com/f0uriest/DESC_freeb/blob/1d071712af6844bc00b4c4071f31777907b0cb50/desc/singularities.py#L278)
  • don't precompute K
  • make field assume fourier series + linear secular terms instead of custom functions (maybe keep the general one?)
  • Maybe add change_potential_resolution methods to class (to allow changing of the basis)
  • add Phi to data index
  • pass in modenumbers instead of basis for fourier currentpotenitalfield
  • maybe remove AD from general currentpotential field
    • this is because a lot of stuff is sort of annoying when dealing with the AD
  • add basis as public attribute of fourier currentpotentialfield
  • change Fouriercurrentpoetnal field to change the potential methods when basis resolution is changed, as right now it will fail if the basis resolution is changed bc the potential functions were defined with the initial basis, and are not changed when change_resolution is called

@dpanici dpanici marked this pull request as ready for review August 9, 2023 18:44
@dpanici dpanici requested review from f0uriest, ddudt and unalmis August 9, 2023 18:45
@dpanici
Copy link
Collaborator Author

dpanici commented Aug 9, 2023

This class is mainly to be used with current potentials of form $I\theta/2/\pi + G\zeta/2/\pi + \sum \Phi_{mn} sin(m\theta - n\zeta$
Unsure if it is worth subclassing for this specific implementation,or leaving it general and then having the user define their functions for the potential. Probably leaving it general is fine though I think

desc/magnetic_fields.py Outdated Show resolved Hide resolved
Copy link
Member

@f0uriest f0uriest left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Are there examples of realistic use cases where the potential would NOT be given as a fourier series + secular terms? I'm wondering if it would be better to have this assume that parameterization, so the user only needs to supply the fourier coefficients rather than functions.

  • I'm not crazy about precomputing K on a fixed grid, since I'm not sure that will work if we want to optimize K in the loop. Is it that much more expensive to evaluate on the fly? IE, the user can pass in the surface grid to compute_magnetic_field. We could also then add compute functions for Phi, K etc to desc.compute._surface.py (I have something like this on the free boundary branch)

@dpanici
Copy link
Collaborator Author

dpanici commented Aug 10, 2023

  • Are there examples of realistic use cases where the potential would NOT be given as a fourier series + secular terms? I'm wondering if it would be better to have this assume that parameterization, so the user only needs to supply the fourier coefficients rather than functions.

My only use case for this right now is the REGCOIL and LM stellarator stuff. Both of those we use fourier series and secular terms (linear in theta and zeta). It might be that for the LM stell stuff we want secular terms that are not just linear in theta and zeta, but at least right now I thin that this can be made to only be this parametrization

  • I'm not crazy about precomputing K on a fixed grid, since I'm not sure that will work if we want to optimize K in the loop. Is it that much more expensive to evaluate on the fly? IE, the user can pass in the surface grid to compute_magnetic_field. We could also then add compute functions for Phi, K etc to desc.compute._surface.py (I have something like this on the free boundary branch)

Yea this I had in mind mainly for a fixed surface, since if the surface changes at all K would need to be precomputed. K can be kind of expensive but I dont think as expensive as the biot savart loop is. I can change so it is not precomputed, and add the compute._surface functions

@f0uriest
Copy link
Member

My only use case for this right now is the REGCOIL and LM stellarator stuff. Both of those we use fourier series and secular terms (linear in theta and zeta). It might be that for the LM stell stuff we want secular terms that are not just linear in theta and zeta, but at least right now I thin that this can be made to only be this parametrization

ok yeah that sounds cleaner. And maybe rename the class to FourierCurrentPotentialField or something, to leave room for other parameterizations in the future

@f0uriest
Copy link
Member

Yea this I had in mind mainly for a fixed surface, since if the surface changes at all K would need to be precomputed. K can be kind of expensive but I dont think as expensive as the biot savart loop is. I can change so it is not precomputed, and add the compute._surface functions

This might also be a place where doing FFTs in the poloidal direction could help speed things up. Can probably be a separate PR though.

@dpanici dpanici marked this pull request as draft August 11, 2023 19:11
@codecov
Copy link

codecov bot commented Aug 14, 2023

Codecov Report

Merging #592 (f756622) into master (3bb4a37) will increase coverage by 0.05%.
The diff coverage is 99.05%.

Additional details and impacted files
@@            Coverage Diff             @@
##           master     #592      +/-   ##
==========================================
+ Coverage   94.89%   94.95%   +0.05%     
==========================================
  Files          79       79              
  Lines       19087    19294     +207     
==========================================
+ Hits        18113    18320     +207     
  Misses        974      974              
Files Coverage Δ
desc/compute/_surface.py 95.54% <100.00%> (+0.47%) ⬆️
desc/compute/data_index.py 90.00% <ø> (ø)
desc/compute/utils.py 95.93% <100.00%> (+0.07%) ⬆️
desc/geometry/surface.py 96.07% <ø> (ø)
desc/magnetic_fields.py 95.81% <98.87%> (+1.60%) ⬆️

... and 1 file with indirect coverage changes

desc/magnetic_fields.py Outdated Show resolved Hide resolved
desc/magnetic_fields.py Outdated Show resolved Hide resolved
desc/magnetic_fields.py Outdated Show resolved Hide resolved
Copy link
Member

@f0uriest f0uriest left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A few minor fixes noted. Should also try to get the coverage up a bit higher, it looks like there are a lot of conditionals etc that aren't tested.

Should also add the new classes to test_compute_everything

@dpanici dpanici requested a review from f0uriest October 12, 2023 18:18
@github-actions
Copy link
Contributor

|             benchmark_name             |         dt(%)          |         dt(s)          |        t_new(s)        |        t_old(s)        | 
| -------------------------------------- | ---------------------- | ---------------------- | ---------------------- | ---------------------- |
 test_build_transform_fft_lowres         |     -0.65 +/- 2.05     | -1.27e-04 +/- 3.98e-04 |  1.93e-02 +/- 2.6e-04  |  1.94e-02 +/- 3.0e-04  |
 test_build_transform_fft_midres         |     -1.26 +/- 0.70     | -1.42e-03 +/- 7.85e-04 |  1.11e-01 +/- 5.9e-04  |  1.13e-01 +/- 5.2e-04  |
 test_build_transform_fft_highres        |     -1.46 +/- 0.79     | -7.55e-03 +/- 4.11e-03 |  5.11e-01 +/- 3.0e-03  |  5.19e-01 +/- 2.8e-03  |
 test_equilibrium_init_lowres            |     -2.13 +/- 1.24     | -9.47e-03 +/- 5.53e-03 |  4.35e-01 +/- 3.4e-03  |  4.45e-01 +/- 4.4e-03  |
 test_equilibrium_init_medres            |     -1.21 +/- 1.77     | -6.83e-03 +/- 9.94e-03 |  5.56e-01 +/- 6.7e-03  |  5.63e-01 +/- 7.4e-03  |
 test_equilibrium_init_highres           |     -0.94 +/- 0.88     | -9.41e-03 +/- 8.81e-03 |  9.88e-01 +/- 6.0e-03  |  9.97e-01 +/- 6.4e-03  |
 test_objective_compile_dshape_current   |     -2.71 +/- 7.65     | -1.45e-01 +/- 4.10e-01 |  5.21e+00 +/- 2.8e-01  |  5.35e+00 +/- 3.0e-01  |
 test_objective_compile_atf              |     +0.62 +/- 4.65     | +9.64e-02 +/- 7.22e-01 |  1.56e+01 +/- 5.8e-01  |  1.55e+01 +/- 4.3e-01  |
 test_objective_compute_dshape_current   |     -3.12 +/- 2.28     | -1.08e-04 +/- 7.85e-05 |  3.34e-03 +/- 3.0e-05  |  3.45e-03 +/- 7.3e-05  |
 test_objective_compute_atf              |     +0.63 +/- 4.44     | +7.32e-05 +/- 5.12e-04 |  1.16e-02 +/- 4.6e-04  |  1.15e-02 +/- 2.4e-04  |
 test_objective_jac_dshape_current       |     +3.16 +/- 9.21     | +4.02e-03 +/- 1.17e-02 |  1.31e-01 +/- 1.1e-02  |  1.27e-01 +/- 4.3e-03  |
 test_objective_jac_atf                  |     -1.17 +/- 2.09     | -7.95e-02 +/- 1.42e-01 |  6.72e+00 +/- 9.1e-02  |  6.80e+00 +/- 1.1e-01  |
 test_perturb_1                          |     -0.18 +/- 14.21    | -1.61e-02 +/- 1.24e+00 |  8.74e+00 +/- 8.9e-01  |  8.76e+00 +/- 8.7e-01  |
 test_perturb_2                          |     -0.27 +/- 5.14     | -4.46e-02 +/- 8.59e-01 |  1.67e+01 +/- 3.9e-01  |  1.67e+01 +/- 7.6e-01  |

@github-actions
Copy link
Contributor

|             benchmark_name             |         dt(%)          |         dt(s)          |        t_new(s)        |        t_old(s)        | 
| -------------------------------------- | ---------------------- | ---------------------- | ---------------------- | ---------------------- |
 test_build_transform_fft_lowres         |     -2.02 +/- 2.46     | -3.93e-04 +/- 4.80e-04 |  1.91e-02 +/- 2.8e-04  |  1.95e-02 +/- 3.9e-04  |
 test_build_transform_fft_midres         |     +0.52 +/- 0.92     | +5.79e-04 +/- 1.03e-03 |  1.13e-01 +/- 7.0e-04  |  1.12e-01 +/- 7.5e-04  |
 test_build_transform_fft_highres        |     +1.15 +/- 0.72     | +5.93e-03 +/- 3.68e-03 |  5.20e-01 +/- 2.1e-03  |  5.14e-01 +/- 3.0e-03  |
 test_equilibrium_init_lowres            |     -0.75 +/- 1.50     | -3.30e-03 +/- 6.59e-03 |  4.36e-01 +/- 5.1e-03  |  4.39e-01 +/- 4.2e-03  |
 test_equilibrium_init_medres            |     -1.06 +/- 1.13     | -5.98e-03 +/- 6.39e-03 |  5.58e-01 +/- 5.0e-03  |  5.64e-01 +/- 4.0e-03  |
 test_equilibrium_init_highres           |     -0.32 +/- 0.92     | -3.20e-03 +/- 9.11e-03 |  9.88e-01 +/- 6.9e-03  |  9.91e-01 +/- 6.0e-03  |
 test_objective_compile_dshape_current   |     -0.47 +/- 6.09     | -2.51e-02 +/- 3.23e-01 |  5.27e+00 +/- 2.0e-01  |  5.30e+00 +/- 2.5e-01  |
 test_objective_compile_atf              |     +1.01 +/- 5.03     | +1.54e-01 +/- 7.67e-01 |  1.54e+01 +/- 5.7e-01  |  1.52e+01 +/- 5.1e-01  |
 test_objective_compute_dshape_current   |     -4.15 +/- 3.27     | -1.45e-04 +/- 1.14e-04 |  3.34e-03 +/- 5.8e-05  |  3.49e-03 +/- 9.8e-05  |
 test_objective_compute_atf              |     +4.35 +/- 2.17     | +4.88e-04 +/- 2.43e-04 |  1.17e-02 +/- 2.1e-04  |  1.12e-02 +/- 1.3e-04  |
 test_objective_jac_dshape_current       |     -0.71 +/- 8.39     | -9.37e-04 +/- 1.10e-02 |  1.31e-01 +/- 6.3e-03  |  1.32e-01 +/- 9.1e-03  |
 test_objective_jac_atf                  |     +0.46 +/- 1.53     | +2.98e-02 +/- 9.98e-02 |  6.54e+00 +/- 6.9e-02  |  6.51e+00 +/- 7.2e-02  |
 test_perturb_1                          |     -0.68 +/- 13.04    | -5.93e-02 +/- 1.14e+00 |  8.65e+00 +/- 7.7e-01  |  8.71e+00 +/- 8.3e-01  |
 test_perturb_2                          |     -0.57 +/- 4.90     | -9.44e-02 +/- 8.12e-01 |  1.65e+01 +/- 3.9e-01  |  1.66e+01 +/- 7.1e-01  |

desc/magnetic_fields.py Outdated Show resolved Hide resolved
@github-actions
Copy link
Contributor

|             benchmark_name             |         dt(%)          |         dt(s)          |        t_new(s)        |        t_old(s)        | 
| -------------------------------------- | ---------------------- | ---------------------- | ---------------------- | ---------------------- |
 test_build_transform_fft_lowres         |     +0.14 +/- 9.95     | +4.04e-05 +/- 2.79e-03 |  2.81e-02 +/- 2.2e-03  |  2.80e-02 +/- 1.7e-03  |
 test_build_transform_fft_midres         |     -2.48 +/- 11.80    | -3.76e-03 +/- 1.79e-02 |  1.48e-01 +/- 1.7e-02  |  1.52e-01 +/- 6.8e-03  |
 test_build_transform_fft_highres        |     -4.17 +/- 9.84     | -2.86e-02 +/- 6.76e-02 |  6.58e-01 +/- 4.9e-02  |  6.87e-01 +/- 4.6e-02  |
 test_equilibrium_init_lowres            |     +3.95 +/- 9.77     | +2.29e-02 +/- 5.65e-02 |  6.01e-01 +/- 4.7e-02  |  5.78e-01 +/- 3.1e-02  |
 test_equilibrium_init_medres            |     -2.46 +/- 10.16    | -1.90e-02 +/- 7.82e-02 |  7.51e-01 +/- 5.9e-02  |  7.70e-01 +/- 5.2e-02  |
 test_equilibrium_init_highres           |     +0.60 +/- 10.64    | +7.05e-03 +/- 1.26e-01 |  1.19e+00 +/- 6.7e-02  |  1.18e+00 +/- 1.1e-01  |
 test_objective_compile_dshape_current   |     +2.21 +/- 11.75    | +1.51e-01 +/- 8.01e-01 |  6.97e+00 +/- 7.3e-01  |  6.81e+00 +/- 3.4e-01  |
 test_objective_compile_atf              |     +1.70 +/- 10.44    | +3.37e-01 +/- 2.07e+00 |  2.02e+01 +/- 5.2e-01  |  1.98e+01 +/- 2.0e+00  |
 test_objective_compute_dshape_current   |     +4.96 +/- 10.58    | +1.85e-04 +/- 3.95e-04 |  3.92e-03 +/- 1.9e-04  |  3.74e-03 +/- 3.5e-04  |
 test_objective_compute_atf              |    +10.62 +/- 16.23    | +1.22e-03 +/- 1.86e-03 |  1.27e-02 +/- 1.2e-03  |  1.14e-02 +/- 1.4e-03  |
 test_objective_jac_dshape_current       |     +3.15 +/- 12.54    | +3.62e-03 +/- 1.44e-02 |  1.19e-01 +/- 1.3e-02  |  1.15e-01 +/- 7.0e-03  |
 test_objective_jac_atf                  |     +8.87 +/- 5.79     | +6.54e-01 +/- 4.27e-01 |  8.03e+00 +/- 2.7e-01  |  7.38e+00 +/- 3.3e-01  |
 test_perturb_1                          |     -1.59 +/- 15.38    | -1.74e-01 +/- 1.68e+00 |  1.07e+01 +/- 1.1e+00  |  1.09e+01 +/- 1.3e+00  |
 test_perturb_2                          |    +10.69 +/- 7.71     | +2.09e+00 +/- 1.51e+00 |  2.17e+01 +/- 3.4e-01  |  1.96e+01 +/- 1.5e+00  |

@github-actions
Copy link
Contributor

|             benchmark_name             |         dt(%)          |         dt(s)          |        t_new(s)        |        t_old(s)        | 
| -------------------------------------- | ---------------------- | ---------------------- | ---------------------- | ---------------------- |
 test_build_transform_fft_lowres         |     -1.95 +/- 1.79     | -3.85e-04 +/- 3.55e-04 |  1.94e-02 +/- 3.0e-04  |  1.98e-02 +/- 1.9e-04  |
 test_build_transform_fft_midres         |     -0.48 +/- 0.70     | -5.39e-04 +/- 7.97e-04 |  1.13e-01 +/- 6.4e-04  |  1.13e-01 +/- 4.7e-04  |
 test_build_transform_fft_highres        |     -0.57 +/- 1.02     | -2.97e-03 +/- 5.29e-03 |  5.16e-01 +/- 4.1e-03  |  5.19e-01 +/- 3.4e-03  |
 test_equilibrium_init_lowres            |     -0.73 +/- 2.61     | -3.32e-03 +/- 1.19e-02 |  4.51e-01 +/- 8.8e-03  |  4.55e-01 +/- 8.0e-03  |
 test_equilibrium_init_medres            |     -2.00 +/- 1.23     | -1.16e-02 +/- 7.14e-03 |  5.68e-01 +/- 4.9e-03  |  5.79e-01 +/- 5.1e-03  |
 test_equilibrium_init_highres           |     -0.60 +/- 0.99     | -6.08e-03 +/- 9.97e-03 |  1.00e+00 +/- 5.0e-03  |  1.01e+00 +/- 8.6e-03  |
 test_objective_compile_dshape_current   |     -2.02 +/- 6.29     | -1.11e-01 +/- 3.45e-01 |  5.37e+00 +/- 2.1e-01  |  5.48e+00 +/- 2.7e-01  |
 test_objective_compile_atf              |     +3.54 +/- 6.85     | +6.03e-01 +/- 1.17e+00 |  1.77e+01 +/- 8.6e-01  |  1.70e+01 +/- 7.9e-01  |
 test_objective_compute_dshape_current   |     +4.35 +/- 6.22     | +1.48e-04 +/- 2.12e-04 |  3.56e-03 +/- 1.6e-04  |  3.41e-03 +/- 1.3e-04  |
 test_objective_compute_atf              |     +5.06 +/- 3.92     | +5.64e-04 +/- 4.38e-04 |  1.17e-02 +/- 3.3e-04  |  1.12e-02 +/- 2.9e-04  |
 test_objective_jac_dshape_current       |     -0.12 +/- 10.54    | -1.66e-04 +/- 1.46e-02 |  1.38e-01 +/- 6.3e-03  |  1.38e-01 +/- 1.3e-02  |
 test_objective_jac_atf                  |     +1.18 +/- 3.28     | +9.54e-02 +/- 2.65e-01 |  8.17e+00 +/- 2.1e-01  |  8.07e+00 +/- 1.6e-01  |
 test_perturb_1                          |     +2.02 +/- 17.65    | +1.81e-01 +/- 1.58e+00 |  9.14e+00 +/- 1.1e+00  |  8.96e+00 +/- 1.1e+00  |
 test_perturb_2                          |     -2.11 +/- 6.15     | -3.69e-01 +/- 1.07e+00 |  1.71e+01 +/- 4.5e-01  |  1.75e+01 +/- 9.7e-01  |

@dpanici dpanici merged commit b0f8691 into master Oct 14, 2023
@dpanici dpanici deleted the dp/CurrentPotentialField branch October 14, 2023 14:51
Copy link
Collaborator

@unalmis unalmis left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

approving

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants