Post-Processing
The tdgl.Solution
class provides a convenient container for the results of a TDGL simulation,
including methods for post-processing and visualization results. Calls to tdgl.solve()
return an
instance of tdgl.Solution
, which can be used for post-processing. A tdgl.Solution
can be serialized to and deserialized from disk.
For each instance solution
of tdgl.Solution
, the raw data from the TDGL simulation (in dimensionless units)
are stored in solution.tdgl_data
, which is an instance of tdgl.solution.data.TDGLData
. Any data
that is measured at each time step in the simulation, i.e., the measured voltage and phase difference between
the tdgl.Device
’s probe_points
, are stored in solution.dynamics
, which is an instance of
tdgl.solution.data.DynamicsData
.
Overview
Post-processing methods:
Visualization methods:
I/O methods:
Solution
- class tdgl.Solution(*, device, options, path, applied_vector_potential, terminal_currents, disorder_epsilon, total_seconds, _solve_step=-1)[source]
A container for the results of a TDGL simulation.
- Parameters:
device (
Device
) – Thetdgl.Device
that was solvedoptions (
SolverOptions
) – Atdgl.SolverOptions
instance.path (
str
) – Path to the HDF5 file containing the raw output data.applied_vector_potential (
Parameter
) – TheParameter
defining the applied vector potential.terminal_currents (
Union
[Dict
[str
,float
],Callable
]) – A dict of{terminal_name: current}
or a callable with signaturefunc(time) -> {terminal_name: current}
, wherecurrent
is a float in units ofcurrent_units
.disorder_epsilon (
Union
[float
,Callable
]) – The disorder parameter \(\epsilon\). If \(\epsilon(\mathbf{r}) < 1\) weakens the order parameter at position \(\mathbf{r}\), which can be used to model inhomogeneity.total_seconds (
float
) – The total wall time in seconds._solve_step (int) –
-
dynamics:
Optional
[DynamicsData
] A container for the time dynamics of the solution (in dimensionless units).
- property solve_step: int
The solver iteration corresponding to the current
tdgl.solution.data.TDGLData
.Setting
solve_step
automatically loads data for the specitied step.
- load_tdgl_data(solve_step=-1, h5file=None)[source]
Loads the TDGL results from file for a given solve step.
- property vorticity: Optional[ndarray]
The current vorticity, \(\omega=(\nabla\times\mathbf{K})\cdot\hat{\mathbf{z}}\)
- property current_density: ndarray
The total sheet current density, \(\mathbf{K}=\mathbf{K}_s+\mathbf{K}_n\).
- magnetic_moment(units=None, with_units=True)[source]
Computes the \(z\)-component of the magnetic dipole moment of the film.
\[m_z = \hat{z}\cdot\frac{1}{2} \int_\mathrm{film}\mathbf{r}\times\mathbf{K}(\mathbf{r})\,\mathrm{d}^2r,\]where \(\mathbf{r}\) is the position in the film relative to the film center of mass.
- Parameters:
- Return type:
- Returns:
The magnetic dipole moment of the film as either a float or a pint.Quantity
- grid_current_density(*, dataset=None, grid_shape=(200, 200), method='linear', units=None, with_units=False, **kwargs)[source]
Interpolates the sheet current density to a rectangular grid.
Keyword arguments are passed to
scipy.interpolate.griddata()
.- Parameters:
dataset (
Optional
[str
]) – The dataset to interpolate. One of"supercurrent"
,"normal_current"
, orNone
. IfNone
, then the total sheet current density is used.grid_shape (
Union
[int
,Tuple
[int
,int
]]) – Shape of the desired rectangular grid. If a single integer N is given, then the grid will be square, shape = (N, N).method (
str
) – Interpolation method to use (seescipy.interpolate.griddata()
).units (
Optional
[str
]) – The desired units for the current density. Defaults toself.current_units / self.device.length_units
.with_units (
bool
) – Whether to return apint.Quantity
array with units attached.
- Return type:
- Returns:
x grid, y grid, interpolated current density
- interp_current_density(positions, *, dataset=None, method='linear', units=None, with_units=False)[source]
Interpolates the sheet current density at unstructured coordinates.
See also
- Parameters:
positions (
ndarray
) – Shape(m, 2)
array of x, y coordinates at which to evaluate the current density.dataset (
Optional
[str
]) – The dataset to interpolate. One of"supercurrent"
,"normal_current"
, orNone
. IfNone
, then the total sheet current density is used.method (
Literal
['linear'
,'cubic'
]) – Interpolation method to use,"linear"
or"cubic"
.units (
Optional
[str
]) – The desired units for the current density. Defaults toself.current_units / self.device.length_units
.with_units (
bool
) – Whether to return apint.Quantity
array with units attached.
- Return type:
- Returns:
The interpolated current density as an array of floats or a
pint.Quantity
array.
- interp_order_parameter(positions, method='linear')[source]
Interpolates the order parameter at unstructured coordinates.
- polygon_fluxoid(polygon_points, interp_method='linear', units='Phi_0', with_units=True)[source]
Computes the
tdgl.Fluxoid
(flux + supercurrent) for a given polygonal region.The fluxoid for a closed region \(S\) with boundary \(\partial S\) is defined as:
\[\begin{split}\begin{split} \Phi^f_S &= \Phi^f_{S,\text{ flux}} + \Phi^f_{S,\text{ supercurrent}} \\&=\int_S \mu_0 H_z(\mathbf{r})\,\mathrm{d}^2r + \oint_{\partial S} \mu_0\Lambda(\mathbf{r})\mathbf{K}_s(\mathbf{r})\cdot\mathrm{d}\mathbf{r} \\&=\oint_{\partial S} \mathbf{A}(\mathbf{r})\cdot\mathrm{d}\mathbf{r} + \oint_{\partial S} \mu_0\Lambda(\mathbf{r})\mathbf{K}_s(\mathbf{r})\cdot\mathrm{d}\mathbf{r} \end{split}\end{split}\]See also
- Parameters:
polygon_points (
Union
[ndarray
,Polygon
]) – A shape(n, 2)
array of(x, y)
coordinates of polygon vertices defining the closed region \(S\).interp_method (
Literal
['linear'
,'cubic'
]) – Interpolation method to use,"linear"
or"cubic"
.units (
str
) – The desired units for the fluxoid.with_units (
bool
) – Whether to return values aspint.Quantity
instances with units attached.
- Return type:
- Returns:
The polygon’s
Fluxoid
.
- hole_fluxoid(hole_name, points=None, interp_method='linear', units='Phi_0', with_units=True)[source]
Calculcates the fluxoid for a polygon enclosing the specified hole.
- Parameters:
hole_name (
str
) – The name of the hole for which to calculate the fluxoid.points (
Optional
[ndarray
]) – The vertices of the polygon enclosing the hole. If None is given, a polygon is generated usingtdgl.make_fluxoid_polygons()
.interp_method (
Literal
['linear'
,'cubic'
]) – Interpolation method to use,"linear"
or"cubic"
.units (
str
) – The desired units for the fluxoid.with_units (
bool
) – Whether to return values aspint.Quantity
instances with units attached.
- Return type:
- Returns:
The hole’s
tdgl.Fluxoid
.
- boundary_phases(delta=False)[source]
Returns a dict of
{polygon_name: (boundary_indices, boundary_phases)}
.(boundary_phases[-1] - boundary_phases[0]) / (2 * np.pi)
gives the winding number for the polygon, i.e., the fluxoid in units ofPhi_0
.See also
- current_through_path(path_coords, dataset=None, method='linear', units=None, with_units=True)[source]
Calculates the total current crossing a given path.
- Parameters:
path_coords (
ndarray
) – An(n, 2)
array of(x, y)
coordinates defining the path.dataset (
Optional
[str
]) –None
,"supercurrent"
, or"normal_current"
.None
indicates the total current.method (
Literal
['linear'
,'cubic'
]) – Interpolation method: either “linear” or “cubic”.with_units (
bool
) – Whether to return apint.Quantity
with units attached.
- Return type:
- Returns:
The total current crossing the path as either a float or a
pint.Quantity
.
- field_at_position(positions, *, zs=None, vector=False, units=None, with_units=True, return_sum=True)[source]
Calculates the field due to currents in the device at any point(s) in space.
See also
- Parameters:
positions (
ndarray
) – Shape (m, 2) array of (x, y) coordinates, or (m, 3) array of (x, y, z) coordinates at which to calculate the magnetic field. A single sequence like [x, y] or [x, y, z] is also allowed.zs (
Union
[float
,ndarray
,None
]) – z coordinates at which to calculate the field. If positions has shape (m, 3), then this argument is not allowed. If zs is a scalar, then the fields are calculated in a plane parallel to the x-y plane. If zs is any array, then it must be same length as positions.vector (
bool
) – Whether to return the full vector magnetic field or just the z component.units (
Optional
[str
]) – Units to which to convert the fields (can be either magnetic field H or magnetic flux density B = mu0 * H). If not given, then the fields are returned in units ofself.field_units
.with_units (
bool
) – Whether to return the fields aspint.Quantity
with units attached.return_sum (
bool
) – IfFalse
, this method will return atdgl.BiotSavartField
instance, where the field from the supercurrent and normal current are identified separately.
- Return type:
- Returns:
An np.ndarray if
return_sum
isTrue
, otherwise an instance oftdgl.BiotSavartField
. Ifwith_units
isTrue
, then the array(s) will be of typepint.Quantity
. The array(s) will have shape(m, )
if vector is False, or shape(m, 3)
ifvector
is True.
- vector_potential_at_position(positions, *, zs=None, units=None, with_units=True, return_sum=True)[source]
Calculates the vector potential due to currents in the device at any point(s) in space, plus the applied vector potential.
The vector potential \(\mathbf{A}\) at position \(\mathbf{r}\) due to sheet current density \(\mathbf{K}(\mathbf{r}')\) flowing in a film with lateral geometry \(S\) is:
\[\mathbf{A}(\mathbf{r}) = \frac{\mu_0}{4\pi} \int_S\frac{\mathbf{K}(\mathbf{r}')}{|\mathbf{r}-\mathbf{r}'|}\mathrm{d}^2r'.\]- Parameters:
positions (
ndarray
) – Shape (m, 2) array of (x, y) coordinates, or (m, 3) array of (x, y, z) coordinates at which to calculate the vector potential. A single list like [x, y] or [x, y, z] is also allowed.zs (
Union
[float
,ndarray
,None
]) – z coordinates at which to calculate the potential. If positions has shape (m, 3), then this argument is not allowed. If zs is a scalar, then the fields are calculated in a plane parallel to the x-y plane. If zs is any array, then it must be same length as positions.units (
Optional
[str
]) – Units to which to convert the vector potential.with_units (
bool
) – Whether to return the vector potential as apint.Quantity
with units attached.return_sum (
bool
) – Whether to return the total potential or a dict with keys("applied", "supercurrent", "normal_current")
.
- Return type:
- Returns:
An np.ndarray if
return_sum
isTrue
, otherwise a dict of{source: potential_from_source}
. Ifwith_units
isTrue
, then the array(s) will be of typepint.Quantity
.potential_from_source
will have shape(m, 3)
.
- to_hdf5(h5path=None, save_mesh=True)[source]
Save the Solution to the existing output HDF5 file or to a new HDF5 file.
- Parameters:
h5path (
Optional
[str
]) – Path to an HDF5 file. IfNone
is given, thetdgl.Solution
will be saved to the existing HDF5 output file located atself.path
.save_mesh (
bool
) – Whether to save the Device’s mesh.
- Return type:
- static from_hdf5(path, solve_step=-1)[source]
Loads a
tdgl.Solution
from file.- Parameters:
path (
str
) – Path to the HDF5 file containing a serializedtdgl.Solution
.solve_step (
int
) – The solve step to load.
- Return type:
- Returns:
The loaded Solution instance.
- delete_hdf5()[source]
Delete the HDF5 file accompanying the
tdgl.Solution
.- Return type:
- equals(other, require_same_timestamp=False)[source]
Checks whether two solutions are equal.
- Parameters:
other (
Any
) – Thetdgl.Solution
to compare for equality.require_same_timestamp (
bool
) – If True, two solutions are only considered equal if they have the exact same time_created.
- Return type:
- Returns:
A boolean indicating whether the two solutions are equal
- plot_currents(**kwargs)[source]
An alias for
tdgl.plot_currents()
.Plots the sheet current density for a given
tdgl.Solution
.Additional keyword arguments are passed to
plt.subplots()
.- Parameters:
dataset – The dataset to plot, either
"supercurrent"
or"normal_current"
.None
indicates the total current density.ax – Matplotlib axes on which to plot.
units – Units in which to plot the current density. Defaults to
solution.current_units / solution.device.length_units
.cmap – Name of the matplotlib colormap to use.
colorbar – Whether to add a colorbar to each subplot.
auto_range_cutoff – Cutoff percentile for
tdgl.solution.plot_solution.auto_range_iqr()
.symmetric_color_scale – Whether to use a symmetric color scale (vmin = -vmax).
vmin – Color scale minimum to use for all layers
vmax – Color scale maximum to use for all layers
streamplot – Whether to overlay current streamlines on the plot.
min_stream_amp – Streamlines will not be drawn anywhere the current density is less than min_stream_amp * max(current_density). This avoids streamlines being drawn where there is no current flowing.
cross_section_coords – Shape (m, 2) array of (x, y) coordinates for a cross-section (or a list of such arrays).
- Return type:
- Returns:
matplotlib figure and axes
- plot_order_parameter(**kwargs)[source]
An alias for
tdgl.plot_order_parameter()
.Plots the magnitude (or the magnitude squared) and phase of the complex order parameter, \(\psi=|\psi|e^{i\theta}\).
- Parameters:
squared – Whether to plot the magnitude squared, \(|\psi|^2\).
mag_cmap – Name of the colormap to use for the magnitude.
phase_cmap – Name of the colormap to use for the phase.
shading – May be
"flat"
or"gouraud"
. The latter does some interpolation.
- Return type:
- Returns:
matplotlib Figure and an array of two Axes objects.
- plot_field_at_positions(positions, **kwargs)[source]
An alias for
tdgl.plot_field_at_positions()
.Plots the Biot-Savart field (either all three components or just the z component) at a given set of positions (x, y, z) outside of the device.
Note
This function plots only the field due to currents flowing in the device. It does not include the applied field.
Additional keyword arguments are passed to
plt.subplots()
. This function first evaluates the field atpositions
, then interpolates the resulting fields to a rectangular grid for plotting.- Parameters:
positions (
ndarray
) – Shape (m, 2) array of (x, y) coordinates, or (m, 3) array of (x, y, z) coordinates at which to calculate the magnetic field.zs – z coordinates at which to calculate the field. If positions has shape (m, 3), then this argument is not allowed. If zs is a scalar, then the fields are calculated in a plane parallel to the x-y plane. If zs is an array, then it must be same length as positions.
vector – Whether to plot the full vector magnetic field or just the z component.
units – Units in which to plot the fields. Defaults to
solution.field_units
.grid_shape – Shape of the desired rectangular grid. If a single integer
n
is given, then the grid will be square, shape(n, n)
.grid_method – Interpolation method to use (see
scipy.interpolate.griddata()
).max_cols – Maximum number of columns in the grid of subplots.
cmap – Name of the matplotlib colormap to use.
colorbar – Whether to add a colorbar to each subplot.
auto_range_cutoff – Cutoff percentile for
tdgl.solution.plot_solution.auto_range_iqr()
.share_color_scale – Whether to force all layers to use the same color scale.
symmetric_color_scale – Whether to use a symmetric color scale (vmin = -vmax).
vmin – Color scale minimum to use for all layers
vmax – Color scale maximum to use for all layers
cross_section_coords – Shape (m, 2) array of (x, y) coordinates for a cross-section (or a list of such arrays).
- Return type:
- Returns:
matplotlib figure and axes
- plot_vorticity(**kwargs)[source]
An alias for
tdgl.plot_vorticity()
.Plots the vorticity in the film: \(\mathbf{\omega}=\mathbf{\nabla}\times\mathbf{K}\).
- Parameters:
ax – Matplotlib axes on which to plot.
cmap – Name of the matplotlib colormap to use.
units – The units in which to plot the vorticity. Must have dimensions of [current] / [length]^2.
auto_range_cutoff – Cutoff percentile for
tdgl.solution.plot_solution.auto_range_iqr()
.symmetric_color_scale – Whether to use a symmetric color scale (vmin = -vmax).
vmin – Color scale minimum.
vmax – Color scale maximum.
shading – May be
"flat"
or"gouraud"
. The latter does some interpolation.
- Return type:
- Returns:
matplotlib Figure and and Axes.
- plot_scalar_potential(**kwargs)[source]
An alias for
tdgl.plot_scalar_potential()
.Plots the scalar potential \(\mu(\mathbf{r})\) in the film.
- Parameters:
ax – Matplotlib axes on which to plot.
cmap – Name of the matplotlib colormap to use.
auto_range_cutoff – Cutoff percentile for
tdgl.solution.plot_solution.auto_range_iqr()
.vmin – Color scale minimum.
vmax – Color scale maximum.
shading – May be
"flat"
or"gouraud"
. The latter does some interpolation.
- Return type:
- Returns:
matplotlib Figure and and Axes.
- class tdgl.solution.data.TDGLData(step, epsilon, psi, mu, applied_vector_potential, induced_vector_potential, supercurrent, normal_current, state)[source]
A container for raw data from the TDGL solver at a single solve step.
- Parameters:
step (
int
) – The solver iteration.epsilon (
ndarray
) – The disorder parameter. \(\epsilon<1\) weakens the order parameter.psi (
ndarray
) – The complex order parameter at each site in the mesh.mu (
ndarray
) – The scalar potential at each site in the mesh.applied_vector_potential (
ndarray
) – The applied vector potential at each edge in the mesh.induced_vector_potential (
ndarray
) – The induced vector potential at each edge in the mesh.supercurrent (
ndarray
) – The supercurrent density at each edge in the mesh.normal_current (
ndarray
) – The normal density at each edge in the mesh.state (
Dict
[str
,Any
]) – The solver state for the current iteration.
- class tdgl.solution.data.DynamicsData(dt, mu=None, theta=None, screening_iterations=None)[source]
A container for the measured dynamics of a TDGL solution, measured at each time step in the simulation.
- Parameters:
dt (
ndarray
) – The solver time step, \(\Delta t^{n}\).time – The solver time, a derived attribute which is equal to the cumulative sum of the time step.
theta (
Optional
[ndarray
]) – The phase of the order parameter, \(\theta=\arg\psi\)screening_iterations (
Optional
[ndarray
]) – The number of screening iterations performed at each time step.
- time_slice(tmin=-inf, tmax=inf)[source]
Returns the integer indices corresponding to the specified time window.
- voltage(i=0, j=1)[source]
Returns the voltage, i.e., the electric potential difference between probe points
i
andj
, as a function of time.
- phase_difference(i=0, j=1)[source]
Returns the phase difference between probe points
i
andj
as a function of time.
- mean_voltage(i=0, j=1, tmin=-inf, tmax=inf)[source]
Returns the time-averaged voltage \(\langle \Delta\mu \rangle\) over the specified time interval.
\[\langle V_{i,j} \rangle = \frac{\sum_n V_{i,j}^{n}\cdot\Delta t^{n}}{\sum_n\Delta t^{n}}\]- Parameters:
- Return type:
- Returns:
The time-averaged voltage over the specified time window.
- resample(num_points=None)[source]
Re-sample the dynamics to a uniform grid using linear interpolation.
- Parameters:
num_points (
Optional
[int
]) – The number of points to interpolate to.- Return type:
- Returns:
A new
DynamicsData
instance with the re-sampled data.
- plot(i=0, j=1, tmin=-inf, tmax=inf, grid=True, mean_voltage=True, labels=True, legend=False)[source]
Plot the voltage and phase difference over the specified time window.
- Parameters:
i (
int
) – Index for the first probe point.j (
int
) – Index for the second probe point.tmin (
float
) – The minimum of the time window to plot.tmax (
float
) – The maximum of the time window to plot.grid (
bool
) – Whether to add grid lines to the plots.mean_voltage (
bool
) – Whether to plot a horizontal line at the mean voltage.labels (
bool
) – Whether to include axis labels.legend (
bool
) – Whether to include a legend.
- Return type:
- Returns:
matplotlib figure and axes.
- plot_dt(tmin=-inf, tmax=inf, grid=True, labels=True, **histogram_kwargs)[source]
Plots the time step \(\Delta t^{n}\) vs. time and a histogram of \(\Delta t^{n}\).
- Parameters:
- Return type:
- Returns:
matplotlib figure and two axes.
- static from_hdf5(h5file, step_min=None, step_max=None)[source]
Load a
DynamicsData
instance from an outputh5py.File
.- Parameters:
- Return type:
- Returns:
A new
DynamicsData
instance.
- to_hdf5(h5group)[source]
Save a
DynamicsData
instance to anh5py.Group
.- Parameters:
h5group (
Group
) – An openh5py.Group
in which to save the data.- Return type:
- static from_solution(solution_path, probe_points=None, progress_bar=False)[source]
Load
DynamicsData
from the saved time steps of atdgl.Solution
.- Parameters:
- Return type:
- Returns:
A new
DynamicsData
instance
- class tdgl.BiotSavartField(supercurrent, normal_current)[source]
Bases:
NamedTuple
The magnetic field due to a current distribution, with the field due to the supercurrent and normal current labeled separately.
- tdgl.get_current_through_paths(solution_path, paths, dataset=None, interp_method='linear', units=None, with_units=True, progress_bar=True)[source]
Calculates the current through one or more paths for each saved time step.
- Parameters:
solution_path (
str
) – Path to the solution HDF5 file.paths (
Union
[ndarray
,List
[ndarray
]]) – A list of(n, 2)
arrays of(x, y)
coordinates defining the paths. A single(n, 2)
array is also allowed.dataset (
Optional
[str
]) –None
,"supercurrent"
, or"normal_current"
.None
indicates the total current.interp_method (
Literal
['linear'
,'cubic'
]) – Interpolation method: either “linear” or “cubic”.with_units (
bool
) – Whether to return apint.Quantity
with units attached.progress_bar (
bool
) – Whether to display a progress bar.
- Return type:
- Returns:
(times, currents)
, wherecurrents
is a list of arrays of the time-dependent current through each path. Ifpaths
is given as a single array,currents
will be returned as a single array.
Fluxoid Quantization
See also
tdgl.Solution.polygon_fluxoid()
, tdgl.Solution.hole_fluxoid()
, tdgl.Solution.boundary_phases()
- class tdgl.Fluxoid(flux_part, supercurrent_part)[source]
Bases:
NamedTuple
The fluxoid for a closed region \(S\) with boundary \(\partial S\) is defined as:
\[\begin{split}\begin{split} \Phi^f_S &= \Phi^f_{S,\text{ flux}} + \Phi^f_{S,\text{ supercurrent}} \\&=\int_S \mu_0 H_z(\mathbf{r})\,\mathrm{d}^2r + \oint_{\partial S} \mu_0\Lambda(\mathbf{r})\mathbf{K}_s(\mathbf{r})\cdot\mathrm{d}\mathbf{r} \\&=\oint_{\partial S} \mathbf{A}(\mathbf{r})\cdot\mathrm{d}\mathbf{r} + \oint_{\partial S} \mu_0\Lambda(\mathbf{r})\mathbf{K}_s(\mathbf{r})\cdot\mathrm{d}\mathbf{r} \end{split}\end{split}\]- Parameters:
- tdgl.make_fluxoid_polygons(device, holes=None, join_style='mitre', interp_points=None)[source]
Generates polygons enclosing the given holes to calculate the fluxoid.
- Parameters:
device (
Device
) – The Device for which to generate polygons.holes (
Union
[str
,List
[str
],None
]) – Name(s) of the hole(s) in the device for which to generate polygons. Defaults to all holes in the device.join_style (
str
) – Seetdgl.Polygon.buffer()
.interp_points (
Optional
[int
]) – If provided, the resulting polygons will be interpolated tointerp_points
vertices.
- Return type:
- Returns:
A dict of
{hole_name: fluxoid_polygon}
.