Device Interface
The tdgl.device
subpackage provides the following functionalities:
Definition of material properties in instances of
tdgl.Layer
.Definition of device geometry in terms of
tdgl.Polygon
instances, which can be created from simple geometric primitives using constructive solid geometry.Mesh generation for
tdgl.Polygon
andtdgl.Device
instances.Translation between physical units (e.g., microns and microamperes) and the dimensionless units used in TDGL.
Visualization and serialization of
tdgl.Device
instances.
Overview
Here is a quick overview of the most useful functions and methods for creating and working with polygons and devices. For a demonstration of how all of these pieces work together in practice, see Working with polygons.
- Geometric primitives:
tdgl.Polygon
constructive solid geometry methods:
tdgl.Polygon
geometrical transformation methods:
tdgl.Device
geometrical transformation methods:
- Visualization methods:
- Mesh generation methods:
Layer
- class tdgl.Layer(*, london_lambda, coherence_length, thickness, conductivity=None, u=5.79, gamma=10.0, z0=0)[source]
A superconducting thin film.
- Parameters:
london_lambda (
float
) – The London penetration depth of the film.coherence_length (
float
) – The superconducting coherence length of the film.thickness (
float
) – The thickness of the film.conductivity (
Optional
[float
]) – The normal state conductivity of the superconductor in Siemens / length_unit.u (
float
) – The ratio of the relaxation times for the order parameter amplitude and phase. This value is 5.79 for dirty superconductors.gamma (
float
) – This parameter quantifies the effect of inelastic phonon-electron scattering. \(\gamma\) is proportional to the inelastic scattering time and the size of the superconducting gap.z0 (
float
) – Vertical location of the film.
- copy()[source]
Create a deep copy of the
tdgl.Layer
.- Return type:
- to_hdf5(h5_group)[source]
Save the
tdgl.Layer
to anh5py.Group
.- Parameters:
h5_group (
Group
) – An openh5py.Group
to which to save the layer.- Return type:
- static from_hdf5(h5_group)[source]
Load a
tdgl.Layer
from anh5py.Group
.- Parameters:
h5_group (
Group
) – An openh5py.Group
from which to load the layer.- Return type:
- Returns:
A new
tdgl.Layer
instance.
Polygon
- class tdgl.Polygon(name=None, *, points, mesh=True)[source]
A polygonal region located in a
tdgl.Layer
.- Parameters:
points (
Union
[Polygon
,ndarray
,LineString
,LinearRing
,Polygon
]) – The polygon vertices. This can be a shape(n, 2)
array of x, y coordinates or a shapelyLineString
,LinearRing
, orPolygon
.mesh (
bool
) – Whether to include this polygon when computing a mesh.
- property bbox: Tuple[Tuple[float, float], Tuple[float, float]]
Returns the coordinates of the lower left and upper right corners of the polygon’s bounding box.
- property extents: Tuple[float, float]
Returns the total x, y extent of the polygon, (Delta_x, Delta_y).
- contains_points(points, index=False, radius=0)[source]
Determines whether
points
lie within the polygon.- Parameters:
points (
ndarray
) – Shape(n, 2)
array of x, y coordinates.index (
bool
) – If True, then return the indices of the points inpoints
that lie within the polygon. Otherwise, returns a shape(n, )
boolean array.radius (
float
) – An additional margin onself.path
. Seematplotlib.path.Path.contains_points()
.
- Return type:
- Returns:
If index is True, returns the indices of the points in
points
that lie within the polygon. Otherwise, returns a shape(n, )
boolean array indicating whether each point lies within the polygon.
- on_boundary(points, radius=0.001, index=False)[source]
Determines whether
points
lie within a given radius of the Polygon boundary.- Parameters:
- Returns:
If index is True, returns the indices of the points in
points
that lie within the polygon. Otherwise, returns a shape(n, )
boolean array indicating whether each point lies within the polygon.
- make_mesh(min_points=None, smooth=0, **meshpy_kwargs)[source]
Returns the vertices and triangles of a Delaunay mesh covering the Polygon.
- Parameters:
min_points (
Optional
[int
]) – Minimum number of vertices in the mesh. If None, then the number of vertices will be determined by meshpy_kwargs and the number of vertices in the underlying polygons.smooth (
int
) – Number of Laplacian smoothing steps to perform.**meshpy_kwargs – Passed to meshpy.triangle.build().
- Return type:
- Returns:
Mesh vertex coordinates and triangle indices
- rotate(degrees, origin=(0.0, 0.0), inplace=False)[source]
Rotates the polygon counterclockwise by a given angle.
- Parameters:
degrees (
float
) – The amount by which to rotate the polygon.origin (
Union
[str
,Tuple
[float
,float
]]) – (x, y) coorindates about which to rotate, or the strings “center” (for the bounding box center) or “centroid” (for the polygon center of mass).inplace (
bool
) – If True, modify the polygon in place. Otherwise, return a modified copy.
- Return type:
- Returns:
The rotated polygon.
- scale(xfact=1.0, yfact=1.0, origin=(0, 0), inplace=False)[source]
Scales the polygon horizontally by
xfact
and vertically byyfact
.Negative
xfact
(yfact
) can be used to reflect the polygon horizontally (vertically) about theorigin
.- Parameters:
xfact (
float
) – Distance by which to translate along the x-axis.yfact (
float
) – Distance by which to translate along the y-axis.origin (
Union
[str
,Tuple
[float
,float
]]) – (x, y) coorindates for the scaling origin, or the strings “center” (for the bounding box center) or “centroid” (for the polygon center of mass).inplace (
bool
) – If True, modify the polygon in place. Otherwise, return a modified copy.
- Return type:
- Returns:
The scaled polygon.
- union(*others, name=None)[source]
Returns the union of the polygon and zero or more other polygons.
- Parameters:
others (
Union
[Polygon
,ndarray
,LineString
,LinearRing
,Polygon
]) – One or more objects with which to join the polygon.name (
Optional
[str
]) – A name for the resulting joined Polygon (defaults toself.name
.)
- Return type:
- Returns:
A new
Polygon
instance representing the union ofself
andothers
.
- intersection(*others, name=None)[source]
Returns the intersection of the polygon and zero or more other polygons.
- Parameters:
others (
Union
[Polygon
,ndarray
,LineString
,LinearRing
,Polygon
]) – One or more objects with which to join the polygon.name (
Optional
[str
]) – A name for the resulting joined Polygon (defaults toself.name
.)
- Return type:
- Returns:
A new
Polygon
instance representing the intersection ofself
andothers
.
- difference(*others, name=None)[source]
Returns the difference of the polygon and zero more other polygons.
- Parameters:
others (
Union
[Polygon
,ndarray
,LineString
,LinearRing
,Polygon
]) – One or more objects with which to join the polygon.name (
Optional
[str
]) – A name for the resulting joined Polygon (defaults toself.name
.)
- Return type:
- Returns:
A new
Polygon
instance representing the difference ofself
andothers
.
- buffer(distance, join_style='mitre', mitre_limit=5.0, single_sided=True, as_polygon=True)[source]
Returns polygon points or a new Polygon object with vertices offset from
self.points
by a givendistance
. Ifdistance > 0
this “inflates” the polygon, and ifdistance < 0
this shrinks the polygon.- Parameters:
distance (
float
) – The amount by which to inflate or deflate the polygon.join_style (
Union
[str
,int
]) – One of “round” (1), “mitre” (2), or “bevel” (3). See the shapely documentation.mitre_limit (
float
) – See the shapely documentation.single_sided (
bool
) – See the shapely documentation.as_polygon (
bool
) – If True, returns a newPolygon
instance, otherwise returns a shape(n, 2)
array of polygon vertices.
- Return type:
- Returns:
A new
Polygon
or an array of vertices offset bydistance
.
- resample(num_points=None, degree=1, smooth=0)[source]
Resample vertices so that they are approximately uniformly distributed along the polygon boundary.
- Parameters:
num_points (
Optional
[int
]) – Number of points to interpolate to. Ifnum_points
is None, the polygon is resampled tolen(self.points)
points. Ifnum_points
is not None and has a boolean value of False, then an unaltered copy of the polygon is returned.degree (
int
) – The degree of the spline with which to interpolate. Defaults to 1 (linear spline).smooth (
float
) – Smoothing condition.
- Return type:
- classmethod from_union(items, *, name=None, mesh=True)[source]
Creates a new
Polygon
from the union of a sequence of polygons.
- classmethod from_intersection(items, *, name=None, mesh=True)[source]
Creates a new
Polygon
from the intersection of a sequence of polygons.
- classmethod from_difference(items, *, name=None, mesh=True)[source]
Creates a new
Polygon
from the difference of a sequence of polygons.
- to_hdf5(h5_group)[source]
Save the
Polygon
to anh5py.Group
.- Return type:
- Parameters:
h5_group (Group) –
Device
- class tdgl.Device(name, *, layer, film, holes=None, terminals=None, probe_points=None, length_units='um')[source]
An object representing a thin film superconducting device.
- Parameters:
name (
str
) – Name of the device.layer (
Layer
) – The superconductingLayer
.film (
Polygon
) – ThePolygon
representing the superconducting film.holes (
Optional
[List
[Polygon
]]) –Polygons
representing holes in the superconducting film.terminals (
Optional
[List
[Polygon
]]) – A sequence ofPolygons
representing the current terminals. Any points that are on the boundary of the mesh and lie inside a terminal will have current source/sink boundary conditions.probe_points (
Optional
[Sequence
[Tuple
[float
,float
]]]) – A shape(n, 2)
sequence of floats, with each row representing the(x, y)
position of a voltage probe.length_units (
str
) – Distance units for the coordinate system.
- property coherence_length: Quantity
Ginzburg-Landau coherence length, \(\xi\)
- property london_lambda: Quantity
London penetration depth, \(\lambda\)
- property thickness: Quantity
Film thickness, \(d\)
- property Lambda: Quantity
Effective magnetic penetration depth, \(\Lambda=\lambda^2/d\).
- property Bc2: Quantity
Upper critical field, \(B_{c2}=\Phi_0/(2\pi\xi^2)\).
- property A0: Quantity
Scale for the magnetic vector potential, \(A_0=\xi B_{c2}\).
- property K0: Quantity
Sheet current density scale (dimensions of current / length), \(K_0=4\xi B_{c2}/(\mu_0\Lambda)\).
- V0(conductivity=None)[source]
Electric potential scale, \(\V_0=\xi J_0/\sigma\).
- Parameters:
conductivity (
Optional
[Quantity
]) – The normal state conductivity of the film, which defaults todevice.layer.conductivity
.- Return type:
Quantity
- Returns:
The electric potential scale, \(\V_0=\xi J_0/\sigma\)
- property triangulation: Optional[Triangulation]
Matplotlib triangulation of the mesh.
- terminal_info()[source]
Returns a tuple of
TerminalInfo
objects, one for each current terminal in the device.
- property points: Optional[ndarray]
The mesh vertex coordinates in
length_units
(shape(n, 2)
, typefloat
).
- property probe_point_indices: Optional[List[int]]
A list of the mesh site indices for the probe points.
- boundary_sites()[source]
Returns a dict of
{polygon_name: boundary_indices}
, whereboundary_indices
is an integer array of site indices for mesh sites on the boundary of each polygon.The length of the returned dictionary will be the number of holes in the device plus one.
- contains_points(points, index=False, radius=0)[source]
Determines whether
points
lie within the Device.- Parameters:
points (
ndarray
) – Shape(n, 2)
array of x, y coordinates.index (
bool
) – If True, then return the indices of the points inpoints
that lie within the polygon. Otherwise, returns a shape(n, )
boolean array.radius (
float
) – An additional margin onself.path
. Seematplotlib.path.Path.contains_points()
.
- Return type:
- Returns:
If index is True, returns the indices of the points in
points
that lie within the polygon. Otherwise, returns a shape(n, )
boolean array indicating whether each point lies within the polygon.
- scale(xfact=1, yfact=1, origin=(0, 0))[source]
Returns a new device with polygons scaled horizontally and/or vertically.
Negative
xfact
(yfact
) can be used to reflect the device horizontally (vertically) about theorigin
.
- rotate(degrees, origin=(0, 0))[source]
Returns a new device with polygons rotated a given amount counterclockwise about specified origin.
- translate(dx=0, dy=0, dz=0, inplace=False)[source]
Translates the device polygons, layer, and mesh in space by a given amount.
- Parameters:
dx (
float
) – Distance by which to translate along the x-axis.dy (
float
) – Distance by which to translate along the y-axis.dz (
float
) – Distance by which to translate layer along the z-axis.inplace (
bool
) – If True, modifies the device (self
) in-place and returns None, otherwise, creates a new device, translates it, and returns it.
- Return type:
- Returns:
The translated device.
- translation(dx, dy, dz=0)[source]
A context manager that temporarily translates a device in-place, then returns it to its original position.
- make_mesh(max_edge_length=None, min_points=None, smooth=0, **meshpy_kwargs)[source]
Generates and optimizes the triangular mesh.
- Parameters:
min_points (
Optional
[float
]) – Minimum number of vertices in the mesh. If None, then the number of vertices will be determined by meshpy_kwargs, the number of vertices in the underlying polygons, andmax_edge_length
.max_edge_length (
Optional
[float
]) – The maximum distance between vertices in the mesh. Passing a value <= 0 means that the number of mesh points will be determined solely by the density of points in the Device’s film and holes. Defaults to 1.0 * self.coherence_length.smooth (
int
) – Number of Laplacian smoothing iterations to perform.**meshpy_kwargs – Passed to meshpy.triangle.build().
- Return type:
- mesh_stats(precision=3)[source]
When called with in Jupyter notebook, displays a table of information about the mesh.
- Parameters:
precision (
int
) – Number of digits after the decimal for float values.- Return type:
HTML
- Returns:
An HTML table of mesh statistics.
- plot(ax=None, legend=True, figsize=None, mesh=False, mesh_kwargs={'color': 'k', 'lw': 0.5}, **kwargs)[source]
Plot all of the device’s polygons.
- Parameters:
ax (
Optional
[Axes
]) – matplotlib axis on which to plot. If None, a new figure is created.subplots – If True, plots each layer on a different subplot.
legend (
bool
) – Whether to add a legend.figsize (
Optional
[Tuple
[float
,float
]]) – matplotlib figsize, only used if ax is None.mesh (
bool
) – If True, plot the mesh.mesh_kwargs (
Dict
[str
,Any
]) – Keyword arguments passed toax.triplot()
ifmesh
is True.kwargs – Passed to
ax.plot()
for the polygon boundaries.
- Return type:
- Returns:
Matplotlib Figure and Axes
- draw(ax=None, legend=True, figsize=None, alpha=0.5, exclude=None)[source]
Draws all polygons in the device as matplotlib patches.
- Parameters:
ax (
Optional
[Axes
]) – matplotlib axis on which to plot. If None, a new figure is created.legend (
bool
) – Whether to add a legend.figsize (
Optional
[Tuple
[float
,float
]]) – matplotlib figsize, only used if ax is None.alpha (
float
) – The alpha (opacity) value for the patches (0 <= alpha <= 1).exclude (
Union
[str
,List
[str
],None
]) – A polygon name or list of polygon names to exclude from the figure.
- Return type:
- Returns:
Matplotlib Figre and Axes.
Geometry
The tdgl.geometry
module contains functions for creating polygons approximating simple
shapes, i.e., rectangles and ellipses, which can be combined using
constructive solid geometry
to create complex shapes.
- tdgl.geometry.box(width, height=None, points=101, center=(0, 0), angle=0)[source]
Returns the coordinates for a rectangle with a given width and height, centered at the specified center.
- Parameters:
width (
float
) – Width of the rectangle (in the x direction).height (
Optional
[float
]) – Height of the rectangle (in the y direction). If None is given, then height is set to width and the function returns a square.points (
int
) – The target number of points making up the box. The actual number of points may be slightly different than this value.center (
Tuple
[float
,float
]) – Coordinates of the center of the rectangle.angle (
float
) – Angle (in degrees) by which to rotate counterclockwise about (0, 0) after translating to the specified center.
- Return type:
- Returns:
A shape
(m, 2)
or array of (x, y) coordinates.
- tdgl.geometry.circle(radius, points=100, center=(0, 0))[source]
Returns the coordinates for a circle with a given radius, centered at the specified center.
- tdgl.geometry.ellipse(a, b, points=100, center=(0, 0), angle=0)[source]
Returns the coordinates for an ellipse with major axis a and semimajor axis b, rotated by the specified angle about (0, 0), then translated to the specified center.
- Parameters:
a (
float
) – Major axis lengthb (
float
) – Semi-major axis lengthpoints (
int
) – Number of points in the circlecenter (
Tuple
[float
,float
]) – Coordinates of the center of the circleangle (
float
) – Angle (in degrees) by which to rotate counterclockwise about (0, 0) after translating to the specified center.
- Returns:
A shape
(points, 2)
array of (x, y) coordinates
Mesh Generation
- tdgl.generate_mesh(poly_coords, hole_coords=None, min_points=None, max_edge_length=None, convex_hull=False, boundary=None, min_angle=32.5, **kwargs)[source]
Generates a Delaunay mesh for a given set of polygon vertex coordinates.
Additional keyword arguments are passed to
triangle.build()
.- Parameters:
poly_coords (
ndarray
) – Shape(n, 2)
array of polygon(x, y)
coordinates.hole_coords (
Optional
[List
[ndarray
]]) – A list of arrays of hole boundary coordinates.min_points (
Optional
[int
]) – The minimimum number of vertices in the resulting mesh.max_edge_length (
Optional
[float
]) – The maximum distance between vertices in the resulting mesh.convex_hull (
bool
) – If True, then the entire convex hull of the polygon (minus holes) will be meshed. Otherwise, only the polygon interior is meshed.boundary (
Optional
[ndarray
]) – Shape(m, 2)
(wherem <= n
) array of(x, y)
coordinates for points on the boundary of the polygon.min_angle (
float
) – The minimum angle in the mesh’s triangles. Setting a larger value will make the triangles closer to equilateral, but the mesh generation may fail if the value is too large.
- Return type:
- Returns:
Mesh vertex coordinates and triangle indices.