# 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:
data_range: Optional[Tuple[int, int]]

A tuple of (min_step, max_step).

supercurrent_density: Optional[ndarray]

Sheet supercurrent density, $$\mathbf{K}_s$$

normal_current_density: Optional[ndarray]

Sheet normal density, $$\mathbf{K}_n$$

tdgl_data: Optional[TDGLData]

A container for the raw TDGL data (in dimensionless units).

dynamics: Optional[DynamicsData]

A container for the time dynamics of the solution (in dimensionless units).

property saved_on_disk: bool

Returns True if the underlying HDF5 file exists on disk.

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.

property times: Optional[ndarray]

The time associated with each solve step.

closest_solve_step(time)[source]

Returns the index of the saved step closest in time to time.

Parameters:

time (float) – The time for which to find the closest index.

Return type:

int

Returns:

The index of the saved solve step whose time is closest to time

Loads the TDGL results from file for a given solve step.

Parameters:
• solve_step (int) – The step index for which to load data. Defaults to -1, i.e. the final step.

• h5file (Optional[File]) –

Return type:

None

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$$.

property field_units: str

The units in which magnetic fields are specified.

property current_units: str

The units in which currents are specified.

property time_created: datetime

The time at which the solution was originally created.

property version_info: Dict[str, str]

A dictionary of dependency versions.

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:
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.

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", or None. If None, then the total sheet current density is used.

• method (str) – Interpolation method to use, "nearest", "linear", or "cubic".

• units (Optional[str]) – The desired units for the current density. Defaults to self.current_units / self.device.length_units.

• with_units (bool) – Whether to return a pint.Quantity array with units attached.

Return type:

ndarray

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.

Parameters:
• positions (ndarray) – Shape (m, 2) array of x, y coordinates at which to evaluate the order parameter.

• method (str) – Interpolation method to use, "nearest", "linear", or "cubic".

Return type:

ndarray

Returns:

The interpolated order parameter.

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}$
Parameters:
Return type:

Fluxoid

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:
Return type:

Fluxoid

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 of Phi_0.

Parameters:

delta (bool) – If True, boundary_phases[0] will be subtracted for each polygon.

Return type:
Returns:

{polygon_name: (boundary_indices, boundary_phases)}

current_through_path(path_coords, dataset=None, method='linear', units=None, with_units=True)[source]

Calculates the total current crossing a given path.

Parameters:
Return type:

Union[float, Quantity]

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.

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 of self.field_units.

• with_units (bool) – Whether to return the fields as pint.Quantity with units attached.

• return_sum (bool) – If False, this method will return a tdgl.BiotSavartField instance, where the field from the supercurrent and normal current are identified separately.

Return type:
Returns:

An np.ndarray if return_sum is True, otherwise an instance of tdgl.BiotSavartField. If with_units is True, then the array(s) will be of type pint.Quantity. The array(s) will have shape (m, ) if vector is False, or shape (m, 3) if vector 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 a pint.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 is True, otherwise a dict of {source: potential_from_source}. If with_units is True, then the array(s) will be of type pint.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:
Return type:

None

static from_hdf5(path, solve_step=-1)[source]

Loads a tdgl.Solution from file.

Parameters:
Return type:

Solution

Returns:

delete_hdf5()[source]

Delete the HDF5 file accompanying the tdgl.Solution.

Return type:

None

equals(other, require_same_timestamp=False)[source]

Checks whether two solutions are equal.

Parameters:
Return type:

bool

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 at positions, 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:
static from_hdf5(h5file, step)[source]

Load a TDGLData instance from an output h5py.File.

Parameters:
Return type:

TDGLData

Returns:

A TDGLData instance containing data for the requested solve step.

to_hdf5(h5group)[source]

Save a TDGLData instance to an h5py.Group.

Parameters:

h5group (Group) – An open h5py.Group in which to save the data.

Return type:

None

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:
time_slice(tmin=-inf, tmax=inf)[source]

Returns the integer indices corresponding to the specified time window.

Parameters:
Return type:

ndarray

Returns:

An array of indices for the time window.

closest_time(time)[source]

Returns the index of the time step closest to time.

Parameters:

time (float) – The time for which to find the index.

Return type:

int

Returns:

The index of the time step closest to time.

voltage(i=0, j=1)[source]

Returns the voltage, i.e., the electric potential difference between probe points i and j, as a function of time.

Parameters:
Return type:

ndarray

Returns:

The voltage $$V_{ij}(t)=\mu_i(t)-\mu_j(t)$$

phase_difference(i=0, j=1)[source]

Returns the phase difference between probe points i and j as a function of time.

Parameters:
Return type:

ndarray

Returns:

The phase difference $$\Delta\theta_{ij}(t)=\theta_i(t)-\theta_j(t)$$, where $$\theta=\arg\psi$$.

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:

float

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:

DynamicsData

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:
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 output h5py.File.

Parameters:
Return type:

DynamicsData

Returns:

A new DynamicsData instance.

to_hdf5(h5group)[source]

Save a DynamicsData instance to an h5py.Group.

Parameters:

h5group (Group) – An open h5py.Group in which to save the data.

Return type:

None

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.

Parameters:
• supercurrent (ndarray) – An array of fields due to the supercurrent.

• normal_current (ndarray) – An array of fields due to the normal current.

## Fluxoid Quantization

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:
• flux_part (Union[float, Quantity]) – $$\int_S \mu_0 H_z(\mathbf{r})\,\mathrm{d}^2r=\oint_{\partial S}\mathbf{A}(\mathbf{r})\cdot\mathrm{d}\mathbf{r}$$.

• supercurrent_part (Union[float, Quantity]) – $$\oint_{\partial S}\mu_0\Lambda(\mathbf{r})\mathbf{K}_s(\mathbf{r})\cdot\mathrm{d}\mathbf{r}$$.

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:
Return type:
Returns:

A dict of {hole_name: fluxoid_polygon}.