Open matt-pharr opened 3 months ago
The OMFIT gEQDSK has a built-in method for flux-surface averages: https://github.com/gafusion/OMFIT-source/blob/df9536b5967d89b967b7d30ec9a71aa6681d0332/omfit/omfit_classes/omfit_eqdsk.py#L389
Can we use these directly?
I have been trying to avoid using OMFIT as omfit-classes both requires a ton of dependencies that would be for just reading eqdsk files, and does not work past python 3.7. Starting in python 3.8, attempting from omfit_classes import omfit_eqdsk
yields ModuleNotFoundError: No module named 'xarray.plot.plot'
Additionally, one of the priorities with this refactor was to speed up the code, and numba will likely have issues on python 3.7 for those running on new (>2021) Macs, as python 3.7 was never updated to run natively on them, and is emulated. This doesn't play nicely with compile-at-runtime libraries like numba.
We want to get the current profile for ohmic heating from the gfile. I do not believe that $|\bf{J}|$ has symmetry about a flux surface, since in equilibrium $$ \bf{J}\times\bf{B} = \nabla p $$ and $|bf{B}|$ is known to vary about a flux surface; however, I have tried to get average current via a line integral about each flux surface but am having trouble working through the interpolation step. This needs to be revisited, as I have approximated the current profile by using $R_0$ in place of $R$ in the equation derived from grad-shafranov, taking $F$, $FF'$, and $P'$ as inputs:
$$ \begin{align} |J| & = \sqrt{ \left( \frac{1}{\mu{0}}\frac{1}{R}F' \right)^2 + \left( R P' + \frac{1}{2\mu{0}}FF' \right)^2 } \ &= \sqrt{ \frac{1}{(\mu{0}R)^2}\left( \frac{\partial}{\partial\psi}F \right)^2 + \left( R P' + \frac{1}{2\mu{0}}FF' \right)^2 } \end{align} $$
Using this, one can obtain 'a' current profile, but upon integration it is clear that this has a wide margin of error. The section that extracts current density from the gfile needs to be rewritten to do an average of the quantity seen in the third plot.