zoidberg package

Module contents

zoidberg.make_maps(grid, magnetic_field, nslice=1, quiet=False, **kwargs)

Make the forward and backward FCI maps

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • nslice (int) – Number of parallel slices in each direction
  • quiet (bool) – Don’t display progress bar
  • kwargs – Optional arguments for field line tracing, etc.
Returns:

Dictionary containing the forward/backward field line maps

Return type:

dict

zoidberg.write_maps(grid, magnetic_field, maps, gridfile='fci.grid.nc', new_names=False, metric2d=True, format='NETCDF3_64BIT', quiet=False)

Write FCI maps to BOUT++ grid file

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • maps (dict) – Dictionary of FCI maps
  • gridfile (str, optional) – Output filename
  • new_names (bool, optional) – Write “g_yy” rather than “g_22”
  • metric2d (bool, optional) – Output only 2D metrics
  • format (str, optional) – Specifies file format to use, passed to boutdata.DataFile
  • quiet (bool, optional) – Don’t warn about 2D metrics
Returns:

Return type:

Writes the following variables to the grid file

Submodules

zoidberg.boundary module

Boundary objects that define an ‘outside’

class zoidberg.boundary.NoBoundary

No boundary, so no points outside

outside(x, y, z)

Returns True if the point is outside the boundary

Parameters:x, y, z (array_like) – Coordinates of the point(s) to check
Returns:True if point is outside boundary
Return type:bool
class zoidberg.boundary.PolygonBoundaryXZ(xarr, zarr)
outside(x, y, z)

Returns true if the given point is outside the domain

Parameters:x, y, z (array_like) – Coordinates of the point(s) to check
Returns:True if point is outside boundary
Return type:bool
class zoidberg.boundary.RectangularBoundaryXZ(xmin, xmax, zmin, zmax)
outside(x, y, z)

Returns true if the given point is outside the domain

Parameters:x, y, z (array_like) – Coordinates of the point(s) to check
Returns:True if point is outside boundary
Return type:bool

zoidberg.field module

class zoidberg.field.CurvedSlab(By=1.0, Bz=0.1, xcentre=0.0, Bzprime=1.0, Rmaj=1.0)

Represents a magnetic field in a curved slab geometry

Magnetic field in z = Bz + (x - xcentre) * Bzprime

x - Distance in radial direction [m] y - Azimuthal (toroidal) angle z - Height [m]

Parameters:
  • By (float) – Magnetic field in y direction
  • Bz (float) – Magnetic field in z at xcentre (float)
  • xcentre (float) – Reference x coordinate
  • Bzprime (float) – Rate of change of Bz with x
  • Rmaj (float) – Major radius of the slab
Bxfunc(x, z, phi)

Magnetic field in x direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:X-component of the magnetic field
Return type:ndarray
Byfunc(x, z, phi)

Magnetic field in y direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Y-component of the magnetic field
Return type:ndarray
Bzfunc(x, z, phi)

Magnetic field in z direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Z-component of the magnetic field
Return type:ndarray
Rfunc(x, z, phi)

Major radius [meters]

Returns None if in Cartesian coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The major radius
Return type:ndarray
class zoidberg.field.GEQDSK(gfile)

Read a EFIT G-Eqdsk file for a toroidal equilibrium

This generates a grid in cylindrical geometry

Parameters:gfile (str) – Name of the file to open
Bxfunc(x, z, phi)

Magnetic field in x direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:X-component of the magnetic field
Return type:ndarray
Byfunc(x, z, phi)

Magnetic field in y direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Y-component of the magnetic field
Return type:ndarray
Bzfunc(x, z, phi)

Magnetic field in z direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Z-component of the magnetic field
Return type:ndarray
Rfunc(x, z, phi)

Major radius [meters]

Returns None if in Cartesian coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The major radius
Return type:ndarray
pressure(x, z, phi)

Pressure [Pascals]

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The plasma pressure
Return type:ndarray
class zoidberg.field.MagneticField

Represents a magnetic field in either Cartesian or cylindrical geometry

This is the base class, you probably don’t want to instantiate one of these directly. Instead, create an instance of one of the subclasses.

Functions which can be overridden

  • Bxfunc = Function for magnetic field in x
  • Bzfunc = Function for magnetic field in z
  • Byfunc = Function for magnetic field in y (default = 1.)
  • Rfunc = Function for major radius. If None, y is in meters
boundary

An object with an “outside” function. See zoidberg.boundary

attributes

Contains attributes to be written to the output

Type:A dictionary of string -> function(x,z,phi)

See also

Slab
A straight field in normal Cartesian coordinates
CurvedSlab
A field in curvilinear coordinates
StraightStellarator
A rotating ellipse stellarator without curvature
VMEC
A numerical field from a VMEC equilibrium file
GEQDSK
A numerical field from an EFIT g-file
Bmag(x, z, phi)

Magnitude of the magnetic field

\[Bmag = \sqrt(B_x^2 + B_y^2 + B_z^2)\]
Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The magnitude of the magnetic field
Return type:ndarray
Bxfunc(x, z, phi)

Magnetic field in x direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:X-component of the magnetic field
Return type:ndarray
Byfunc(x, z, phi)

Magnetic field in y direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Y-component of the magnetic field
Return type:ndarray
Bzfunc(x, z, phi)

Magnetic field in z direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Z-component of the magnetic field
Return type:ndarray
Rfunc(x, z, phi)

Major radius [meters]

Returns None if in Cartesian coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The major radius
Return type:ndarray
attributes = {}
boundary = <zoidberg.boundary.NoBoundary object>
field_direction(pos, ycoord, flatten=False)

Calculate the direction of the magnetic field Returns the change in x with phi and change in z with phi

Parameters:
  • pos (ndarray) – 2-D NumPy array, with the second dimension being [x,z], with x and z in meters
  • ycoord (float) – Toroidal angle in radians if cylindrical coordinates, metres if Cartesian
  • flatten (bool, optional) – If True, return a flattened form of the vector components. This is useful for passing to FieldTracer
Returns:

(dx/dy, dz/dy)

  • = (R*Bx/Bphi, R*Bz/Bphi) if cylindrical
  • = (Bx/By, Bz/By) if Cartesian

Return type:

list of floats or ndarray

pressure(x, z, phi)

Pressure [Pascals]

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The plasma pressure
Return type:ndarray
class zoidberg.field.Slab(By=1.0, Bz=0.1, xcentre=0.0, Bzprime=1.0)

Represents a magnetic field in an infinite flat slab

Magnetic field in z = Bz + (x - xcentre) * Bzprime

Coordinates (x,y,z) assumed to be Cartesian, all in metres

Parameters:
  • By (float, optional) – Magnetic field in y direction
  • Bz (float, optional) – Magnetic field in z at xcentre
  • xcentre (float, optional) – Reference x coordinate
  • Bzprime (float, optional) – Rate of change of Bz with x
Bxfunc(x, z, phi)

Magnetic field in x direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:X-component of the magnetic field
Return type:ndarray
Byfunc(x, z, phi)

Magnetic field in y direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Y-component of the magnetic field
Return type:ndarray
Bzfunc(x, z, phi)

Magnetic field in z direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Z-component of the magnetic field
Return type:ndarray
class zoidberg.field.SmoothedMagneticField(field, grid, xboundary=None, zboundary=None)

Represents a magnetic field which is smoothed so it never leaves the boundaries of a given grid.

Parameters:
  • field (zoidberg.field.MagneticField) – A MagneticField object
  • grid (zoidberg.grid.Grid) – A Grid object
  • xboundary (int, optional) – Number of grid points in x over which the magnetic field is smoothed
  • zboundary (int, optional) – Number of grid points in x over which the magnetic field is smoothed
Bxfunc(x, z, phi)

Magnetic field in x direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:X-component of the magnetic field
Return type:ndarray
Byfunc(x, z, phi)

Not modified by smoothing

Rfunc(x, z, phi)

Major radius [meters]

Returns None if in Cartesian coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:The major radius
Return type:ndarray
smooth_field_line(xa, za)

Linearly damp the field to be parallel to the edges of the box

Should take some parameters to adjust rate of smoothing, etc.

class zoidberg.field.StraightStellarator(*args, **kwargs)

Invalid StraightStellarator, since no Sympy module.

Rather than printing an error on startup, which may be missed or ignored, this raises an exception if StraightStellarator is ever used.

class zoidberg.field.VMEC(vmec_file, ntheta=None, nzeta=None, nr=32, nz=32)

A numerical magnetic field from a VMEC equilibrium file

Parameters:
  • vmec_file (str) – Name of the VMEC file to read
  • ntheta (int, optional) – Number of theta points to use (default: use ‘mpol’ from VMEC file)n
  • nzeta (int, optional) – Number of zeta points to use (default: use ‘ntor’ from VMEC file)
  • nr (int) – Number of R points to use
  • nz (int) – Number of Z points to use
Bxfunc(x, z, phi)

Magnetic field in x direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:X-component of the magnetic field
Return type:ndarray
Byfunc(x, z, phi)

Magnetic field in y direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Y-component of the magnetic field
Return type:ndarray
Bzfunc(x, z, phi)

Magnetic field in z direction at given coordinates

Parameters:x, z, phi (array_like) – X, Z, and toroidal coordinates
Returns:Z-component of the magnetic field
Return type:ndarray
Rfunc(x, z, phi)

Major radius

cfunct(field)

VMEC DCT

read_vmec_file(vmec_file, ntheta=None, nzeta=None)

Read a VMEC equilibrium file

sfunct(field)

VMEC DST

zoidberg.fieldtracer module

class zoidberg.fieldtracer.FieldTracer(field)

A class for following magnetic field lines

Parameters:field (MagneticField) – A Zoidberg MagneticField instance
follow_field_lines(x_values, z_values, y_values, rtol=None)

Uses field_direction to follow the magnetic field from every grid (x,z) point at toroidal angle y through a change in toroidal angle dy

Parameters:
  • x_values (array_like) – Starting x coordinates
  • z_values (array_like) – Starting z coordinates
  • y_values (array_like) – y coordinates to follow the field line to. y_values[0] is the starting position
  • rtol (float, optional) – The relative tolerance to use for the integrator. If None, use the default value
Returns:

result – Field line ending coordinates

The first dimension is y, the last is (x,z). The middle dimensions are the same shape as [x|z]: [0,…] is the initial position […,0] are the x-values […,1] are the z-values If x_values is a scalar and z_values a 1D array, then result has the shape [len(y), len(z), 2], and vice-versa. If x_values and z_values are 1D arrays, then result has the shape [len(y), len(x), 2]. If x_values and z_values are 2D arrays, then result has the shape [len(y), x.shape[0], x.shape[1], 2].

Return type:

numpy.ndarray

class zoidberg.fieldtracer.FieldTracerReversible(field, rtol=1e-08, eps=1e-05, nsteps=20)

Traces magnetic field lines in a reversible way by using trapezoidal integration:

\[pos_{n+1} = pos_n + 0.5*( f(pos_n) + f(pos_{n+1}) )*dy\]

This requires a Newton iteration to solve the nonlinear set of equations for the unknown pos_{n+1}.

Parameters:
  • field (MagneticField) – A Zoidberg MagneticField instance
  • rtol (float, optional) – Tolerance applied to changes in dx**2 + dz**2
  • eps (float, optional) – Change in x,z used to calculate finite differences of magnetic field direction
  • nsteps (int, optional) – Number of sub-steps between outputs
follow_field_lines(x_values, z_values, y_values, rtol=None, eps=None, nsteps=None)

Uses field_direction to follow the magnetic field from every grid (x,z) point at toroidal angle y through a change in toroidal angle dy

Parameters:
  • x_values (array_like) – Starting x coordinates
  • z_values (array_like) – Starting z coordinates
  • y_values (array_like) – y coordinates to follow the field line to. y_values[0] is the starting position
  • rtol (float, optional) – Tolerance applied to changes in dx**2 + dz**2. If None, use the default value
  • eps (float, optional) – Change in x,z used to calculate finite differences of magnetic field direction
  • nsteps (int, optional) – Number of sub-steps between outputs
Returns:

result – Field line ending coordinates

The first dimension is y, the last is (x,z). The middle dimensions are the same shape as [x|z]: [0,…] is the initial position […,0] are the x-values […,1] are the z-values If x_values is a scalar and z_values a 1D array, then result has the shape [len(y), len(z), 2], and vice-versa. If x_values and z_values are 1D arrays, then result has the shape [len(y), len(x), 2]. If x_values and z_values are 2D arrays, then result has the shape [len(y), x.shape[0], x.shape[1], 2].

Return type:

numpy.ndarray

zoidberg.fieldtracer.trace_poincare(magnetic_field, xpos, zpos, yperiod, nplot=3, y_slices=None, revs=20, nover=20)

Trace a Poincare graph of the field lines

Does no plotting, see zoidberg.plot.plot_poincare()

Parameters:
  • magnetic_field (MagneticField) – Magnetic field object
  • xpos, zpos (array_like) – Starting X, Z locations
  • yperiod (float) – Length of period in y domain
  • nplot (int, optional) – Number of equally spaced y-slices to trace to
  • y_slices (list of ints) – List of y-slices to plot; overrides nplot
  • revs (int, optional) – Number of revolutions (times around y)
  • nover (int, optional) – Over-sample. Produced additional points in y then discards. This seems to be needed for accurate results in some cases
Returns:

coords is a Numpy array of data:

[revs, nplot, ..., R/Z]

where the first index is the revolution, second is the y slice, and last is 0 for R, 1 for Z. The middle indices are the shape of the input xpos,zpos

Return type:

coords, y_slices

zoidberg.grid module

class zoidberg.grid.Grid(poloidal_grids, ycoords, Ly, yperiodic=False, name='fci_grid')

Represents a 3D grid, consisting of a collection of poloidal grids

shape

Tuple of grid sizes (nx, ny, nz)

Type:(int, int, int)
Parameters:
  • poloidal_grids (list of PoloidalGrid) – The collection of poloidal grids to group together
  • ycoords (array_like) – The y-coordinate corresponding to each element of poloidal_grids

Examples

>>> poloidal_grids = [RectangularPoloidalGrid(5, 5, 1, 1)]
>>> ycoords = [0.0]
>>> grid = Grid(poloidal_grids, ycoords)

To iterate over the poloidal grids, and get the grids to either side:

>>> for i in range(grid.numberOfPoloidalGrids()):
...     pol, y = grid.getPoloidalGrid(i)
...     pol_next, y_next = grid.getPoloidalGrid(i+1)
...     pol_last, y_last = grid.getPoloidalGrid(i-1)

The getPoloidalGrid method ensures that y_last <= y <= y_next

getPoloidalGrid(yindex)

Returns the poloidal grid and y value at the given y index

This handles negative values and values out of range, if the domain is periodic

Parameters:yindex (int) – The desired index in y
Returns:
  • PoloidalGrid – The poloidal grid at yindex
  • float – The value of the y coordinate at yindex
metric()

Return the metric tensor, dx and dz

Returns:Dictionary containing: - dx, dy, dz: Grid spacing - gxx, gxz, gyy, gzz: Covariant components - g_xx, g_xz, g_yy, g_zz: Contravariant components
Return type:dict
numberOfPoloidalGrids()

Returns the number of poloidal grids i.e. number of points in Y

Returns:Number of poloidal grids
Return type:int
zoidberg.grid.rectangular_grid(nx, ny, nz, Lx=1.0, Ly=10.0, Lz=1.0, xcentre=0.0, zcentre=0.0, yperiodic=False)

Create a rectangular grid in (x,y,z)

Here y is along the magnetic field (typically toroidal angle), and (x,z) are in the poloidal plane

Parameters:
  • nx, ny, nz (int) – Number of points in x, y, z
  • Lx, Ly, Lz (float, optional) – Size of the domain in x, y, z
  • xcentre, zcentre (float, optional) – The middle of the domain
  • yperiodic (bool, optional) – Determines if the y direction is periodic
Returns:

A Grid representing a rectangular domain

Return type:

Grid

zoidberg.plot module

class zoidberg.plot.AnimateVectorField(X, Y, U, V)

Very basic/experimental class for animating vector fields

Transpose U, V to have dimension to animate at index 0, e.g. to animate along y, pass:

>>> AnimateVectorField(X, Y, U.transpose((1,0,2)), V.transpose((1,0,2)))
Parameters:
  • X, Y (array_like) – The X, Y coordinates
  • U, V (ndarray) – Vector components in X, Y respectively

Examples

>>> anim = AnimateVectorField(X, Y, U, V)
>>> anim.animate()
animate()
zoidberg.plot.plot_3d_field_line(magnetic_field, xpos, zpos, yperiod, cycles=20, y_res=50)

Make a 3D plot of field lines

Parameters:
  • magnetic_field (zoidberg.field.MagneticField) – Magnetic field object
  • xpos, zpos (array_like) – Starting X, Z locations
  • yperiod (float) – Length of period in y domain
  • cycles (int, optional) – Number of times to go round in y
  • y_res (int, optional) – Number of points in y in each cycle
Returns:

The matplotlib figure and axis used

Return type:

fig, ax

zoidberg.plot.plot_backward_map(grid, maps, yslice=0)

Plots the backward map from yslice to yslice-1

Parameters:
  • grid (‘zoidberg.grid.Grid`) – Grid generated by Zoidberg
  • maps (dict) – Dictionary containing the backward FCI maps
  • y_slice (int, optional) – Originating y-index to plot map from
zoidberg.plot.plot_forward_map(grid, maps, yslice=0)

Plots the forward map from yslice to yslice + 1

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • maps (dict) – Dictionary containing the forward FCI maps
  • y_slice (int, optional) – Originating y-index to plot map from
zoidberg.plot.plot_poincare(magnetic_field, xpos, zpos, yperiod, nplot=3, y_slices=None, revs=40, nover=20, interactive=False)

Plot a Poincare graph of the field lines.

Parameters:
  • magnetic_field (zoidberg.field.MagneticField) – Magnetic field object
  • xpos, zpos (array_like) – Starting X, Z locations
  • yperiod (float) – Length of period in y domain
  • nplot (int, optional) – Number of equally spaced y-slices to plot
  • y_slices (list of int, optional) – List of y-slices to plot; overrides nplot
  • revs (int, optional) – Number of revolutions (times around phi)
  • interactive (bool, optional) – If True, plots can be interacted with via the mouse: - Left-click on the plot to trace a field-line from that point - Right-click to add an additional trace - Middle-click to clear added traces
Returns:

The matplotlib figure and axis used

Return type:

fig, ax

zoidberg.plot.plot_streamlines(grid, magnetic_field, y_slice=0, width=None, **kwargs)

Plot streamlines of the magnetic field in the poloidal plane

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • y_slice (int, optional) – y-index to plot streamlines at
  • width (float, optional) – If not None, line widths are proportional to the magnitude of the magnetic_field times width
Returns:

The matplotlib figure and axis used

Return type:

fig, ax

zoidberg.poloidal_grid module

Routines for generating structured meshes on poloidal domains

Classes

RectangularPoloidalGrid
Simple rectangles in R-Z
StructuredPoloidalGrid
Curvilinear structured grids in R-Z

Functions

grid_annulus
Create a StructuredPoloidalGrid from inner and outer RZLine objects using a simple algorithm
grid_elliptic
Create a StructuredPoloidalGrid from inner and outer RZLine objects using elliptic meshing method
class zoidberg.poloidal_grid.PoloidalGrid

Represents a poloidal grid

Note: Here the 2D plane (R,Z) is labelled by (x,z) indices

nx, nz

Number of points in x and z

Type:int
R

2D Numpy array of R coordinates

Type:ndarray
Z

2D Numpy array of Z coordinates

Type:ndarray
plot(axis=None, show=True)

Plot grid using matplotlib

Parameters:
  • axis (matplotlib axis, optional) – A matplotlib axis to plot on. By default a new figure is created
  • show (bool, optional) – Calls plt.show() at the end
Returns:

The matplotlib axis that was used

Return type:

axis

class zoidberg.poloidal_grid.RectangularPoloidalGrid(nx, nz, Lx, Lz, Rcentre=0.0, Zcentre=0.0, MXG=2)

Represents a poloidal grid consisting of a rectangular domain

Note: Here the 2D plane (R,Z) is labelled by (x,z) indices

nx, nz

Number of points in x and z

Type:int
R

2D Numpy array of R coordinates

Type:ndarray
Z

2D Numpy array of Z coordinates

Type:ndarray
Parameters:
  • nx (int) – Number of points in major radius (including boundaries)
  • nz (int) – Number of points in height (including boundaries)
  • Lx (float) – Radial domain size [m]
  • Lz (float) – Vertical domain size [m]
  • Rcentre (float, optional) – Coordinate at the middle of the domain
  • Zcentre (float, optional) – Coordinate at the middle of the domain
  • MXG (int, optional) – Number of guard cells in X. The boundary is put half-way between the guard cell and the domain
findIndex(R, Z)

Finds the (x,z) index corresponding to the given (R,Z) coordinate

Parameters:R, Z (array_like) – Locations to find indices for
Returns:x, z – Index as a float, same shape as R,Z
Return type:(ndarray, ndarray)
getCoordinate(xind, zind, dx=0, dz=0)

Get coordinates (R,Z) at given (xind,zind) index

Parameters:
  • xind, zind (array_like) – Indices in X and Z. These should be the same shape
  • dx (int, optional) – Order of x derivative
  • dz (int, optional) – Order of z derivative
Returns:

R, Z – Locations of point or derivatives of R,Z with respect to indices if dx,dz != 0

Return type:

(ndarray, ndarray)

metric()

Return the metric tensor, dx and dz

For this rectangular grid the metric is the identity

Returns:Dictionary containing: - dx, dz: Grid spacing - gxx, gxz, gzz: Covariant components - g_xx, g_xz, g_zz: Contravariant components
Return type:dict
class zoidberg.poloidal_grid.StructuredPoloidalGrid(R, Z)

Represents a structured poloidal grid in R-Z

nx, nz

Number of points in x and z

Type:int
R

2D Numpy array of R coordinates

Type:ndarray
Z

2D Numpy array of Z coordinates

Type:ndarray
Parameters:

R, Z (ndarray) – 2D Numpy arrays of R,Z points

Note

R,Z are not copied, so these arrays should not be modified afterwards

findIndex(R, Z, tol=1e-10, show=False)

Finds the (x, z) index corresponding to the given (R, Z) coordinate

Parameters:
  • R, Z (array_like) – Locations. Can be scalar or array, must be the same shape
  • tol (float, optional) – Maximum tolerance on the square distance
Returns:

x, z – Index as a float, same shape as R, Z

Return type:

(ndarray, ndarray)

getCoordinate(xind, zind, dx=0, dz=0)

Get coordinates (R, Z) at given (xind, zind) index

Parameters:
  • xind, zind (array_like) – Indices in X and Z. These should be the same shape
  • dx (int, optional) – Order of x derivative
  • dz (int, optional) – Order of z derivative
Returns:

R, Z – Locations of point or derivatives of R,Z with respect to indices if dx,dz != 0

Return type:

(ndarray, ndarray)

metric()

Return the metric tensor, dx and dz

Returns:Dictionary containing: - dx, dz: Grid spacing - gxx, gxz, gzz: Covariant components - g_xx, g_xz, g_zz: Contravariant components
Return type:dict
zoidberg.poloidal_grid.grid_annulus(inner, outer, nx, nz, show=True, return_coords=False)

Grid an annular region, given inner and outer boundaries both of which are RZline objects

This is a very simple algorithm which just draws straight lines between inner and outer boundaries.

Parameters:
  • inner, outer (RZline) – Inner and outer boundaries of the domain
  • nx (int) – The required radial resolution, including boundaries
  • nz (int) – The required poloidal resolution
  • show (bool, optional) – If True, plot the resulting grid
  • return_coords (bool, optional) – If True, return the R, Z coordinates of the grid points, instead of a StructuredPoloidalGrid
Returns:

A grid of the region

Return type:

StructuredPoloidalGrid

zoidberg.poloidal_grid.grid_elliptic(inner, outer, nx, nz, show=False, tol=1e-10, align=True, restrict_size=20, restrict_factor=2, return_coords=False)

Create a structured grid between inner and outer boundaries using elliptic method

Coordinates x = x(R, Z) and z = z(R,Z) obey an elliptic equation:

\[ \begin{align}\begin{aligned}d^2x/dR^2 + d^2x/dZ^2 = 0\\d^2z/dR^2 + d^2z/dZ^2 = 0\end{aligned}\end{align} \]

where here x is in in the domain (0, 1) and z in (0, 2pi)

The above equations are inverted, giving:

\[ \begin{align}\begin{aligned}a*R_xx - 2*b*R_xz + c*R_zz = 0\\a*Z_xx - 2*b*Z_xz + c*Z_zz = 0\end{aligned}\end{align} \]

where

\[ \begin{align}\begin{aligned}a &= R_z^2 + Z_z^2\\b &= R_z*R_x + Z_x*Z_z\\c &= R_x^2 + Z_x^2\end{aligned}\end{align} \]

This is a nonlinear system of equations which is solved iteratively.

Parameters:
  • inner, outer (RZline) – Inner and outer boundaries of the domain
  • nx (int) – The required radial resolution, including boundaries
  • nz (int) – The required poloidal resolution
  • show (bool, optional) – Display plots of intermediate results
  • tol (float, optional) – Controls when iteration stops
  • align (bool, optional) – Attempt to align the inner and outer boundaries
  • restrict_size (int, optional) – The size (nx or nz) above which the grid is coarsened
  • restrict_factor (int, optional) – The factor by which the grid is divided if coarsened
  • return_coords (bool, optional) – If True, return the R, Z coordinates of the grid points, instead of a StructuredPoloidalGrid
Returns:

  • If return_coords is true, returns R,Z as arrays.
  • If return_coords is false, returns a StructuredPoloidalGrid object

References

https://www.nada.kth.se/kurser/kth/2D1263/l2.pdf https://en.wikipedia.org/wiki/Principles_of_grid_generation

zoidberg.progress module

zoidberg.progress.update_progress(progress, barLength=10, ascii=False, **kwargs)

Displays or updates a console progress bar

Accepts a float between 0 and 1. Any int will be converted to a float. A value under 0 represents a ‘halt’. A value at 1 or bigger represents 100%

Parameters:
  • progress (float) – Number between 0 and 1
  • barLength (int, optional) – Length of the progress bar
  • ascii (bool, optional) – If True, use ‘#’ as the progress indicator, otherwise use a Unicode character (the default)

zoidberg.rzline module

Routines and classes for representing periodic lines in R-Z poloidal planes

class zoidberg.rzline.RZline(r, z, anticlockwise=True)

Represents (R,Z) coordinates of a periodic line

R

Major radius [m]

Type:array_like
Z

Height [m]

Type:array_like
theta

Angle variable [radians]

R, Z and theta all have the same length

Type:array_like
Parameters:
  • r, z (array_like) – 1D arrays of the major radius (r) and height (z) which are of the same length. A periodic domain is assumed, so the last point connects to the first.
  • anticlockwise (bool, optional) – Ensure that the line goes anticlockwise in the R-Z plane (positive theta)
  • Note that the last point in (r,z) arrays should not be the same
  • as the first point. The (r,z) points are in [0,2pi)
  • The input r,z points will be reordered, so that the
  • theta angle goes anticlockwise in the R-Z plane
Rvalue(theta=None, deriv=0)

Calculate the value of R at given theta locations

Parameters:
  • theta (array_like, optional) – Theta locations to find R at. If None (default), use the values of theta stored in the instance
  • deriv (int, optional) – The order of derivative to compute (default is just the R value)
Returns:

Value of R at each input theta point

Return type:

ndarray

Zvalue(theta=None, deriv=0)

Calculate the value of Z at given theta locations

Parameters:
  • theta (array_like, optional) – Theta locations to find Z at. If None (default), use the values of theta stored in the instance
  • deriv (int, optional) – The order of derivative to compute (default is just the Z value)
Returns:

Value of Z at each input theta point

Return type:

ndarray

closestPoint(R, Z, niter=3, subdivide=20)

Find the closest point on the curve to the given (R,Z) point

Parameters:
  • R, Z (float) – The input R, Z point
  • niter (int, optional) – How many iterations to use
Returns:

The value of theta (angle)

Return type:

float

distance(sample=20)

Integrates the distance along the line.

Parameters:sample (int, optional) – Number of samples to take per point
Returns:
  • An array one longer than theta. The first element is zero,
  • and the last element is the total distance around the loop
equallySpaced(n=None)

Returns a new RZline which has a theta uniform in distance along the line

Parameters:n (int, optional) – Number of points. Default is the same as the current line
Returns:A new RZline based on this instance, but with uniform theta-spacing
Return type:RZline
plot(axis=None, show=True)

Plot the RZline, either on the given axis or a new figure

Parameters:
  • axis (matplotlib axis, optional) – A matplotlib axis to plot on. By default a new figure is created
  • show (bool, optional) – Calls plt.show() at the end
Returns:

The matplotlib axis that was used

Return type:

axis

position(theta=None)

Calculate the value of both R, Z at given theta locations

Parameters:theta (array_like, optional) – Theta locations to find R, Z at. If None (default), use the values of theta stored in the instance
Returns:R, Z – Value of R, Z at each input theta point
Return type:(ndarray, ndarray)
positionPolygon(theta=None)

Calculates (R,Z) position at given theta angle by joining points by straight lines rather than a spline. This avoids the overshoots which can occur with splines.

Parameters:theta (array_like, optional) – Theta locations to find R, Z at. If None (default), use the values of theta stored in the instance
Returns:R, Z – Value of R, Z at each input theta point
Return type:(ndarray, ndarray)
zoidberg.rzline.circle(R0=1.0, r=0.5, n=20)

Creates a pair of RZline objects, for inner and outer boundaries

Parameters:
  • R0 (float, optional) – Centre point of the circle
  • r (float, optional) – Radius of the circle
  • n (int, optional) – Number of points to use in the boundary
Returns:

A circular RZline

Return type:

RZline

zoidberg.rzline.line_from_points(rarray, zarray, show=False)

Find a periodic line which goes through the given (r,z) points

This function starts at a point, and finds the nearest neighbour which is not already in the line

Parameters:rarray, zarray (array_like) – R, Z coordinates. These arrays should be the same length
Returns:An RZline object representing a periodic line
Return type:RZline
zoidberg.rzline.line_from_points_poly(rarray, zarray, show=False)

Find a periodic line which goes through the given (r,z) points

This function starts with a triangle, then adds points one by one, inserting into the polygon along the nearest edge

Parameters:rarray, zarray (array_like) – R, Z coordinates. These arrays should be the same length
Returns:An RZline object representing a periodic line
Return type:RZline
zoidberg.rzline.shaped_line(R0=3.0, a=1.0, elong=0.0, triang=0.0, indent=0.0, n=20)

Parametrisation of plasma shape from J. Manickam, Nucl. Fusion 24 595 (1984)

Parameters:
  • R0 (float, optional) – Major radius
  • a (float, optional) – Minor radius
  • elong (float, optional) – Elongation, 0 for a circle
  • triang (float, optional) – Triangularity, 0 for a circle
  • indent (float, optional) – Indentation, 0 for a circle
Returns:

An RZline matching the given parameterisation

Return type:

RZline

zoidberg.test_field module

zoidberg.test_field.test_curved_slab()
zoidberg.test_field.test_slab()

zoidberg.test_fieldtracer module

zoidberg.test_fieldtracer.test_FieldTracerReversible_slab()
zoidberg.test_fieldtracer.test_poincare()
zoidberg.test_fieldtracer.test_slab()

zoidberg.test_grid module

zoidberg.test_grid.test_getPoloidalGrid()
zoidberg.test_grid.test_getPoloidalGrid_periodic()

zoidberg.test_poloidal_grid module

zoidberg.test_poloidal_grid.test_out_of_domain()

zoidberg.test_rzline module

zoidberg.test_rzline.test_circular_boundaries()
zoidberg.test_rzline.test_distance()
zoidberg.test_rzline.test_line_from_points()
zoidberg.test_rzline.test_order_by_distance()

zoidberg.test_zoidberg module

zoidberg.test_zoidberg.test_make_maps_slab()
zoidberg.test_zoidberg.test_make_maps_straight_stellarator()

zoidberg.zoidberg module

zoidberg.zoidberg.fci_to_vtk(infile, outfile, scale=5)
zoidberg.zoidberg.make_maps(grid, magnetic_field, nslice=1, quiet=False, **kwargs)

Make the forward and backward FCI maps

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • nslice (int) – Number of parallel slices in each direction
  • quiet (bool) – Don’t display progress bar
  • kwargs – Optional arguments for field line tracing, etc.
Returns:

Dictionary containing the forward/backward field line maps

Return type:

dict

zoidberg.zoidberg.make_surfaces(grid, magnetic_field, nsurfaces=10, revs=100)

Essentially interpolate a poincare plot onto the grid mesh

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • nsurfaces (int, optional) – Number of surfaces to interpolate to [10]
  • revs (int, optional) – Number of points on each surface [100]
Returns:

Array of psuedo-psi on the grid mesh

Return type:

surfaces

zoidberg.zoidberg.parallel_slice_field_name(field, offset)

Form a unique, backwards-compatible name for field at a given offset

Parameters:
  • field (str) – Name of the field to convert
  • offset (int) – Parallel slice offset
zoidberg.zoidberg.upscale(field, maps, upscale_factor=4, quiet=True)

Increase the resolution in y of field along the FCI maps.

First, interpolate onto the (forward) field line end points, as in normal FCI technique. Then interpolate between start and end points. We also need to interpolate the xt_primes and zt_primes. This gives a cloud of points along the field lines, which we can finally interpolate back onto a regular grid.

Parameters:
  • field (array_like) – 3D field to be upscaled
  • maps (dict) – Zoidberg field line maps
  • upscale_factor (int, optional) – Factor to increase resolution by [4]
  • quiet (bool, optional) – Don’t show progress bar [True]
Returns:

  • Field with y-resolution increased *upscale_factor times. Shape is*
  • (nx, upscale_factor*ny, nz).

zoidberg.zoidberg.write_Bfield_to_vtk(grid, magnetic_field, scale=5, vtkfile='fci_zoidberg', psi=True)

Write the magnetic field to a VTK file

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • scale (int, optional) – Factor to scale x, z dimensions by [5]
  • vtkfile (str, optional) – Output filename without extension [“fci_zoidberg”]
  • psi (bool, optional) – Write psi?
Returns:

Return type:

path - Full path to vtkfile

zoidberg.zoidberg.write_maps(grid, magnetic_field, maps, gridfile='fci.grid.nc', new_names=False, metric2d=True, format='NETCDF3_64BIT', quiet=False)

Write FCI maps to BOUT++ grid file

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • maps (dict) – Dictionary of FCI maps
  • gridfile (str, optional) – Output filename
  • new_names (bool, optional) – Write “g_yy” rather than “g_22”
  • metric2d (bool, optional) – Output only 2D metrics
  • format (str, optional) – Specifies file format to use, passed to boutdata.DataFile
  • quiet (bool, optional) – Don’t warn about 2D metrics
Returns:

Return type:

Writes the following variables to the grid file

Module contents

zoidberg.make_maps(grid, magnetic_field, nslice=1, quiet=False, **kwargs)

Make the forward and backward FCI maps

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • nslice (int) – Number of parallel slices in each direction
  • quiet (bool) – Don’t display progress bar
  • kwargs – Optional arguments for field line tracing, etc.
Returns:

Dictionary containing the forward/backward field line maps

Return type:

dict

zoidberg.write_maps(grid, magnetic_field, maps, gridfile='fci.grid.nc', new_names=False, metric2d=True, format='NETCDF3_64BIT', quiet=False)

Write FCI maps to BOUT++ grid file

Parameters:
  • grid (zoidberg.grid.Grid) – Grid generated by Zoidberg
  • magnetic_field (zoidberg.field.MagneticField) – Zoidberg magnetic field object
  • maps (dict) – Dictionary of FCI maps
  • gridfile (str, optional) – Output filename
  • new_names (bool, optional) – Write “g_yy” rather than “g_22”
  • metric2d (bool, optional) – Output only 2D metrics
  • format (str, optional) – Specifies file format to use, passed to boutdata.DataFile
  • quiet (bool, optional) – Don’t warn about 2D metrics
Returns:

Return type:

Writes the following variables to the grid file