simsopt.util package

class simsopt.util.FocusData(filename, keep_Ic_zeros=False, downsample=1)

Bases: object

Class object for reading in FOCUS-style data files.

Parameters:

  • filename – a FOCUS file

  • keep_Ic_zeros – (optional) if True, all dipoles in the file will be retained in the data structure irrespective of the value of their Ic parameter. If False, only magnets with Ic == 1 will be imported. Default is False.

  • downsample – (optional) if set to integer value n, the initialization method will only load every nth magnet from the file. If 1, all magnets will be loaded. Default is 1.

adjust_rho(q_new)

Adjust rho according to the desired momentq

Parameters:

q_new – New exponent to use for FAMUS representation of dipole magnitudes.

flip_negative_magnets()

For magnets with negative rho values, adjust the orientation angles (mt and mp) so that the rho value can be negated to a positive value

float_inds = [3, 4, 5, 7, 8, 10, 11]
init_pol_vecs(n_pol)

Initializes arrays for the x, y, and z components of allowable polarization vectors, as an array giving the ID of the polarization vector actually used for the respective magnet (0 means magnet is off)

Parameters:

n_pol – int Number of allowed polarization vectors, typically 3 or 12.

perp_vector(inds)

Returns the x, y, and z components of a Cartesian unit vector perpendicular to the polarizaton direction of the nth magnet

print_to_file(filename)

Write the FocusData class object information to a FOCUS style file.

Parameters:

filename – string denoting the name of the output file.

propNames = ['type', 'symm', 'coilname', 'ox', 'oy', 'oz', 'Ic', 'M_0', 'pho', 'Lc', 'mp', 'mt', 'op']
read_from_file(filename, keep_Ic_zeros, downsample)
repeat_hp_to_fp(nfp, magnet_sector=1)

Duplicates the magnets to the adjacent half-period, implementing the appropriate transformations to uphold stellarator symmetry.

Notes

  • At the present time, duplicate magnets will have the same pol_id as their corresponding originals

  • At the present time, duplicate magnets will have the same name their corresponding originals.

Parameters:

  • nfp (integer) – Number of field periods in the configuration

  • magnet_sector (integer (optional)) – Sector (half-period) of the torus in which the magnets are assumed to be located. Must be between 1 and 2*nfp, inclusive. Sector 1 starts at toroidal angle 0.

unit_vector(inds)

Returns the x, y, and z components of a Cartesian unit vector in the polarizaton direction of the nth magnet

class simsopt.util.FocusPlasmaBnormal(fname)

Bases: object

bnormal_grid(nphi, ntheta, range)

Calculates the normal component of the magnetic field on a regularly-spaced 2D grid of points on the surface with a specified resolution.

Parameters:

  • nphi (int) – Number of grid points in the toroidal dimension per field period.

  • ntheta (int) – Number of grid points in the poloidal dimension.

  • range (str) – Extent of the grid in the toroidal dimension. Choices include half period, field period, or full torus.

Returns:

bnormal – Normal component of the magnetic field evaluated at each of the grid points.

Return type:

2D array

simsopt.util.MpiPartition(ngroups=None, comm_world=None)

This module contains functions related to dividing up the set of MPI processes into groups, each of which can work together. For more information, see mpi.

Parameters:

  • ngroups – The number of worker groups desired. If None, a worker group will be created for each MPI process. If a value is supplied that is larger than the number of processes, the number will be lowered to the number of processes.

  • comm_world – The MPI communicator containing all processes to split into worker groups.

simsopt.util.ngroups

The number of worker groups.

Type:

int

simsopt.util.group

Index giving the worker group that this MPI process belongs to.

Type:

int

simsopt.util.comm_world

The MPI communicator that includes all processes together.

Type:

mpi4py.MPI.Intracomm

simsopt.util.comm_groups

The MPI communicator representing the worker groups.

Type:

mpi4py.MPI.Intracomm

simsopt.util.comm_leaders

The MPI communicator that includes only leaders of worker groups.

Type:

mpi4py.MPI.Intracomm

simsopt.util.rank_world

The MPI rank in the comm_world communicator.

Type:

int

simsopt.util.rank_groups

The MPI rank in the comm_groups communicator.

Type:

int

simsopt.util.rank_leaders

The MPI rank in the comm_leaders communicator, or -1 for processes that are not group leaders.

Type:

int

simsopt.util.nprocs_world

The number of MPI processes in the comm_world communicator.

Type:

int

simsopt.util.nprocs_groups

The number of MPI processes in this process’s worker group.

Type:

int

simsopt.util.nprocs_leaders

The number of group leaders, if this process is a group leader, otherwise -1.

Type:

int

simsopt.util.proc0_world

Whether this MPI process has rank 0 in comm_world.

Type:

bool

simsopt.util.proc0_groups

Whether this MPI process has rank 0 in comm_groups, i.e. whether this process is a group leader.

Type:

bool

.

simsopt.util.build_stage_II_data_array(INPUT_DIR='./', margin_up: float = 1.5, margin_low: float = 0.5, **kwargs)

Load, filter, and compute Pareto front for coil optimization results. Filtering is done based on the following engineering constraints. Max denotes the maximum over all coils, max denotes the maximum over a single coil, mean denotes an average over the plasma surface. These numbers will need to be adjusted for other stellarator configurations and rescaled if the major radius is scaled from the 1 m baseline:

Parameters:

  • (str (INPUT_DIR) – Directory containing optimization result folders with results.json files.

  • optional) – Directory containing optimization result folders with results.json files.

  • (float (margin_low) – The upper bound margin to consider for tolerable engineering constraints.

  • 1.5) (default) – The upper bound margin to consider for tolerable engineering constraints.

  • (float – The lower bound margin to consider for tolerable engineering constraints.

  • 0.5) (default) – The lower bound margin to consider for tolerable engineering constraints.

  • (dict) (**kwargs) – max_coil_length (float, default 5): Maximum coil length. max_max_kappa (float, default 12): Maximum curvature over all coils. max_mean_squared_curvature (float, default 6): Maximum mean squared curvature. max_coil_coil_distance (float, default 0.083): Maximum coil-coil distance. max_coil_surface_distance (float, default 0.166): Maximum coil-surface distance. min_coil_coil_distance (float, default 0.083): Minimum coil-coil distance. min_coil_surface_distance (float, default 0.166): Minimum coil-surface distance. mean_abs_B (float, default 0.22): Mean absolute magnetic field. max_arclength_variance (float, default 1e-2): Maximum arclength variance. max_normalized_BdotN (float, default 4e-2): Maximum normalized BdotN. max_max_force (float, default 50000): Maximum maximum force. min_coil_length (float, default 3): Minimum coil length.

Returns:

  • df (list) (List of dictionaries containing all raw optimization results.)

  • df_filtered (list) (List of dictionaries containing filtered optimization results.)

  • df_pareto (list) (List of dictionaries containing Pareto-optimal results (minimizing normalized_BdotN and max_max_force).) – If paretoset is not available, returns an empty list.

simsopt.util.calculate_modB_on_major_radius(bs, s)

Check the average magnetic field strength along the major radius (m=n=0 mode of a SurfaceRZFourier object) to make sure the configuration is scaled correctly. For highly shaped stellarators, this can deviate a bit from the on-axis B field strength.

Parameters:
Returns:

Average magnetic field strength along

the major radius (m=n=0 mode of a SurfaceRZFourier object) of the device.

Return type:

B0avg (float)

simsopt.util.coil_optimization(s, bs, base_curves, curves, **kwargs)

Basic stellarator coil optimization wrapper function that can be reused in many different scripts and different configurations. The objective function is given by:

JF = Jf + LENGTH_WEIGHT * QuadraticPenalty(sum([CurveLength(c) for c in base_curves]), LENGTH_TARGET, “max”) + CC_WEIGHT * CurveCurveDistance(curves, CC_THRESHOLD, num_basecurves=len(base_curves)) + CS_WEIGHT * CurveSurfaceDistance(curves, s, CS_THRESHOLD) + CURVATURE_WEIGHT * sum([LpCurveCurvature(c, 2, CURVATURE_THRESHOLD) for c in base_curves]) + MSC_WEIGHT * sum([QuadraticPenalty(MeanSquaredCurvature(c), MSC_THRESHOLD) for c in base_curves]) + LINKING_NUMBER_WEIGHT * LinkingNumber(curves) + FORCE_WEIGHT * Jforce

Parameters:

  • s (Surface) – Surface representing the plasma boundary

  • bs (BiotSavart) – Biot Savart class object, presumably representing the magnetic fields generated by the coils.

  • base_curves (list[Curve]) – List of Curve class objects.

  • curves (list[Curve]) – List of Curve class objects obtained after applying symmetries to the base_curves.

  • **kwargs (dict) – Keyword arguments for the optimization. LENGTH_WEIGHT (float, default=1): Weight on the curve lengths in the objective function. LENGTH_THRESHOLD (float, default=18.0): Threshold for the curve lengths in the objective function. CC_THRESHOLD (float, default=0.1): Threshold for the coil-to-coil distance penalty in the objective function. CC_WEIGHT (float, default=1): Weight for the coil-to-coil distance penalty in the objective function. CS_THRESHOLD (float, default=1.5): Threshold for the coil-to-surface distance penalty in the objective function. CS_WEIGHT (float, default=1e-2): Weight for the coil-to-surface distance penalty in the objective function. CURVATURE_THRESHOLD (float, default=0.1): Threshold for the curvature penalty in the objective function. CURVATURE_WEIGHT (float, default=1e-9): Weight for the curvature penalty in the objective function. MSC_THRESHOLD (float, default=0.1): Threshold for the mean squared curvature penalty in the objective function. MSC_WEIGHT (float, default=1e-9): Weight for the mean squared curvature penalty in the objective function. LINKING_NUMBER_WEIGHT (float, default=0): Weight for the linking number penalty in the objective function. FORCE_THRESHOLD (float, default=0): Threshold for the force penalty in the objective function. FORCE_WEIGHT (float, default=0): Weight for the force penalty in the objective function. MAXITER (int, default=500): Maximum number of iterations for the optimization.

Returns:

Biot Savart class object, representing the

OPTIMIZED magnetic fields generated by the coils.

Return type:

bs (BiotSavart)

simsopt.util.continuation_vacuum_stage_II_optimizations(N=10000, dx=0.05, config='QA', INPUT_DIR=None, OUTPUT_DIR=None, INPUT_FILE=None, FORCE_OBJ=None, debug=False, MAXITER=14000)

Perform a continuation method on a set of previous vacuum stage II optimizations, perturbing parameters and rerunning optimizations.

Parameters:

  • N (int) – Number of new optimizations to run.

  • dx (float) – Relative perturbation size for parameters.

  • INPUT_DIR (str) – Directory containing previous optimization results (with results.json).

  • OUTPUT_DIR (str) – Directory to save new optimization results.

  • INPUT_FILE (str) – Path to VMEC input file for the target surface.

  • MAXITER (int) – Maximum number of optimizer iterations per job.

simsopt.util.discretize_polarizations(mag_data, orientation_phi, pol_axes, pol_type)

Adjusts the polarization axes of a set of dipole moments to the closest from a set of allowable axes.

Also populates the pol_x, pol_y, pol_z, and pol_id arrays within the input instance of the mag_data class, providing all of the allowable polarization vectors for each magnet in the lab frame.

Parameters:

  • mag_data – mag_data object Instance of the mag_data class characterizing the dipole moments of a set of magnets, whose orientations are to be updated.

  • orientation_phi – double array Orientation angle defining a reference axes for the allowable axes for each respective magnet (as described in the documentation for orientation_phi().) Number of elements must agree with the number of magnets in mag_data.

  • pol_axes – double array with 3 columns Set of Cartesian unit vectors defining the allowable polarization vectors for each magnet. The first, second, and third columns contain the x, y, and z components, respectively.

  • pol_type – integer array Type IDs corresponding to each row in pol_axes.

simsopt.util.edge_triplet(theta_fe, theta_fc)

Generates a set of vectors corresponding to the “edge-triplet” set of polarization types: edge, face-edge, and face-corner.

Parameters:

  • theta_fe – double Angle with respect to the nearest face-normal for the face-edge type

  • theta_fc – double Angle with respect to the nearest face-normal for face-corner type

Returns:

double array

3-column array in which each row is one unit vector in the full set

Return type:

vectors

simsopt.util.face_triplet(theta_fe, theta_fc)

Generates a set of vectors corresponding to the “face-triplet” set of polarization types: face, face-edge, and face-corner.

Parameters:

  • theta_fe – double Angle with respect to the nearest face-normal for the face-edge type

  • theta_fc – double Angle with respect to the nearest face-normal for face-corner type

Returns:

double array

3-column array in which each row is one unit vector in the full set

Return type:

vectors

simsopt.util.initial_vacuum_stage_II_optimizations(N=10000, FORCE_OBJ=None, with_force=False, debug=False, MAXITER=14000, config='QA', OUTPUT_DIR=None, INPUT_FILE=None, ncoils=5)

Perform a batch of initial random parameter optimizations for coil design for the Landreman-Paul 2021 QA or QH configurations.

Parameters:

  • N (int) – Number of optimizations to run.

  • OUTPUT_DIR (str) – Directory to save optimization results.

  • INPUT_FILE (str) – Path to VMEC input file for the target surface.

  • debug (bool) – If True, print extra diagnostics.

  • ncoils (int) – Number of unique coil shapes.

  • (dict) (**kwargs) – See vacuum_stage_II_optimization for more details.

simsopt.util.initialize_coils_for_pm_optimization(config_flag, TEST_DIR, s, out_dir='')

Initializes coils for each of the target configurations that are used for permanent magnet optimization. The total current is chosen to achieve a desired average field strength along the major radius.

Parameters:

  • config_flag – String denoting the stellarator configuration being initialized.

  • TEST_DIR – String denoting where to find the input files.

  • out_dir – Path or string for the output directory for saved files.

  • s – plasma boundary surface.

Returns:

List of CurveXYZ class objects. curves: List of Curve class objects. coils: List of Coil class objects.

Return type:

base_curves

simsopt.util.initialize_default_kwargs(algorithm='RS')

Keywords to the permanent magnet optimizers are now passed by the kwargs dictionary. Default dictionaries are initialized here.

Parameters:

algorithm – String denoting which algorithm is being used for permanent magnet optimization. Options are ‘RS’ (relax and split), ‘GPMO’ (greedy placement), ‘backtracking’ (GPMO with backtracking), ‘multi’ (GPMO, placing multiple magnets each iteration), ‘ArbVec’ (GPMO with arbitrary dipole orientiations), ‘ArbVec_backtracking’ (GPMO with backtracking and arbitrary dipole orientations).

Returns:

Dictionary of keywords to pass to a permanent magnet

optimizer.

Return type:

kwargs

simsopt.util.initialize_logging(filename: str = None, level: str = None, mpi: bool = False) → None

Initializes logging in a simple way for both serial and MPI jobs. The MPI logging uses MPILogger package.

Parameters:

  • filename – Name of file to store the logging info

  • level – Logging level. Could be ‘INFO’, ‘DEBUG’, ‘WARNING’, etc.

  • mpi – If True MPI logging is used provided mpi4py is installed.

simsopt.util.log(level: int = 20)

Turn on logging. If MPI is available, the processor number will be added to all logging entries.

Parameters:

level – Typically logging.INFO for regular output, or logging.DEBUG for more extensive output.

simsopt.util.make_optimization_plots(RS_history, m_history, m_proxy_history, pm_opt, out_dir='')

Make line plots of the algorithm convergence and make histograms of m, m_proxy, and if available, the FAMUS solution.

Parameters:

  • RS_history – List of the relax-and-split optimization progress.

  • m_history – List of the permanent magnet solutions generated as the relax-and-split optimization progresses.

  • m_proxy_history – Same but for the ‘proxy’ variable in the relax-and-split optimization.

  • pm_opt – PermanentMagnetGrid class object that was optimized.

  • out_dir – Path or string for the output directory for saved files.

simsopt.util.make_qfm(s, Bfield, n_iters=200)

Given some Bfield, compute a quadratic flux-minimizing surface (QFMS) for the field configuration on the surface s.

Parameters:

  • s (Surface) – plasma boundary surface.

  • Bfield (MagneticField or InterpolatedField) – MagneticField or InterpolatedField class object.

  • n_iters (int) – Maximum number of iterations for the optimization.

Returns:

The identified QfmSurface class object.

Return type:

qfm_surface (QfmSurface)

simsopt.util.make_stage_II_pareto_plots(df: list, df_filtered: list, OUTPUT_DIR: str = './')

Generate and save histograms comparing distributions of key metrics before and after filtering for vacuum stage II coil optimizations.

Parameters:

  • (list) (df_filtered)

  • (list)

  • (str) (OUTPUT_DIR)

simsopt.util.orientation_phi(corners_fname)

Determines the azimuthal (phi) angle that sets the orientations of rectangular prisms from an arrangement based on data from the corners file

Phi is the azimuthal angle of the normal vector of the face whose normal vector is in the “radial-like” direction with a polar angle of pi/2 radians (i.e., parallel to the x-y plane). The remaining faces will be assumed to have normal vectors in either:

  1. the z-direction, or

  2. parallel to the x-y plane with an azimuthal angle equal to the orientation phi plus pi/2 radians (the “phi-like” direction).

Parameters:

corners_fname (string) – Corners file

Returns:

Orientation phi angle (radians) of each prism in the arrangement

Return type:

orientation_phi (double array)

simsopt.util.polarization_axes(polarizations)

Populates an array of allowable polarization vectors based on an input list of names. Also outputs a 1D array with an integer corresponding to the polarization type in the order it appears in the list of names.

Parameters:

polarizations – list of strings List of names of polarization types. Each type should not be listed more than once

Returns:

2d double array

Allowable axes for the magnets expressed as unit vectors in local cylindrical coordinates. The first, second, and third columns contain the r, phi, and z components, respectively.

pol_type: 1d integer array

ID number for the polarization type associated with each row in pol_axes. The ID corresponds to the order in which the type appears in the input list (polarizations). The first type is given an ID of 1.

Return type:

pol_axes

simsopt.util.proc0_print(*args, **kwargs)
simsopt.util.read_focus_coils(filename)

Reads in the coils from a FOCUS file. For instance, this is used for loading in the MUSE phased TF coils.

Parameters:

filename – String denoting the name of the coils file.

Returns:

List of CurveXYZFourier class objects. base_currents: List of Current class objects. ncoils: Integer representing the number of coils.

Return type:

coils

simsopt.util.run_Poincare_plots_with_permanent_magnets(s_plot, bs, b_dipole, comm, filename_poincare, out_dir='')

Wrapper function for making Poincare plots.

Parameters:

  • s_plot – plasma boundary surface with range = ‘full torus’.

  • bs – MagneticField class object.

  • b_dipole – DipoleField class object.

  • comm – MPI COMM_WORLD object for using MPI for tracing.

  • filename_poincare – Filename for the output poincare picture.

  • out_dir – Path or string for the output directory for saved files.

simsopt.util.stell_point_transform(mode, phi, x_in, y_in, z_in)
Transforms a point in one of two ways, depending on the mode selected:
reflect: Reflects a point in a given poloidal plane presumed to form

the boundary between two stellarator-symmetryc half-periods. The output point will have the equivalent location associated with the adjacent half-period.

translate: Translates a point in the toroidal direction (i.e., along

a circle with a fixed radius about the z axis) by a given angle.

Parameters:

  • mode (string) – ‘translate’ or ‘reflect’ (see description above)

  • phi (double) – ‘reflect’ mode: toroidal angle (rad.) of the symmetry plane ‘translate’ mode: toroidal interval (rad.) along which to translate

  • x_in (double (possibly arrays of equal dimensions)) – x, y, z coordinates of the point(s) to be transformed

  • y_in (double (possibly arrays of equal dimensions)) – x, y, z coordinates of the point(s) to be transformed

  • z_in (double (possibly arrays of equal dimensions)) – x, y, z coordinates of the point(s) to be transformed

Returns:

x_out, y_out, z_out – x, y, z coordinates of the transformed point(s)

Return type:

double (possibly arrays)

simsopt.util.stell_vector_transform(mode, phi, vx_in, vy_in, vz_in)
Transforms a vector in one of two ways, depending on the mode selected:
reflect: Reflects a vector, with a defined origin point, in a given

poloidal plane that is presumed to form the boundary between two stellarator-symmetryc half-periods. The output vector will have the equivalent origin point and direction associated with the target half-period. Vector magnitude is preserved.

translate: Translates a vector, with a defined origin point, in the

toroidal direction. The origin thus moves along a circle of fixed radius about the z axis by a given angle. The azimuthal components of the vector change in order to preserve the radial and toroidal components. The vertical component of the vector, as well as its length, are held fixed.

Parameters:

  • mode (string) – ‘translate’ or ‘reflect’ (see description above)

  • phi (double) – ‘reflect’ mode: toroidal angle (rad.) of the symmetry plane ‘translate’ mode: toroidal interval (rad.) along which to translate

  • vx_in (double (possibly arrays of equal dimensions)) – x, y, z components of the input vectors

  • vy_in (double (possibly arrays of equal dimensions)) – x, y, z components of the input vectors

  • vz_in (double (possibly arrays of equal dimensions)) – x, y, z components of the input vectors

Returns:

x, y, z components of the output transformed vector(s)

Return type:

vx_out, vy_out, vz_out

simsopt.util.trace_fieldlines(bfield, label, s, comm, out_dir='', nfieldlines=4, tmax_fl=10000, num_iters=20000)

Make Poincare plots on a surface as in the trace_fieldlines example in the examples/1_Simple/ directory.

Parameters:

  • bfield (MagneticField or InterpolatedField) – MagneticField or InterpolatedField class object.

  • label (str) – Name of the file to write to.

  • s (Surface) – plasma boundary surface.

  • comm (MPI COMM_WORLD) – MPI COMM_WORLD object for using MPI for tracing.

  • out_dir (str) – Path or string for the output directory for saved files.

  • nfieldlines (int) – Number of fieldlines to trace.

  • tmax_fl (float) – Maximum time to trace the fieldlines.

  • num_iters (int) – Maximum number of iterations for the fieldline tracing.

simsopt.util.vacuum_stage_II_optimization(config='QA', OUTPUT_DIR=None, INPUT_FILE=None, R1=0.5, order=5, ncoils=5, UUID_init_from=None, **kwargs)

Perform a vacuum stage II force-based coil optimization with specified parameters and objectives.

Parameters:

  • (str) (config)

  • (float) (R1) – Minor radius for initial coil geometry.

  • (int) (ncoils)

  • (int)

  • None) (UUID_init_from (str or) – Universally Unique Identifiers (UUIDs) are standardized 128-bit identifiers that provide a practical way to ensure uniqueness across systems and time. UUID of previous optimization to initialize from, or None for random.

  • (dict) (**kwargs) – LENGTH_TARGET (float): Target length for coil optimization. LENGTH_WEIGHT (float): Weight for length penalty. CURVATURE_THRESHOLD (float): Threshold for curvature penalty. CURVATURE_WEIGHT (float): Weight for curvature penalty. MSC_THRESHOLD (float): Threshold for mean squared curvature penalty. MSC_WEIGHT (float): Weight for mean squared curvature penalty. CC_THRESHOLD (float): Threshold for coil-coil distance penalty. CC_WEIGHT (float): Weight for coil-coil distance penalty. CS_THRESHOLD (float): Threshold for coil-surface distance penalty. Parameters and weights for objective terms. FORCE_OBJ (class): Force objective class to use (e.g., LpCurveForce). ARCLENGTH_WEIGHT (float): Weight for arclength variation penalty. dx (float or None): Optional perturbation for initialization. with_force (bool): Whether to include force terms in the objective. debug (bool): If True, print extra diagnostics. MAXITER (int): Maximum number of optimizer iterations.

Returns:

  • res (OptimizeResult) – Result object from scipy.optimize.minimize.

  • results (dict) – Dictionary of exported optimization results and metrics.

  • coils (list) – List of optimized coil objects.