ascotpy

Python interface to libascot.so.

class a5py.ascotpy.Ascotpy

Class with methods to initialize and access the data via Python.

This class will be inherited by Ascot class. Here we hide all the dirty implementation of the libascot.so interface in private methods and attributes so that the Ascot front is clean for the users.

Attributes:

_sim

Simulation data struct.

_nmrk

Number of markers currently in the marker array.

_inistate

Marker input array for interactive simulations.

_endstate

Marker output array for interactive simulations.

_diag_occupied

Flag indicating if diagnostics output contains data.

_mute

Mute output from libascot.so: “yes” - stdout and stderr both muted, “no” - output is not muted, “err” - stderr is printed.

DUMMY_QID = b''

Flag to use in _sim.qid_* to indicate that the input is not initialized.

get_plasmaquantities()

Return species present in plasma input.

input_check_bfield_seam(qnt, r=None, z=None, axes=None, centremostpercentage=10, diff=True, forceplot=False, n_centremost_elements_for_fit=None, plot_spline_interpolated=False)

Tries to check if an (R,z) plane is missing near the “seam” of the bfield input.

This function is intended to be used to verify that the correct number of (R,z) planes is provided when using a stellarator bfield input. To recap, psi and Bfield components should be given at the points np.linspace(phimin, phimax, nphi+1)[:-1]. (It is possible to independently specify phimin, phimax and nphi for both psi and the bfield components.)

If an (R,z) plane might be missing, plots qnt as a function of phi, evaluated only at the grid points where the input is given.

Parameters:

qntstr

One of the following: {“psi”, “br”, “bphi”, “bz”}

rfloat, optional

R coordinate where data is evaluated. (The functions uses the R value of the input that is closest to this requested one) [m]

zfloat, optional

z coordinate where data is evaluated. (The functions uses the z value of the input that is closest to this requested one) [m]

axesAxes, optional

The axes where figure is plotted or otherwise new figure is created.

centremostpercentagefloat, optional

Percentage of the grid closest to the requested r and z. (Default is 10, i.e. 5 % of the grid on either side)

diffbool, optional

If True, plots np.diff(qnt)/np.diff(phi) instead of qnt. A possible kink may be more clearly visible in the diff(qnt) plot.

forceplotbool, optional

If True, plots the requested quantity even if no (R,z) plane is deemed missing.

n_centremost_elements_for_fitint, optional

Number of elements from the centre of the quantity array to use for least square sum fitting a polynomial. Default is None, in which case the whole array, consisting of ‘centremostpercentage’, is used.

plot_spline_interpolatedbool, optional

If True, plots the spline interpolated version of qnt.

input_eval(r, phi, z, t, *qnt, grid=False)

Evaluate input quantities at given coordinates.

This method uses ASCOT5 C-routines for evaluating inputs exactly as they are evaluated during a simulation. See input_eval_list() for a complete list of available input and derived quantities.

Parameters:

rarray_like

R coordinates where data is evaluated.

phiarray_like

phi coordinates where data is evaluated.

zarray_like

z coordinates where data is evaluated.

tarray_like

Time coordinates where data is evaluated.

*qntstr

Names of evaluated quantities.

gridbool, optional

Treat input coordinates as abscissae and return the evaluated quantities on a grid instead.

Returns:

valarray_like (n,) or (nr,nphi,nz,nt)

Quantity evaluated in each queried point.

NaN is returned at those points where the evaluation failed.

4D matrix is returned when grid=``True``.

Raises:

AssertionError

If required data has not been initialized.

RuntimeError

If evaluation in libascot.so failed.

input_eval_list(show=True)

Return quantities that can be evaluated from inputs.

The returned quantities are those that can be evaluated from the inputs with input_eval method. Units are not listed as that function returns quantities with units included.

The number ion species in plasma input is not fixed so “ni” refers to i’th plasma species as listed by input_getplasmaspecies (index starts from 1).

Parameters:

showbool, optional

Print output.

Returns:

quantitiesdict [str, str]

Name of each quantity and a short description.

input_eval_orbitresonance(rhogrid, xigrid, egrid, species, plot=False, n=1, p=0, omega=0.0, axes=None, cax=None)

Evaluate the resonance condition between particles and perturbations.

Orbit resonance refers to the phenomenon where the frequency of a particle’s orbit in a magnetic field matches the frequency of a perturbation, such as a wave or an external perturbation. This resonance can lead to increased particle transport.

The resonance condition for orbit resonance is given by:

f_{tor} / f_{pol} = m / n

where: - f_{tor} is the toroidal orbit frequency of the particle, - f_{pol} is the poloidal orbit frequency of the particle, - n is the toroidal mode number of the perturbation. - m is an arbitrary integer.

This function can evaluate and plot the resonance in 1D or 2D space. Note that one of the (‘rhogrid’, ‘xigrid’, ‘egrid’) must be a scalar.

Parameters:

rhogridarray_like, (nrho)

Value of the radial coordinate of radial grid for which the resonance is evaluated.

xigridfloat or array_like, (nxi)

Value of the pitch of pitch grid for which the resonance is evaluated.

egridfloat or array_like, (ne)

Value of the energy of energy grid for which the resonance is evaluated.

speciesstr

Name of the particle species.

plotbool, optional

If True, the results are not returned but visualized instead.

resonancefloat or array_like, optional

Resonance(s) to be indicated on the plot e.g. m/n=2/1.

omegafloat, optional

Perturbation mode frequency.

The default value assumes stationary perturbation.

axesAxes, optional

The axes where figure is plotted or otherwise new figure is created.

caxAxes, optional

The color bar axes or otherwise taken from the main axes.

Returns:

torfreqarray_like (nrho, nxi) or (nrho, ne) or (nxi, ne)

Evaluated toroidal orbit frequencies if ‘plot’ = False.

polfreqarray_like (nrho, nxi) or (nrho, ne) or (nxi, ne)

Evaluated poloidal orbit frequencies if ‘plot’ = False.

qeffarray_like (nrho, nxi) or (nrho, ne) or (nxi, ne)

Evaluated effective qfactor = torfreq / polfreq if ‘plot’ = False.

input_eval_ripple(rgrid, zgrid, rlarmor, rho0=unyt_quantity(1., 'dimensionless'), theta0=unyt_quantity(0, 'rad'), t=unyt_quantity(0, 's'), nphi=362, plot=True, axes1=None, axes2=None)

Evaluate various ripple quantities.

The variation of toroidal field strength as a function of toroidal angle is called the toroidal field (TF) ripple. The TF ripple can lead to ripple-induced losses provided that i) the particle is trapped in a ripple well in which case it is promptly lost or ii) the ripple causes particle motion to become stochastic which leads to radial diffusion.

This function evaluates the quantities relevant for ripple transport. Note that deltacrit and ripplewell are approximate quantities based on analytical estimates.

Parameters:

rgridarray_like

R points where delta, deltacrit, and ripplewell are evaluated.

zgridarray_like

z points where delta, deltacrit, and ripplewell are evaluated.

rlarmorfloat or tuple (str, int, float)

Larmor radius of the test particle or tuple with (species, charge, energy).

rho0float, optional

Radial coordinate where amplitude is evaluated.

theta0float, optional

Poloidal coordinate where amplitude is evaluated.

tfloat, optional

Time instance for time-dependent fields.

nphiint, optional

Number of toroidal grid points in ripple evaluation.

There is probably never a need to change this number.

plotbool, optional

If True, the evaluated quantities are also plotted.

axes1, optional

Axes for plotting amplitude if plot = True.

axes2, optional

Axes for plotting delta, deltacrit, and ripplewell if plot = True.

Returns:

phigridarray_like, (nphi)

Toroidal grid [-2, 362] where amplitude is given.

This grid is not limited to [0,360] as usual so that one can also use it to spot any discontinuities in field if the input was poorly constructed.

amplitudearray_like, (nphi,)

Variation in toroidal field strength, i.e. Bphi(phi) - Bphi_axisymmetric.

deltaarray_like, (nr,nz)

Ripple magnitude defined as (Bmax-Bmin) / (Bmax+Bmin), where Bmax and Bmin are toroidal extrema of Bphi.

deltacritarray_like, (nr,nz)

Critical ripple magnitude divided by delta.

Particles are subject to stochastic ripple diffusion when deltacrit < 1.

ripplewellarray_like, (nr,nz)

Ripple well parameter.

Particles are bound to be ripple-trapped when ripplewell < 1. However, this does not necessarily mean that the particle is lost as it may escape the well when ripple-trapping causes it to drift outwards. As particles that are ripple-trapped drift vertically, a good indication that a trapped particle will be lost if the ripple well region covers everything vertically below/above that point.

input_eval_safetyfactor(rho, nth=10000, return_ftrap=False)

Evaluate safety factor and associated quantities.

This function most likely works for tokamaks only.

Parameters:

rhoarray_like, (n,)

Radial positions where the safety factor is evaluated.

Note that evaluating it near the axis or exactly at the separatrix may yield mess.

return_bminmaxbool, optional

Return also minimum and maximum values of B along the flux surface.

Returns:

qprofarray_like, (n,)

Safety factor.

Iprofarray_like, (n,)

Toroidal current term which, when multiplied with 2pi/mu0, gives the enclosed toroidal current.

gprofarray_like, (n,)

This is just R * Bphi which is constant since Bphi ~ 1/R.

ftraparray_like, (n,)

Trapped particle fraction if return_ftrap is true.

input_free_rfof()

input_init_rfof()

input_initialized()

Get inputs that are currently initialized.

Returns:

outputsdict [str, str]

Name of the input parent groups that are initialized and QID of the initialized input.

input_plotalongaxis(qnt, nphi=None, phimin=0, phimax=360, t=unyt_quantity(0., 's'), axes=None)

Plot qnt on the magnetic axis as a function of phi.

By default plots qnt in the phigrid of the axis inputs. However, if nphi is given, qnt is plotted against np.linspace(phimin, phimax, nphi)[:-1]

Parameters:

qntstr

Name of the quantity to be plotted along the magnetic axis.

nphiint, optional

number of phi points for plotting

phiminfloat, optional

Minimum phi value at which qnt is plotted.

phimaxfloat, optional

“Maximum” phivalua for the plot. See the comment above parameters.

input_plotcontoursurface3d(qnt, isovalue, xmin, xmax, ymin, ymax, zmin, zmax, pointspermeterx=10, pointspermetery=10, pointspermeterz=10, t=unyt_quantity(0, 's'), max_rho=1.0, ax=None, c=0.135, delta=None, tri_lc=None, tri_lw=0.4, opacity=0.8, reltol=0.01, axscatter=None, surfcolor='blue')

Plot a 3D contour surface of the specified quantity at a given isovalue.

If there are holes in you surface, try lowering input parameter c or increasing reltol.

Parameters:

qntstr

Name of the quantity whose contour surface is drawn (e.g., “bnorm”).

isovaluefloat

Value of the quantity at which to generate the isosurface.

xmin, xmaxfloat

Minimum and maximum bounds along the x-axis.

ymin, ymaxfloat

Minimum and maximum bounds along the y-axis.

zmin, zmaxfloat

Minimum and maximum bounds along the z-axis.

pointspermeterx, pointspermetery, pointspermeterzint, optional

Number of grid points along the x, y, and z directions per meter. Ideally, use uniform grid spacing in all directions.

tfloat or unyt_quantity, optional

Time coordinate at which the data is evaluated. Default is 0 seconds.

max_rhofloat, optional

Maximum rho until (surface outside this is not drawn)

axmatplotlib.axes.Axes, optional

Matplotlib 3D axes to plot on. If None, a new figure and axes are created.

cfloat, optional

Scaling parameter for the grid spacing or contour thresholding. Default is 0.135.

deltafloat or None, optional

Characteristic grid spacing. Ideally, use uniform grid of points.

tri_lcstr or None, optional

Color of the mesh triangle edges. If None, edges are invisible.

tri_lwfloat, optional

Width of the mesh triangle edges. Default is 0.4.

opacityfloat, optional

Opacity of the drawn surface. (0: transparent, 1: opaque)

reltolfloat, optional

Relative tolerance for isosurface extraction. Default is 0.01.

axscatter3D axes, optional

If given, plots the grid points of the surface on these axes. This can be useful especially when verifying that the the plotted surface is sensible.

input_plotfluxsurface3d(rho, phi, theta, qnt=None, t=unyt_quantity(0, 's'), clim=None, cmap=None, axes=None, cax=None, meshalpha=0.6, tri_lc=None, tri_lw=0.4, opacity=0.8)

Plot a flux surface or, given qnt, a quantity on a flux surface in 3D.

To plot quantity on a logarithmic scale (base 10), add “log” in the name of the quantity e.g. “log ne”.

Parameters:

rhofloat

Rho co-ordinate of the flux surface.

phiarray_like (nphi,)

phi abscissa where data is evaluated and plotted.

thetaarray_like (ntheta,)

theta abscissa where data is evaluated and plotted.

qntstr, optional

Name of the plotted quantity. E.g. “bnorm”

tfloat, optional

Time coordinate where data is evaluated.

climlist [float, float], optional

Minimum and maximum color range values.

cmapstr, optional

Colormap.

axesAxes, optional

The axes where figure is plotted or otherwise new figure is created.

caxAxes, optional

The color bar axes or otherwise taken from the main axes.

meshalpha: float, optional

A parameter used when generating the triangles for the surface. For more info, see “alphashape” and “alpha parameter”.

tri_lcstr, optional

Color of the mesh triangle edges. By default, edges are invisible.

tri_lwfloat, optional

Width of the mesh triangle edges. To enable edges, set also tri_lc.

opacityfloat, optional

Opacity of the drawn surface (0: transparent, 1: opaque).

input_plotphitheta(phi, theta, rho, qnt, t=unyt_quantity(0, 's'), clim=None, cmap=None, axes=None, cax=None, drawcontours=False, contourvalues=None, contourcolors=None, contourlinestyles='-', contourlinewidths=None, labelfontsize=10, drawbfielddir=False, bfieldlinecolor='black', bfieldlinedensity=0.5, blinewidth=1.0, blinelength=50.0)

Plot input quantity on a (phi, theta) surface at given rho.

With drawcontours=True, draws also the contourlines. To plot quantity on a logarithmic scale (base 10), add “log” in the name of the quantity e.g. “log ne”.

Parameters:

phiarray_like (nphi,)

phi abscissa where data is evaluated and plotted.

thetaarray_like (ntheta,)

theta abscissa where data is evaluated and plotted.

rhofloat

Rho co-ordinate of the flux surface.

qntstr

Name of the plotted quantity.

tfloat, optional

Time coordinate where data is evaluated.

climlist [float, float], optional

Minimum and maximum color range values.

cmapstr, optional

Colormap.

axesAxes, optional

The axes where figure is plotted or otherwise new figure is created.

caxAxes, optional

The color bar axes or otherwise taken from the main axes.

drawcontoursbool, optional

If true, draws contourlines on top of the heatmap.

contourvalueslist [float], optional

Which contour lines to draw. Leave empty for automatic selection

contourcolorsstr or list[str], optional

Colors of the contour lines. By giving a list, one can match the colors with the contourvalues (see above).

contourlinestylesstr or list[str], optional

Linestyle of the contour lines – as the reader might have guessed. By giving a list, one can match the styles with the contourvalues ( see above).

contourlinewidthsfloat or list[float], optional

You can probably guess. By giving a list, one can match the widths with the contourvalues (see above).

labelfontsizefloat, optional

Font size of the labels that indicate contour values.

drawbfielddirbool, optional

Draw the direction of the bfield on top of the figure as streamlines. This makes sense only in a tokamak geometry since the poloidal plane is a contant phi plane. If drawfielddir true, one must have theta and phi go from 0 to 360.

bfieldlinecolorstr, optional

Color of the bfield lines if drawbfielddir=True.

bfieldlinedensityfloat, optional

Density of the bfield lines.

blinewidthfloat, optional

Width of the bfield lines.

blinelengthfloat, optional

Length of the bfield lines.

input_plotrhocontour(rho=1.0, phi=unyt_quantity(0, 'degree'), t=unyt_quantity(0, 's'), ntheta=180, axes=None, **kwargs)

Plot contours of constant rho.

Calling this without arguments plots the separatrix.

Parameters:

rhofloat or array_like, optional

Value(s) of square root of normalized poloidal flux whose contours are plotted.

phifloat, optional

Toroidal angle at which contours are plotted.

tfloat, optional

Time instance when the contours are plotted.

nthetaint, optional

Poloidal resolution of the contours.

axesAxes, optional

The axes where figure is plotted or otherwise new figure is created.

**kwargs

Arguments passed to matplotlib.pyplot.plot().

input_plotrz(r, z, qnt, phi=unyt_quantity(0, 'degree'), t=unyt_quantity(0, 's'), rhomax=None, clim=None, cmap=None, axes=None, cax=None, drawcontours=False, contourvalues=None, contourcolors=None, contourlinestyles='-', contourlinewidths=None, labelfontsize=10)

Plot input quantity on a (R, z) plane at given coordinates.

With drawcontours=True, draws also the contourlines. To plot quantity on a logarithmic scale (base 10), add “log” in the name of the quantity e.g. “log ne”.

Parameters:

rarray_like (nr,)

R abscissa where data is evaluated and plotted.

zarray_like (nz,)

z abscissa where data is evaluated and plotted.

qntstr

Name of the plotted quantity.

phifloat, optional

Toroidal coordinate where data is evaluated.

tfloat, optional

Time coordinate where data is evaluated.

rhomaxfloat, optional

Do not plot (R, z) points for which rho > rhomax

climlist [float, float], optional

Minimum and maximum color range values.

cmapstr, optional

Colormap.

axesAxes, optional

The axes where figure is plotted or otherwise new figure is created.

caxAxes, optional

The color bar axes or otherwise taken from the main axes.

drawcontoursbool, optional

If true, draws contourlines on top of the heatmap.

contourvalueslist [float], optional

Which contour lines to draw. Leave empty for automatic selection

contourcolorsstr or list[str], optional

Colors of the contour lines. By giving a list, one can match the colors with the contourvalues (see above).

contourlinestylesstr or list[str], optional

Linestyle of the contour lines – as the reader might have guessed. By giving a list, one can match the styles with the contourvalues ( see above).

contourlinewidthsfloat or list[float], optional

You can probably guess. By giving a list, one can match the widths with the contourvalues (see above).

labelfontsizefloat, optional

Font size of the labels that indicate contour values.

input_rhovolume(nrho=10, ntheta=360, nphi=10, method='mc', tol=0.1, t=unyt_quantity(0, 's'), return_area=False, return_coords=False, minrho=0, maxrho=1, minphi=0, maxphi=360, mintheta=0, maxtheta=360)

Evaluate flux surface volumes.

Parameters:

nrhoint, optional

Number of radial grid edges between [0, 1].

nthetaint, optional

Number of poloidal grid edges.

nphiint, optional

Number of toroidal grid edges.

method{“prism”, “mc”}, optional

Method used to evaluate the volumes.

“prism” divides the volume into structured 3D grid consisting of triangular prims where the triangular plane is poloidal and one of the vertices is at magnetic axis while the other two are along a given flux surface. The volume is the estimated as a sum of these prisms.

“mc” divides the volume into toroidal sectors of same size, and on each sector fins the separatrix bounding box. Bounding box is assumed to be fixed along the sector and markers are randomly sampled within it. The volumes are then estimated via the Monte Carlo method.

tolfloat, optional

If relative difference in volume in subsequent iterations falls below this value, the algorithm finishes.

tfloat, optional

Time slice when the volumes are computed.

return_areabool, optional

Return also the area on (R, z) plane.

return_coordsbool, optional

Return also the (R, phi, z) coordinates for grid center points.

minrhofloat, optional

minimum value of rho

maxrhofloat, optional

maximum value of rho

minphifloat, optional

minimum value of phi (in degrees)

maxphifloat, optional

maximum value of phi (in degrees)

minthetafloat, optional

minimum value of theta (in degrees)

maxthetafloat, optional

maximum value of theta (in degrees)

Returns:

volumearray_like (nrho-1, ntheta-1, nphi-1)

Volume of each bin.

areaarray_like (nrho-1, ntheta-1, nphi-1)

Area of the poloidal cross section of each bin if return_area=True.

rarray_like (nrho-1, ntheta-1, nphi-1)

R coordinates of bin centers if return_coords=True.

parray_like (nrho-1, ntheta-1, nphi-1)

phi coordinates of bin centers if return_coords=True.

zarray_like (nrho-1, ntheta-1, nphi-1)

z coordinates of bin centers if return_coords=True.