simsopt.mhd package¶
Submodules¶
simsopt.mhd.boozer module¶
This module provides a class that handles the transformation to Boozer coordinates, and an optimization target for quasisymmetry.
- class simsopt.mhd.boozer.Boozer(equil: simsopt.mhd.vmec.Vmec, mpol: int = 32, ntor: int = 32)¶
Bases:
simsopt._core.optimizable.Optimizable
This class handles the transformation to Boozer coordinates.
A Boozer instance maintains a set “s”, which is a registry of the surfaces on which other objects want Boozer-coordinate data. When the run() method is called, the Boozer transformation is carried out on all these surfaces. The registry can be cleared at any time by setting the s attribute to {}.
- __init__(equil: simsopt.mhd.vmec.Vmec, mpol: int = 32, ntor: int = 32) → None¶
Constructor
- get_dofs()¶
This base Optimizable object has no degrees of freedom, so return an empty array
- register(s: Union[float, Iterable[float]]) → None¶
This function is called by objects that depend on this Boozer object, to indicate that they will want Boozer data on the given set of surfaces.
- Parameters
s – 1 or more surfaces on which Boozer data will be requested.
- run()¶
Run booz_xform on all the surfaces that have been registered.
- set_dofs(x)¶
This base Optimizable object has no degrees of freedom, so do nothing.
- class simsopt.mhd.boozer.Quasisymmetry(boozer: simsopt.mhd.boozer.Boozer, s: Union[float, Iterable[float]], m: int, n: int, normalization: str = 'B00', weight: str = 'even')¶
Bases:
simsopt._core.optimizable.Optimizable
This class is used to compute the departure from quasisymmetry on a given flux surface based on the Boozer spectrum.
- Parameters
boozer – A Boozer object on which the calculation will be based.
s – The normalized toroidal magnetic flux for the flux surface to analyze. Should be in the range [0, 1].
m – The poloidal mode number of the symmetry you want to achive. The departure from symmetry B(m * theta - nfp * n * zeta) will be reported.
n – The toroidal mode number of the symmetry you want to achieve. The departure from symmetry B(m * theta - nfp * n * zeta) will be reported.
normalization – A uniform normalization applied to all bmnc harmonics. If
"B00"
, the symmetry-breaking modes will be divided by the m=n=0 mode amplitude on the same surface. If"symmetric"
, the symmetry-breaking modes will be divided by the square root of the sum of the squares of all the symmetric modes on the same surface. This is the normalization used in stellopt.weight – An option for a m- or n-dependent weight to be applied to the bmnc amplitudes.
- J() → numpy.ndarray¶
Carry out the calculation of the quasisymmetry error.
- Returns
1D numpy array listing all the normalized mode amplitudes of symmetry-breaking Fourier modes of
|B|
.
- __init__(boozer: simsopt.mhd.boozer.Boozer, s: Union[float, Iterable[float]], m: int, n: int, normalization: str = 'B00', weight: str = 'even') → None¶
Constructor
- get_dofs()¶
This base Optimizable object has no degrees of freedom, so return an empty array
- set_dofs(x)¶
This base Optimizable object has no degrees of freedom, so do nothing.
simsopt.mhd.spec module¶
This module provides a class that handles the SPEC equilibrium code.
- class simsopt.mhd.spec.Residue(spec, pp, qq, vol=1, theta=0, s_guess=None, s_min=- 1.0, s_max=1.0, rtol=1e-09)¶
Bases:
simsopt._core.optimizable.Optimizable
Greene’s residue, evaluated from a Spec equilibrum
- J()¶
Run Spec if needed, find the periodic field line, and return the residue
- __init__(spec, pp, qq, vol=1, theta=0, s_guess=None, s_min=- 1.0, s_max=1.0, rtol=1e-09)¶
spec: a Spec object pp, qq: Numerator and denominator for the resonant iota = pp / qq vol: Index of the Spec volume to consider theta: Spec’s theta coordinate at the periodic field line s_guess: Guess for the value of Spec’s s coordinate at the periodic
field line
s_min, s_max: bounds on s for the search rtol: the relative tolerance of the integrator
- get_dofs()¶
This base Optimizable object has no degrees of freedom, so return an empty array
- set_dofs(x)¶
This base Optimizable object has no degrees of freedom, so do nothing.
- class simsopt.mhd.spec.Spec(filename: Optional[str] = None, mpi: Optional[simsopt.util.mpi.MpiPartition] = None, verbose: bool = True)¶
Bases:
simsopt._core.optimizable.Optimizable
This class represents the SPEC equilibrium code.
Philosophy regarding mpol and ntor: The Spec object keeps track of mpol and ntor values that are independent of those for the boundary Surface object. If the Surface object has different mpol/ntor values, the Surface’s rbc/zbs arrays are first copied, then truncated or expanded to fit the mpol/ntor values of the Spec object to before Spec is run. Therefore, you may sometimes need to manually change the mpol and ntor values for the Spec object.
- Parameters
filename – SPEC input file to use for initialization. It should end in
.sp
. Or, if None, default values will be used.mpi – A
simsopt.util.mpi.MpiPartition
instance, from which the worker groups will be used for SPEC calculations. IfNone
, each MPI process will run SPEC independently.verbose – Whether to print SPEC output to stdout.
- get_dofs()¶
This base Optimizable object has no degrees of freedom, so return an empty array
- init(filename: str)¶
Initialize SPEC fortran state from an input file.
- Parameters
filename – Name of the file to load. It should end in
.sp
.
- iota()¶
Return the rotational transform in the middle of the volume.
- run()¶
Run SPEC, if needed.
- set_dofs(x)¶
This base Optimizable object has no degrees of freedom, so do nothing.
- volume()¶
Return the volume inside the boundary flux surface.
simsopt.mhd.vmec module¶
This module provides a class that handles the VMEC equilibrium code.
- class simsopt.mhd.vmec.Vmec(filename: Optional[str] = None, mpi: Optional[simsopt.util.mpi.MpiPartition] = None)¶
Bases:
simsopt._core.optimizable.Optimizable
This class represents the VMEC equilibrium code.
The input parameters to VMEC are all accessible as attributes of the
indata
attribute. For example, ifvmec
is an instance ofVmec
, then you can read or write the input resolution parameters usingvmec.indata.mpol
,vmec.indata.ntor
,vmec.indata.ns_array
, etc. However, the boundary surface is different:rbc
,rbs
,zbc
, andzbs
from theindata
attribute are always ignored, and these arrays are instead taken from the simsopt surface object associated to theboundary
attribute. Ifboundary
is a surface based on some other representation than VMEC’s Fourier representation, the surface will automatically be converted to VMEC’s representation (SurfaceRZFourier
) before each run of VMEC. You can replaceboundary
with a new surface object, of any type that implements the conversion function ``to_RZFourier()`.VMEC is run either when the
run()
function is called, or when any of the output functions likeaspect()
oriota_axis()
are called.A caching mechanism is implemented, using the attribute
need_to_run_code
. Whenever VMEC is run, this attribute is set toFalse
. Subsequent calls torun()
or output functions likeaspect()
will not actually run VMEC again, untilneed_to_run_code
is changed toTrue
. The attributeneed_to_run_code
is automatically set toTrue
wheneverset_dofs()
is called. However,need_to_run_code
is not automatically set toTrue
when entries ofindata
are modified, or whenboundary
is modified.Once VMEC has run at least once, all of the quantities in the
wout
output file are available as attributes of thewout
attribute. For example, ifvmec
is an instance ofVmec
, then the flux surface shapes can be obtained fromvmec.wout.rmnc
andvmec.wout.zmns
.Since the underlying fortran implementation of VMEC uses global module variables, it is not possible to have more than one python Vmec object with different parameters; changing the parameters of one would change the parameters of the other.
An instance of this class owns just a few optimizable degrees of freedom, particularly
phiedge
andcurtor
. The optimizable degrees of freedom associated with the boundary surface are owned by that surface object.- Parameters
filename – Name of a VMEC input file to use for loading the initial parameters. If
None
, default parameters will be used.mpi – A
simsopt.util.mpi.MpiPartition
instance, from which the worker groups will be used for VMEC calculations. IfNone
, each MPI process will run VMEC independently.
- iter¶
Number of times VMEC has run.
- s_full_grid¶
The “full” grid in the radial coordinate s (normalized toroidal flux), including points at s=0 and s=1. Used for the output arrays and
zmns
.
- s_half_grid¶
The “half” grid in the radial coordinate s, used for
bmnc
,lmns
, and other output arrays. In contrast to wout files, this array has only ns-1 entries, so there is no leading 0.
- ds¶
The spacing between grid points for the radial coordinate s.
- __repr__()¶
Print the object in an informative way.
- aspect()¶
Return the plasma aspect ratio.
- get_dofs()¶
This base Optimizable object has no degrees of freedom, so return an empty array
- get_max_mn()¶
Look through the rbc and zbs data in fortran to determine the largest m and n for which rbc or zbs is nonzero.
- iota_axis()¶
Return the rotational transform on axis
- iota_edge()¶
Return the rotational transform at the boundary
- load_wout()¶
Read in the most recent
wout
file created, and store all the data in awout
attribute of this Vmec object.
- mean_iota()¶
Return the mean rotational transform. The average is taken over the normalized toroidal flux s.
- mean_shear()¶
Return an average magnetic shear, d(iota)/ds, where s is the normalized toroidal flux. This is computed by fitting the rotational transform to a linear (plus constant) function in s. The slope of this fit function is returned.
- run()¶
Run VMEC, if
need_to_run_code
isTrue
.
- set_dofs(x)¶
This base Optimizable object has no degrees of freedom, so do nothing.
- vacuum_well()¶
Compute a single number W that summarizes the vacuum magnetic well, given by the formula
W = (dV/ds(s=0) - dV/ds(s=1)) / (dV/ds(s=0)
where dVds is the derivative of the flux surface volume with respect to the radial coordinate s. Positive values of W are favorable for stability to interchange modes. This formula for W is motivated by the fact that
d^2 V / d s^2 < 0
is favorable for stability. Integrating over s from 0 to 1 and normalizing gives the above formula for W. Notice that W is dimensionless, and it scales as the square of the minor radius. To compute dV/ds, we use
dV/ds = 4 * pi**2 * abs(sqrt(g)_{0,0})
where sqrt(g) is the Jacobian of (s, theta, phi) coordinates, computed by VMEC in the gmnc array, and _{0,0} indicates the m=n=0 Fourier component. Since gmnc is reported by VMEC on the half mesh, we extrapolate by half of a radial grid point to s = 0 and 1.
- volume()¶
Return the volume inside the VMEC last closed flux surface.