Base classes

Base classes for mapped domains (single patch).

class struphy.geometry.base.Domain(Nel: tuple[int] | None = None, p: tuple[int] | None = None, spl_kind: tuple[bool] | None = None)

Bases: object

Base class for mapped domains (single patch).

The (physical) domain \(\Omega \subset \mathbb R^3\) is an open subset of \(\mathbb R^3\), defined by a diffeomorphism

\[F:(0, 1)^3 \to \Omega\,,\qquad \boldsymbol{\eta} \mapsto F(\boldsymbol \eta) = \mathbf x\,,\]

mapping points \(\boldsymbol{\eta} \in (0, 1)^3 = \hat\Omega\) of the (logical) unit cube to physical points \(\mathbf x \in \Omega\). The corresponding Jacobain matrix \(DF:\hat\Omega \to \mathbb R^{3\times 3}\), its volume element \(\sqrt g: \hat\Omega \to \mathbb R\) and the metric tensor \(G:\hat\Omega \to \mathbb R^{3\times 3}\) are defined by

\[DF_{i,j} = \frac{\partial F_i}{\partial \eta_j}\,,\qquad \sqrt g = |\textnormal{det}(DF)|\,,\qquad G = DF^\top DF\,.\]

Only right-handed mappings (\(\textnormal{det}(DF) > 0\)) are admitted.

property kind_map: int

Integer defining the mapping: * <=9: spline mappings * >=10 and <=19: analytical mappings with cubic domain boundary * >=20 and <=29: analytical cylinder and torus mappings * >=30 and <=39: Shafranov mappings (cylinder)

property params: dict

Mapping parameters passed to __init__() of the class in domains.py, as dictionary.

property params_numpy: ndarray

Mapping parameters as numpy array (can be empty).

property pole: bool

Bool; True if mapping has one polar point.

property periodic_eta3: bool

Bool; True if mapping is periodic in \(\eta_3\) coordinate.

property cx

3d array of control points for first mapping component \(F_x\).

property cy

3d array of control points for second mapping component \(F_y\).

property cz

3d array of control points for third mapping component \(F_z\).

property Nel

List of number of elements in each direction.

property p

List of spline degrees in each direction.

property spl_kind

List of spline type (True=periodic, False=clamped) in each direction.

property NbaseN

List of number of basis functions for N-splines in each direction.

property T

List of knot vectors for N-splines in each direction.

property indN

Global indices of non-vanishing splines in each element. Can be accessed via (element index, local spline index).

property args_domain

Object for all parameters needed for evaluation of metric coefficients.

property dict_transformations

Dictionary of str->int for pull, push and transformation functions.

__call__(*etas, change_out_order=False, squeeze_out=False, remove_outside=True, identity_map=False)

Evaluates the mapping \(F : (0, 1)^3 \to \mathbb R^3,\, \boldsymbol \eta \mapsto \mathbf x\).

Logical coordinates outside of \((0, 1)^3\) are evaluated to -1. The type of evaluation depends on the shape of the input etas.

Parameters:

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• change_out_order (bool) – If True, the axis corresponding to x, y, z coordinates in the output array is the last one, otherwise the first one.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

• identity_map (bool) – If True, not the mapping F, but the identity map (0, 1)^3 –> (0, 1)^3 is evaluated

Returns:

out – The Cartesian coordinates corresponding to the given logical ones.

Return type:

ndarray | float

jacobian(*etas, transposed=False, change_out_order=False, squeeze_out=False, remove_outside=True)

Evaluates the Jacobian matrix \(DF : (0, 1)^3 \to \mathbb R^{3 \times 3}\). Logical coordinates outside of \((0, 1)^3\) are evaluated to -1.

Parameters:

• transposed (bool) – If True, the transposed Jacobian matrix is evaluated.

• change_out_order (bool) – If True, the axes corresponding to the 3x3 entries in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

Returns:

out – The Jacobian matrix evaluated at given logical coordinates.

Return type:

ndarray | float

jacobian_det(*etas, squeeze_out=False, remove_outside=True)

Evaluates the Jacobian determinant \(\sqrt g : (0, 1)^3 \to \mathbb R^+\) (only right-handed mappings allowed). Logical coordinates outside of \((0, 1)^3\) are evaluated to -1.

Parameters:

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

Returns:

out – The Jacobian determinant evaluated at given logical coordinates.

Return type:

ndarray | float

jacobian_inv(*etas, transposed=False, change_out_order=False, squeeze_out=False, remove_outside=True, avoid_round_off=True)

Evaluates the inverse Jacobian matrix \(DF^{-1} : (0, 1)^3 \to \mathbb R^{3 \times 3}\). Logical coordinates outside of \((0, 1)^3\) are evaluated to -1.

Parameters:

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• transposed (bool) – If True, the transposed Jacobian matrix is evaluated.

• change_out_order (bool) – If True, the axes corresponding to the 3x3 entries in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

• avoid_round_off (bool) – Whether to manually set exact zeros in arrays.

Returns:

out – The inverse Jacobian matrix evaluated at given logical coordinates.

Return type:

ndarray | float

metric(*etas, transposed=False, change_out_order=False, squeeze_out=False, remove_outside=True, avoid_round_off=True)

Evaluates the metric tensor \(G: (0, 1)^3 \to \mathbb R^{3\times 3}\). Logical coordinates outside of \((0, 1)^3\) are evaluated to -1.

Parameters:

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• transposed (bool) – If True, the transposed Jacobian matrix is evaluated.

• change_out_order (bool) – If True, the axes corresponding to the 3x3 entries in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

• avoid_round_off (bool) – Whether to manually set exact zeros in arrays.

Returns:

out – The metric tensor evaluated at given logical coordinates.

Return type:

ndarray | float

metric_inv(*etas, transposed=False, change_out_order=False, squeeze_out=False, remove_outside=True, avoid_round_off=True)

Evaluates the inverse metric tensor \(G^{-1}: (0, 1)^3 \to \mathbb R^{3\times 3}\). Logical coordinates outside of \((0, 1)^3\) are evaluated to -1.

Parameters:

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• transposed (bool) – If True, the transposed Jacobian matrix is evaluated.

• change_out_order (bool) – If True, the axes corresponding to the 3x3 entries in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

• avoid_round_off (bool) – Whether to manually set exact zeros in arrays.

Returns:

out – The inverse metric tensor evaluated at given logical coordinates.

Return type:

ndarray | float

pull(a, *etas, flat_eval=False, kind='0', a_kwargs={}, change_out_order=False, squeeze_out=False, remove_outside=True, coordinates='physical')

Pull-back of a Cartesian scalar/vector field to a differential p-form.

Parameters:

• a (callable | list | tuple | array-like) – The function a(x, y, z) resp. [a_x(x, y, z), a_y(x, y, z), a_z(x, y, z)] to be pulled.

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• flat_eval (bool) – Allows to perform flat evaluation when len(etas) == 3 with 1D arrays of same size.

• kind (str) – Which pull-back to apply, ‘0’, ‘1’, ‘2’, ‘3’ or ‘v’.

• a_kwargs (dict) – Keyword arguments passed to parameter “a” if “a” is a callable: is called as a(*etas, **a_kwargs).

• change_out_order (bool) – If True, the axes corresponding to the 3 components in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

• coordinates (str) – In which coordinates the input “a” is given (in case of callables). “physical” : a = a(x, y, z). “logical” : a = a(eta1, eta2, eta3).

Returns:

out – Pullback of Cartesian vector/scalar field to p-form evaluated at given logical coordinates.

Return type:

ndarray | float

push(a, *etas, flat_eval=False, kind='0', a_kwargs={}, change_out_order=False, squeeze_out=False, remove_outside=True)

Pushforward of a differential p-form to a Cartesian scalar/vector field.

Parameters:

• a (callable | list | tuple | array-like) – The function a(e1, e2, e3) resp. [a_1(e1, e2, e3), a_2(e1, e2, e3), a_3(e1, e2, e3)] to be pushed.

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• flat_eval (bool) – Allows to perform flat evaluation when len(etas) == 3 with 1D arrays of same size.

• kind (str) – Which pushforward to apply, ‘0’, ‘1’, ‘2’, ‘3’ or ‘v’.

• a_kwargs (dict) – Keyword arguments passed to parameter “a” if “a” is a callable: is called as a(*etas, **a_kwargs).

• change_out_order (bool) – If True, the axes corresponding to the 3 components in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

Returns:

out – Pushforward of p-form to Cartesian vector/scalar field evaluated at given logical coordinates.

Return type:

ndarray | float

transform(a, *etas, flat_eval=False, kind='0_to_3', a_kwargs={}, change_out_order=False, squeeze_out=False, remove_outside=True)

Transformation between different differential p-forms and/or vector fields.

Parameters:

• a (callable | list | tuple | array-like) – The function a(e1, e2, e3) resp. [a_1(e1, e2, e3), a_2(e1, e2, e3), a_3(e1, e2, e3)] to be transformed.

• *etas (array-like | tuple) –

  Logical coordinates at which to evaluate. Two cases are possible:

  1. 2d numpy array, where coordinates are taken from eta1 = etas[:, 0], eta2 = etas[:, 1], etc. (like markers).

  2. list/tuple (eta1, eta2, …), where eta1, eta2, … can be float or array-like of various shapes.

• flat_eval (bool) – Allows to perform flat evaluation when len(etas) == 3 with 1D arrays of same size.

• kind (str) – Which transformation to apply, such as ‘0_to_3’ for example, see dict_transformations[‘tran’] for all options.

• a_kwargs (dict) – Keyword arguments passed to parameter “a” if “a” is a callable: is called as a(*etas, **a_kwargs).

• change_out_order (bool) – If True, the axes corresponding to the 3 components in the output array are the last two, otherwise the first two.

• squeeze_out (bool) – Whether to remove singleton dimensions in output array.

• remove_outside (bool) – If True, logical coordinates outside of (0, 1)^3 are NOT evaluated to -1 and are removed in the output array.

Returns:

out – Transformed p-form evaluated at given logical coordinates.

Return type:

ndarray | float

Notes

Possible choices for kind are ‘0_to_3’, ‘3_to_0’, ‘1_to_2’, ‘2_to_1’, ‘norm_to_v’, ‘norm_to_1’, ‘norm_to_2’, ‘v_to_1’, ‘v_to_2’, ‘1_to_v’ and ‘2_to_v’.

static prepare_eval_pts(x, y, z, flat_eval=False)

Broadcasts evaluation point sets to 3d arrays of correct shape.

Parameters:

• x (float | int | list | array-like) – Evaluation point sets.

• y (float | int | list | array-like) – Evaluation point sets.

• z (float | int | list | array-like) – Evaluation point sets.

• flat_eval (bool) – Whether to do a flat evaluation, i.e. f([e11, e12], [e21, e22]) = [f(e11, e21), f(e12, e22)].

Returns:

• E1, E2, E3 (array-like) – 3d arrays of correct shape for evaluation.

• is_sparse_meshgrid (bool) – Whether arguments fit sparse_meshgrid shape.

static prepare_arg(a_in, *Xs, is_sparse_meshgrid=False, a_kwargs={})

Broadcasts argument to be pulled, pushed or transformed to array of correct shape (2d for markers, 4d else).

Parameters:

• a_in (callable | list | tuple | array-like) – The argument to be pulled, pushed or transformed.

• *Xs (array-like | tuple) – The evaluation point sets. Obtained from prepare_eval_pts function.

• is_sparse_meshgrid (bool) – Whether arguments fit sparse_meshgrid shape. Obtained from prepare_eval_pts function.

• a_kwargs (dict) – Keyword arguments passed to parameter “a_in” if “a_in” is a callable: is called as a_in(*etas, **a_kwargs).

Returns:

a_out – The 2d/4d array suitable for evaluation kernels.

Return type:

array-like

get_params_numpy() → ndarray

Convert parameter dict into numpy array.

show(logical=False, grid_info=None, markers=None, marker_coords='logical', show_control_pts=False, figsize=(12, 5), save_dir=None)

Plots isolines (and control point in case on spline mappings) of the 2D physical domain for eta3 = 0. Markers can be plotted as well (optional).

Parameters:

• logical (bool) – Whether to plot the physical domain (False) or logical domain (True).

• plane (str) – Which physical coordinates to plot (xy, xz or yz) in case of logical=False.

• grid_info (array-like) –

  Information about the grid. If not given, the domain is shown with high resolution. If given, can be either

  • a list of # of elements [Nel1, Nel2, (Nel3)], OR

  • a 2d array with information about MPI decomposition.

• markers (array-like) – Markers to be plotted. Can be of shape (Np, 3) or (:, Np, 3). For the former, all markers are plotted with the same color. For the latter, with different colors (are interpreted as orbits in time).

• marker_coords (bool) – Whether the marker coordinates are logical or physical.

• save_dir (str) – If given, the figure is saved according the given directory.

__weakref__

list of weak references to the object (if defined)

class struphy.geometry.base.Spline(Nel: tuple[int] = (8, 24, 6), p: tuple[int] = (2, 3, 1), spl_kind: tuple[bool] = (False, True, True), cx: ndarray | None = None, cy: ndarray | None = None, cz: ndarray | None = None)

Bases: Domain

3D IGA spline mapping.

\[ \begin{align}\begin{aligned}F: (\eta_1, \eta_2, \eta_3) \mapsto (x, y, z) \textnormal{ as } \left\{\begin{aligned} x &= \sum_{ijk} c^x_{ijk} N_i(\eta_1) N_j(\eta_2) N_k(\eta_3)\,,\\y &= \sum_{ijk} c^y_{ijk} N_i(\eta_1) N_j(\eta_2) N_k(\eta_3)\,,\\z &= \sum_{ijk} c^z_{ijk} N_i(\eta_1) N_j(\eta_2) N_k(\eta_3)\,. \end{aligned}\right.\end{aligned}\end{align} \]

class struphy.geometry.base.PoloidalSpline(Nel: tuple[int] = (8, 24), p: tuple[int] = (2, 3), spl_kind: tuple[bool] = (False, True), cx: ndarray | None = None, cy: ndarray | None = None)

Bases: Domain

Base class for all mappings that use a 2D spline representation \(S:(\eta_1, \eta_2) \to (R, Z) \in \mathbb R^2\) in the poloidal plane:

\[ \begin{align}\begin{aligned}S: (\eta_1, \eta_2) \mapsto (R, Z) \textnormal{ as } \left\{\begin{aligned} R &= \sum_{ij} c^R_{ij} N_i(\eta_1) N_j(\eta_2) \,,\\Z &= \sum_{ij} c^Z_{ij} N_i(\eta_1) N_j(\eta_2) \,. \end{aligned}\right.\end{aligned}\end{align} \]

The full map \(F: [0, 1]^3 \to \Omega\) is obtained by defining \((R, Z, \eta_3) \mapsto (x, y, z)\) in the child class.

class struphy.geometry.base.PoloidalSplineStraight(Nel: tuple[int] = (8, 24), p: tuple[int] = (2, 3), spl_kind: tuple[bool] = (False, True), cx: ndarray | None = None, cy: ndarray | None = None, Lz: float = 4.0)

Bases: PoloidalSpline

Cylinder where the poloidal planes are described by a 2D IGA-spline mapping.

\[ \begin{align}\begin{aligned}F: (R, Z, \eta_3) \mapsto (x, y, z) \textnormal{ as } \left\{\begin{aligned} x &= R \,,\\y &= Z \,,\\z &= L_z\eta_3\,. \end{aligned}\right.\end{aligned}\end{align} \]

class struphy.geometry.base.PoloidalSplineTorus(Nel: tuple[int] = (8, 24), p: tuple[int] = (2, 3), spl_kind: tuple[bool] = (False, True), cx: ndarray | None = None, cy: ndarray | None = None, tor_period: int = 3)

Bases: PoloidalSpline

Torus where the poloidal planes are described by a 2D IGA-spline mapping.

\[ \begin{align}\begin{aligned}F: (R, Z, \eta_3) \mapsto (x, y, z) \textnormal{ as } \left\{\begin{aligned} x &= R \cos(2\pi\eta_3) \,,\\y &= R \sin(- 2\pi\eta_3) \,,\\z &= Z \,. \end{aligned}\right.\end{aligned}\end{align} \]

Parameters:

• Nel (tuple[int]) – Number of elements in each poloidal direction.

• p (tuple[int]) – Spline degree in each poloidal direction.

• spl_kind (tuple[bool]) – Kind of spline in each poloidal direction (True=periodic, False=clamped).

• cx (xp.ndarray) – Control points (spline coefficients) of the poloidal mapping. If None, a default square-to-disc mapping of radius 1 centered around (x, y) = (3, 0) is interpolated.

• cy (xp.ndarray) – Control points (spline coefficients) of the poloidal mapping. If None, a default square-to-disc mapping of radius 1 centered around (x, y) = (3, 0) is interpolated.

• tor_period (int) – The toroidal angle is between [0, 2*pi/tor_period).

struphy.geometry.base.interp_mapping(Nel, p, spl_kind, X, Y, Z=None)

Interpolates the mapping \(F: (0, 1)^3 \to \mathbb R^3\) on the given spline space.

Parameters:

• Nel (array-like) – Defining the spline space.

• p (array-like) – Defining the spline space.

• spl_kind (array-like) – Defining the spline space.

• X (callable) – Either X(eta1, eta2) in 2D or X(eta1, eta2, eta3) in 3D.

• Y (callable) – Either X(eta1, eta2) in 2D or X(eta1, eta2, eta3) in 3D.

• Z (callable) – 3rd mapping component Z(eta1, eta2, eta3) in 3D.

Returns:

cx, cy (, cz) – The control points.

Return type:

array-like

struphy.geometry.base.spline_interpolation_nd(p: list, spl_kind: list, grids_1d: list, values: ndarray)

n-dimensional tensor-product spline interpolation with discrete input.

The interpolation points are passed as a list of 1d arrays, each array with increasing entries g[0]=0 < g[1] < … The last element must be g[-1] = 1 for clamped interpolation and g[-1] < 1 for periodic interpolation.

Parameters:

• p (list[int]) – Spline degree.

• grids_1d (list[array]) – Interpolation points in [0, 1].

• spl_kind (list[bool]) – True: periodic splines, False: clamped splines.

• values (array) – Function values at interpolation points. values.shape = (grid1.size, …, gridn.size).

Returns:

• coeffs (xp.array) – spline coefficients as nd array.

• T (list[array]) – Knot vector of spline interpolant.

• indN (list[array]) – Global indices of non-vanishing splines in each element. Can be accessed via (element, local index).