# 2. filter_functions package¶

## 2.2. filter_functions.analytic module¶

This file provides functions for the analytical solutions to some of the dynamical decoupling sequences. Note that the filter functions given here differ by a factor of 1/omega**2 from those defined in this package due to different conventions. See for example [Cyw08]. Depending on the definition of the noise Hamiltonian one might also get different results. The functions here agree for

$B_\alpha\equiv\sigma_z/2.$

### 2.2.1. Functions¶

FID()

Free Induction Decay / Ramsey pulse

SE()

Spin Echo

PDD()

Periodic Dynamical Decoupling

CPMG()

Carr-Purcell-Meiboom-Gill Sequence

CDD()

Concatenated Dynamical Decoupling

UDD()

Uhrig Dynamical Decoupling

References

Cyw08

Cywiński, Ł., Lutchyn, R. M., Nave, C. P., & Das Sarma, S. (2008). How to enhance dephasing time in superconducting qubits. Physical Review B - Condensed Matter and Materials Physics, 77(17), 1–11. https://doi.org/10.1103/PhysRevB.77.174509

CDD(z, g)[source]
CPMG(z, n)[source]
FID(z)[source]
PDD(z, n)[source]
SE(z)[source]
UDD(z, n)[source]

## 2.3. filter_functions.basis module¶

This module defines the Basis class, a subclass of NumPy’s ndarray, to represent operator bases.

### 2.3.1. Classes¶

Basis

The operator basis as an array of shape (d**2, d, d) with d the dimension of the Hilbert space

### 2.3.2. Functions¶

normalize()

Function to normalize a Basis instance

expand()

Function to expand an array of operators in a given basis

ggm_expand()

Fast function to expand an array of operators in a Generalized Gell-Mann basis

class Basis(basis_array: Sequence, traceless: Optional[bool] = None, btype: Optional[str] = None, labels: Optional[Sequence[str]] = None)[source]

Bases: numpy.ndarray

Class for operator bases. There are several ways to instantiate a Basis object:

• by just calling this constructor with a (possibly incomplete) array of basis matrices. No checks regarding orthonormality or hermiticity are performed.

• by calling one of the classes alternative constructors (classmethods):

These bases guarantee the following properties:

• hermitian

• orthonormal

• [traceless] (can be controlled by a flag)

Since Basis is a subclass of NumPy’s ndarray, it inherits all of its attributes, e.g. shape. The following attributes behave slightly differently to a ndarray, however

• A == B is True if all elements evaluate almost equal, i.e. equivalent to np.allclose(A, B).

• basis.T transposes the last two axes of basis. For a full basis, this corresponds to transposing each element individually. For a basis element, it corresponds to normal transposition.

Parameters
basis_array: array_like, shape (n, d, d)

An array or list of square matrices that are elements of an operator basis spanning $$\mathbb{C}^{d\times d}$$. n should be smaller than or equal to d**2.

traceless: bool, optional (default: auto)

Controls whether a traceless basis is forced. Here, traceless means that the first element of the basis is the identity and the remaining elements are matrices of trace zero. If an element of basis_array is neither traceless nor the identity and traceless == True, an exception will be raised. Defaults to True if basis_array is traceless and False if not.

btype: str, optional (default: ’custom’)

A string describing the basis type. For example, a basis created by the factory method pauli() has btype ‘pauli’.

labels: sequence of str, optional

A list of labels for the individual basis elements. Defaults to ‘C_0’, ‘C_1’, …

Attributes
Other than the attributes inherited from ndarray, a Basis
instance has the following attributes:
btype: str

Basis type.

labels: sequence of str

The labels for the basis elements.

d: int

Dimension of the space spanned by the basis.

H: Basis

Hermitian conjugate.

isherm: bool

If the basis is hermitian.

isorthonorm: bool

If the basis is orthonormal.

istraceless: bool

If the basis is traceless except for an identity element

iscomplete: bool

If the basis is complete, ie spans the full space.

sparse: COO, shape (n, d, d)

Representation in the COO format supplied by the sparse package.

four_element_traces: COO, shape (n, n, n, n)

Traces over all possible combinations of four elements of self. This is required for the calculation of the error transfer matrix and thus cached in the Basis instance.

Most of the attributes above are properties which are lazily
evaluated and cached.

Methods

 Other than the methods inherited from ndarray, a Basis instance has the following methods: normalize(b) Normalizes the basis (used internally when creating a basis from elements) tidyup(eps_scale=None) Cleans up floating point errors in-place to make zeros actual zeros. eps_scale is an optional argument multiplied to the data type’s eps to get the absolute tolerance.

Constructor.

property H: filter_functions.basis.Basis

Return the basis hermitian conjugated element-wise.

property T: filter_functions.basis.Basis

Return the basis transposed element-wise.

property four_element_traces: sparse._coo.core.COO

Return all traces of the form $$\mathrm{tr}(C_i C_j C_k C_l)$$ as a sparse COO array for $$i,j,k,l > 0$$ (i.e. excluding the identity).

classmethod from_partial(partial_basis_array: Sequence, traceless: Optional[bool] = None, btype: Optional[str] = None, labels: Optional[Sequence[str]] = None)filter_functions.basis.Basis[source]

Generate complete and orthonormal basis from a partial set.

The basis is completed using singular value decomposition to determine the null space of the expansion coefficients of the partial basis with respect to another complete basis.

Parameters
partial_basis_array: array_like

A sequence of basis elements.

traceless: bool, optional

If a traceless basis should be generated (i.e. the first element is the identity and all the others have trace zero).

btype: str, optional

A custom identifier.

labels: Sequence[str], optional

A list of custom labels for each element. If len(labels) == len(partial_basis_array), the newly created elements get labels ‘C_i’.

Returns
basis: Basis, shape (d**2, d, d)

The orthonormal basis.

Raises
ValueError

If the given elements are not orthonormal.

ValueError

If the given elements are not traceless but traceless==True.

ValueError

If not len(partial_basis_array) or d**2 labels were given.

classmethod ggm(d: int)filter_functions.basis.Basis[source]

Returns a generalized Gell-Mann basis in $$d$$ dimensions [Bert08] where the elements $$\Lambda_i$$ are normalized with respect to the Hilbert-Schmidt inner product,

$\begin{split}\langle\Lambda_i,\Lambda_j\rangle &= \mathrm{Tr}\,\Lambda_i^\dagger\Lambda_j \\ &= \delta_{ij}.\end{split}$
Parameters
d: int

The dimensionality of the space spanned by the basis

Returns
basis: Basis

The Basis object representing the GGM.

References

Bert08

Bertlmann, R. A., & Krammer, P. (2008). Bloch vectors for qudits. Journal of Physics A: Mathematical and Theoretical, 41(23). https://doi.org/10.1088/1751-8113/41/23/235303

property iscomplete: bool

Returns True if basis is complete.

property isherm: bool

Returns True if all basis elements are hermitian.

property isorthonorm: bool

Returns True if basis is orthonormal.

property istraceless: bool

Returns True if basis is traceless except for possibly the identity.

normalize(copy: bool = False)Union[None, filter_functions.basis.Basis][source]

Normalize the basis.

classmethod pauli(n: int)filter_functions.basis.Basis[source]

Returns a Pauli basis for $$n$$ qubits, i.e. the basis spans the space $$\mathbb{C}^{d\times d}$$ with $$d = 2^n$$:

$\mathcal{P} = \{I, X, Y, Z\}^{\otimes n}.$

The elements $$\sigma_i$$ are normalized with respect to the Hilbert-Schmidt inner product,

$\begin{split}\langle\sigma_i,\sigma_j\rangle &= \mathrm{Tr}\,\sigma_i^\dagger\sigma_j \\ &= \delta_{ij}.\end{split}$
Parameters
n: int

The number of qubits.

Returns
basis: Basis

The Basis object representing the Pauli basis.

property sparse: sparse._coo.core.COO

Return the basis as a sparse COO array

tidyup(eps_scale: Optional[float] = None)None[source]

Wraps util.remove_float_errors.

expand(M: , basis: , normalized: bool = True, hermitian: bool = False, tidyup: bool = False)numpy.ndarray[source]

Expand the array M in the basis given by basis.

Parameters
M: array_like

The square matrix (d, d) or array of square matrices (…, d, d) to be expanded in basis

basis: array_like

The basis of shape (m, d, d) in which to expand.

normalized: bool {True}

Wether the basis is normalized.

hermitian: bool (default: False)

If M is hermitian along its last two axes, the result will be real.

tidyup: bool {False}

Whether to set values below the floating point eps to zero.

Returns
coefficients: ndarray

The coefficient array with shape (…, m) or (m,) if M was 2-d

Notes

For an orthogonal matrix basis $$\mathcal{C} = \big\{C_k\in \mathbb{C}^{d\times d}: \langle C_k,C_l\rangle_\mathrm{HS} \propto \delta_{kl}\big\}_{k=0}^{d^2-1}$$ with the Hilbert-Schmidt inner product as implemented by dot_HS() and $$M\in\mathbb{C}^{d\times d}$$, the expansion of $$M$$ in terms of $$\mathcal{C}$$ is given by

$\begin{split}M &= \sum_j c_j C_j, \\ c_j &= \frac{\mathrm{tr}\big(M C_j\big)} {\mathrm{tr}\big(C_j^\dagger C_j\big)}.\end{split}$
ggm_expand(M: , traceless: bool = False, hermitian: bool = False)numpy.ndarray[source]

Expand the matrix M in a Generalized Gell-Mann basis [Bert08]. This function makes use of the explicit construction prescription of the basis and thus makes do without computing the expansion coefficients as the overlap between the matrix and each basis element.

Parameters
M: array_like

The square matrix (d, d) or array of square matrices (…, d, d) to be expanded in a GGM basis.

traceless: bool (default: False)

Include the basis element proportional to the identity in the expansion. If it is known beforehand that M is traceless, the corresponding coefficient is zero and thus doesn’t need to be computed.

hermitian: bool (default: False)

If M is hermitian along its last two axes, the result will be real.

Returns
coefficients: ndarray

The real coefficient array with shape (d**2,) or (…, d**2)

References

Bert08

Bertlmann, R. A., & Krammer, P. (2008). Bloch vectors for qudits. Journal of Physics A: Mathematical and Theoretical, 41(23). https://doi.org/10.1088/1751-8113/41/23/235303

normalize(b: filter_functions.basis.Basis)filter_functions.basis.Basis[source]

Return a copy of the basis b normalized with respect to the Frobenius norm [Gol85]:

$$||A||_F = \left[\sum_{i,j} |a_{i,j}|^2\right]^{1/2}$$

or equivalently, with respect to the Hilbert-Schmidt inner product as implemented by dot_HS().

References

Gol85

G. H. Golub and C. F. Van Loan, Matrix Computations, Baltimore, MD, Johns Hopkins University Press, 1985, pg. 15

## 2.4. filter_functions.numeric module¶

This module defines the functions to calculate everything related to filter functions.

### 2.4.1. Functions¶

calculate_control_matrix_from_atomic()

Calculate the control matrix from those of atomic pulse sequences

calculate_control_matrix_from_scratch()

Calculate the control matrix from scratch

calculate_control_matrix_periodic()

Calculate the control matrix for a periodic Hamiltonian

calculate_noise_operators_from_atomic()

Calculate the interaction picture noise operators from atomic segments. Same calculation as calculate_control_matrix_from_atomic() except in Hilbert space.

calculate_noise_operators_from_scratch()

Calculate the interaction picture noise operators from scratch. Same calculation as calculate_control_matrix_from_scratch() except in Hilbert space.

calculate_cumulant_function()

Calculate the cumulant function for a given PulseSequence object.

calculate_decay_amplitudes()

Calculate the decay amplitudes, corresponding to first order terms of the Magnus expansion

calculate_frequency_shifts()

Calculate the frequency shifts, corresponding to second order terms of the Magnus expansion

calculate_filter_function()

Calculate the filter function from the control matrix

calculate_second_order_filter_function()

Calculate the second order filter function used to compute the frequency shifts.

calculate_pulse_correlation_filter_function()

Calculate the pulse correlation filter function from the control matrix

diagonalize()

Diagonalize a Hamiltonian

error_transfer_matrix()

Calculate the error transfer matrix of a pulse up to a unitary rotation

infidelity()

Function to compute the infidelity of a pulse defined by a PulseSequence instance for a given noise spectral density and frequencies

calculate_control_matrix_from_atomic(phases: numpy.ndarray, control_matrix_atomic: numpy.ndarray, propagators_liouville: numpy.ndarray, show_progressbar: bool = False, which: str = 'total')numpy.ndarray[source]

Calculate the control matrix from the control matrices of atomic segments.

Parameters
phases: array_like, shape (n_dt, n_omega)

The phase factors for $$g\in\{0, 1, \dots, G-1\}$$.

control_matrix_atomic: array_like, shape (n_dt, n_nops, d**2, n_omega)

The pulse control matrices for $$g\in\{1, 2, \dots, G\}$$.

propagators_liouville: array_like, shape (n_dt, n_nops, d**2, d**2)

The transfer matrices of the cumulative propagators for $$g\in\{0, 1, \dots, G-1\}$$.

show_progressbar: bool, optional

Show a progress bar for the calculation.

which: str, (‘total’, ‘correlations’)

Compute the total control matrix (the sum of all time steps) or the correlation control matrix (first axis holds each pulses’ contribution)

Returns
control_matrix: ndarray, shape ([n_pls,] n_nops, d**2, n_omega)

The control matrix $$\tilde{\mathcal{B}}(\omega)$$.

calculate_control_matrix_from_scratch

Control matrix from scratch.

liouville_representation

Liouville representation for a given basis.

Notes

The control matrix is calculated by evaluating the sum

$\tilde{\mathcal{B}}(\omega) = \sum_{g=1}^G e^{i\omega t_{g-1}} \tilde{\mathcal{B}}^{(g)}(\omega)\mathcal{Q}^{(g-1)}.$
calculate_control_matrix_from_scratch(eigvals: numpy.ndarray, eigvecs: numpy.ndarray, propagators: numpy.ndarray, omega: Sequence[float], basis: filter_functions.basis.Basis, n_opers: Sequence[Union[numpy.ndarray, qutip.qobj.Qobj]], n_coeffs: Sequence[Sequence[float]], dt: Sequence[float], t: Optional[Sequence[float]] = None, show_progressbar: bool = False, cache_intermediates: bool = False, out: Optional[numpy.ndarray] = None)Union[numpy.ndarray, Tuple[numpy.ndarray, Dict[str, numpy.ndarray]]][source]

Calculate the control matrix from scratch, i.e. without knowledge of the control matrices of more atomic pulse sequences.

Parameters
eigvals: array_like, shape (n_dt, d)

Eigenvalue vectors for each time pulse segment g with the first axis counting the pulse segment, i.e. eigvals == array([D_0, D_1, ...]).

eigvecs: array_like, shape (n_dt, d, d)

Eigenvector matrices for each time pulse segment g with the first axis counting the pulse segment, i.e. eigvecs == array([V_0, V_1, ...]).

propagators: array_like, shape (n_dt+1, d, d)

The propagators $$Q_g = P_g P_{g-1}\cdots P_0$$ as a (d, d) array with d the dimension of the Hilbert space.

omega: array_like, shape (n_omega,)

Frequencies at which the pulse control matrix is to be evaluated.

basis: Basis, shape (d**2, d, d)

The basis elements in which the pulse control matrix will be expanded.

n_opers: array_like, shape (n_nops, d, d)

Noise operators $$B_\alpha$$.

n_coeffs: array_like, shape (n_nops, n_dt)

The sensitivities of the system to the noise operators given by n_opers at the given time step.

dt: array_like, shape (n_dt)

Sequence duration, i.e. for the $$g$$-th pulse $$t_g - t_{g-1}$$.

t: array_like, shape (n_dt+1), optional

The absolute times of the different segments. Can also be computed from dt.

show_progressbar: bool, optional

Show a progress bar for the calculation.

cache_intermediates: bool, optional

Keep and return intermediate terms of the calculation that can be reused in other computations (second order or gradients). Otherwise the sum is performed in-place. The default is False.

out: ndarray, optional

A location into which the result is stored. See numpy.ufunc().

Returns
control_matrix: ndarray, shape (n_nops, d**2, n_omega)

The control matrix $$\tilde{\mathcal{B}}(\omega)$$

intermediates: dict[str, ndarray]

Intermediate results of the calculation. Only if cache_intermediates is True.

calculate_control_matrix_from_atomic

Control matrix from concatenation.

calculate_control_matrix_periodic

Control matrix for periodic system.

Notes

The control matrix is calculated according to

$\tilde{\mathcal{B}}_{\alpha k}(\omega) = \sum_{g=1}^G e^{i\omega t_{g-1}} s_\alpha^{(g)}\mathrm{tr}\left( [\bar{B}_\alpha^{(g)}\circ I(\omega)] \bar{C}_k^{(g)} \right)$

where

$\begin{split}I^{(g)}_{nm}(\omega) &= \int_0^{t_l - t_{g-1}}\mathrm{d}t\, e^{i(\omega+\omega_n-\omega_m)t} \\ &= \frac{e^{i(\omega+\omega_n-\omega_m) (t_l - t_{g-1})} - 1} {i(\omega+\omega_n-\omega_m)}, \\ \bar{B}_\alpha^{(g)} &= V^{(g)\dagger} B_\alpha V^{(g)}, \\ \bar{C}_k^{(g)} &= V^{(g)\dagger} Q_{g-1} C_k Q_{g-1}^\dagger V^{(g)},\end{split}$

and $$V^{(g)}$$ is the matrix of eigenvectors that diagonalizes $$\tilde{\mathcal{H}}_n^{(g)}$$, $$B_\alpha$$ the $$\alpha$$-th noise operator $$s_\alpha^{(g)}$$ the noise sensitivity during interval $$g$$, and $$C_k$$ the $$k$$-th basis element.

calculate_control_matrix_periodic(phases: numpy.ndarray, control_matrix: numpy.ndarray, total_propagator_liouville: numpy.ndarray, repeats: int, check_invertible: bool = True)numpy.ndarray[source]

Calculate the control matrix of a periodic pulse given the phase factors, control matrix and transfer matrix of the total propagator, total_propagator_liouville, of the atomic pulse.

Parameters
phases: ndarray, shape (n_omega,)

The phase factors $$e^{i\omega T}$$ of the atomic pulse.

control_matrix: ndarray, shape (n_nops, d**2, n_omega)

The control matrix $$\tilde{\mathcal{B}}^{(1)}(\omega)$$ of the atomic pulse.

total_propagator_liouville: ndarray, shape (d**2, d**2)

The transfer matrix $$\mathcal{Q}^{(1)}$$ of the atomic pulse.

repeats: int

The number of repetitions.

check_invertible: bool, optional

Check for all frequencies if the inverse $$\mathbb{I} - e^{i\omega T} \mathcal{Q}^{(1)}$$ exists by calculating the determinant. If it does not exist, the sum is evaluated explicitly for those points. The default is True.

Returns
control_matrix: ndarray, shape (n_nops, d**2, n_omega)

The control matrix $$\tilde{\mathcal{B}}(\omega)$$ of the repeated pulse.

Notes

The control matrix is computed as

$\begin{split}\tilde{\mathcal{B}}(\omega) &= \tilde{\mathcal{B}}^{(1)}(\omega) \sum_{g=0}^{G-1} \left(e^{i\omega T}\right)^g \\ &= \tilde{\mathcal{B}}^{(1)}(\omega)\bigl( \mathbb{I} - e^{i\omega T} \mathcal{Q}^{(1)}\bigr)^{-1}\bigl( \mathbb{I} - \bigl(e^{i\omega T} \mathcal{Q}^{(1)}\bigr)^G\bigr).\end{split}$

with $$G$$ the number of repetitions.

calculate_cumulant_function(pulse: PulseSequence, spectrum: Optional[numpy.ndarray] = None, omega: Optional[Sequence[float]] = None, n_oper_identifiers: Optional[Sequence[str]] = None, which: str = 'total', second_order: bool = False, decay_amplitudes: Optional[numpy.ndarray] = None, frequency_shifts: Optional[numpy.ndarray] = None, show_progressbar: bool = False, memory_parsimonious: bool = False, cache_intermediates: Optional[bool] = None)numpy.ndarray[source]

Calculate the cumulant function $$\mathcal{K}(\tau)$$.

The error transfer matrix is obtained from the cumulant function by exponentiation, $$\langle\tilde{\mathcal{U}}\rangle = \exp\mathcal{K}(\tau)$$.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to compute the cumulant function.

spectrum: array_like, shape ([[n_nops,] n_nops,] n_omega), optional

The noise power spectral density in units of inverse frequencies as an array of shape (n_omega,), (n_nops, n_omega), or (n_nops, n_nops, n_omega). In the first case, the same spectrum is taken for all noise operators, in the second, it is assumed that there are no correlations between different noise sources and thus there is one spectrum for each noise operator. In the third and most general case, there may be a spectrum for each pair of noise operators corresponding to the correlations between them. n_nops is the number of noise operators considered and should be equal to len(n_oper_identifiers).

omega: array_like, shape (n_omega,), optional

The frequencies at which to evaluate the filter functions.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which to evaluate the cumulant function. The default is all.

which: str, optional

Which decay amplitudes should be calculated, may be either ‘total’ (default) or ‘correlations’. See infidelity() and Notes. Note that the latter is not available for the second order terms.

second_order: bool, optional

Also take into account the frequency shifts $$\Delta$$ that correspond to second order Magnus expansion and constitute unitary terms. Default False.

decay_amplitudes, array_like, shape ([[n_pls, n_pls,] n_nops,] n_nops, d**2, d**2), optional

A precomputed cumulant function. If given, spectrum, omega are not required.

frequency_shifts, array_like, shape ([[n_pls, n_pls,] n_nops,] n_nops, d**2, d**2), optional

A precomputed frequency shift. If given, spectrum, omega are not required for second order terms.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

memory_parsimonious: bool, optional

Trade memory footprint for performance. See calculate_decay_amplitudes(). The default is False.

cache_intermediates: bool, optional

Keep and return intermediate terms of the calculation of the control matrix that can be reused in other computations (second order or gradients). Otherwise the sum is performed in-place. Default is True if second_order=True, else False.

Returns
cumulant_function: ndarray, shape ([[n_pls, n_pls,] n_nops,] n_nops, d**2, d**2)

The cumulant function. The individual noise operator contributions chosen by n_oper_identifiers are on the third to last axis / axes, depending on whether the noise is cross-correlated or not. If which == 'correlations', the first two axes correspond to the contributions of the pulses in the sequence.

calculate_decay_amplitudes

Calculate the $$\Gamma_{\alpha\beta,kl}$$

error_transfer_matrix

Calculate the error transfer matrix $$\exp\mathcal{K}$$.

infidelity

Calculate only infidelity of a pulse.

pulse_sequence.concatenate

Concatenate PulseSequence objects.

calculate_pulse_correlation_filter_function

Notes

The cumulant function is given by

$\begin{split}K_{\alpha\beta,ij}(\tau) = -\frac{1}{2} \sum_{kl}\biggl( &\Delta_{\alpha\beta,kl}\left( T_{klji} - T_{lkji} - T_{klij} + T_{lkij} \right) \\ + &\Gamma_{\alpha\beta,kl}\left( T_{klji} - T_{kjli} - T_{kilj} + T_{kijl} \right) \biggr)\end{split}$

Here, $$T_{ijkl} = \mathrm{tr}(C_i C_j C_k C_l)$$ is a trivial function of the basis elements $$C_i$$, and $$\Gamma_{\alpha\beta,kl}$$ and $$\Delta_{\alpha\beta,kl}$$ are the decay amplitudes and frequency shifts which correspond to first and second order in the Magnus expansion, respectively. Since the latter induce coherent errors, we can approximately neglect them if we assume that the pulse has been experimentally calibrated.

For a single qubit and represented in the Pauli basis, the above reduces to

$\begin{split}K_{\alpha\beta,ij}(\tau) = \begin{cases} - \sum_{k\neq i}\Gamma_{\alpha\beta,kk} &\quad\mathrm{if}\: i = j, \\ - \Delta_{\alpha\beta,ij} + \Delta_{\alpha\beta,ji} + \Gamma_{\alpha\beta,ij} &\quad\mathrm{if}\: i\neq j, \end{cases}\end{split}$

for $$i\in\{1,2,3\}$$ and $$K_{0j} = K_{i0} = 0$$.

Lastly, the pulse correlation cumulant function resolves correlations in the cumulant function of a sequence of pulses $$g = 1,\dotsc,G$$ such that the following holds:

$K_{\alpha\beta,ij}(\tau) = \sum_{g,g'=1}^G K_{\alpha\beta,ij}^{(gg')}(\tau).$
calculate_decay_amplitudes(pulse: PulseSequence, spectrum: numpy.ndarray, omega: Sequence[float], n_oper_identifiers: Optional[Sequence[str]] = None, which: str = 'total', show_progressbar: bool = False, cache_intermediates: bool = False, memory_parsimonious: bool = False)numpy.ndarray[source]

Get the decay amplitudes $$\Gamma_{\alpha\beta, kl}$$ for noise sources $$\alpha,\beta$$ and basis elements $$k,l$$.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to compute the decay amplitudes.

spectrum: array_like, shape ([[n_nops,] n_nops,] n_omega)

The noise power spectral density. If 1-d, the same spectrum is used for all noise operators. If 2-d, one (self-) spectrum for each noise operator is expected. If 3-d, should be a matrix of cross-spectra such that spectrum[i, j] == spectrum[j, i].conj().

omega: array_like,

The frequencies at which to calculate the filter functions.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which to calculate the decay amplitudes. The default is all.

which: str, optional

Which decay amplitudes should be calculated, may be either ‘total’ (default) or ‘correlations’. See infidelity() and Notes.

show_progressbar: bool, optional

Show a progress bar for the calculation.

cache_intermediates: bool, optional

Keep and return intermediate terms of the calculation that are useful in other places (if control matrix not already cached).

memory_parsimonious: bool, optional

For large dimensions, the integrand

$\tilde{\mathcal{B}}^\ast_{\alpha k}(\omega) S_{\alpha\beta}(\omega)\tilde{\mathcal{B}}_{\beta l}(\omega)$

can consume quite a large amount of memory if set up for all $$\alpha,\beta,k,l$$ at once. If True, it is only set up and integrated for a single $$k$$ at a time and looped over. This is slower but requires much less memory. The default is False.

Returns
decay_amplitudes: ndarray, shape ([[n_pls, n_pls,] n_nops,] n_nops, d**2, d**2)

The decay amplitudes.

Raises
ValueError

If spectrum has incompatible shape.

infidelity

Compute the infidelity directly.

pulse_sequence.concatenate

Concatenate PulseSequence objects.

calculate_frequency_shifts

Second order (unitary) terms.

calculate_pulse_correlation_filter_function

Notes

The total decay amplitudes are given by

$\Gamma_{\alpha\beta, kl} = \int\frac{\mathrm{d}\omega}{2\pi} \tilde{\mathcal{B}}^\ast_{\alpha k}(\omega) S_{\alpha\beta}(\omega)\tilde{\mathcal{B}}_{\beta l}(\omega).$

If pulse correlations are taken into account, they are given by

$\Gamma_{\alpha\beta, kl}^{(gg')} = \int \frac{\mathrm{d}\omega}{2\pi} S_{\alpha\beta}(\omega) F_{\alpha\beta, kl}^{(gg')}(\omega).$
calculate_filter_function(control_matrix: numpy.ndarray, which: str = 'fidelity')numpy.ndarray[source]

Compute the filter function from the control matrix.

Parameters
control_matrix: array_like, shape (n_nops, d**2, n_omega)

The control matrix.

whichstr, optional

Which filter function to return. Either ‘fidelity’ (default) or ‘generalized’ (see Notes).

Returns
filter_function: ndarray, shape (n_nops, n_nops, [d**2, d**2,] n_omega)

The filter functions for each noise operator correlation. The diagonal corresponds to the filter functions for uncorrelated noise sources.

Notes

The generalized filter function is given by

$F_{\alpha\beta,kl}(\omega) = \tilde{\mathcal{B}}_{\alpha k}^\ast(\omega) \tilde{\mathcal{B}}_{\beta l}(\omega),$

where $$\alpha,\beta$$ are indices counting the noise operators $$B_\alpha$$ and $$k,l$$ indices counting the basis elements $$C_k$$.

The fidelity filter function is obtained by tracing over the basis indices:

$F_{\alpha\beta}(\omega) = \sum_{k} F_{\alpha\beta,kk}(\omega).$
calculate_frequency_shifts(pulse: PulseSequence, spectrum: numpy.ndarray, omega: Sequence[float], n_oper_identifiers: Optional[Sequence[str]] = None, show_progressbar: bool = False)numpy.ndarray[source]

Get the frequency shifts $$\Delta_{\alpha\beta, kl}$$ for noise sources $$\alpha,\beta$$ and basis elements $$k,l$$.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to compute the frequency shifts.

spectrum: array_like, shape ([[n_nops,] n_nops,] n_omega)

The two-sided noise power spectral density. If 1-d, the same spectrum is used for all noise operators. If 2-d, one (self-) spectrum for each noise operator is expected. If 3-d, should be a matrix of cross-spectra such that spectrum[i, j] == spectrum[j, i].conj().

omega: array_like,

The frequencies. Note that the frequencies are assumed to be symmetric about zero.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which to calculate the frequency shifts. The default is all.

show_progressbar: bool, optional

Show a progress bar for the calculation.

Returns
Delta: ndarray, shape ([n_nops,] n_nops, d**2, d**2)

The frequency shifts.

Raises
ValueError

If spectrum has incompatible shape.

calculate_second_order_filter_function

Corresponding filter function.

calculate_decay_amplitudes

First order (dissipative) terms.

infidelity

Compute the infidelity directly.

pulse_sequence.concatenate

Concatenate PulseSequence objects.

calculate_pulse_correlation_filter_function

Notes

The total frequency shifts are given by

$\Delta_{\alpha\beta, kl} = \int_{-\infty}^\infty \frac{\mathrm{d}{\omega}}{2\pi} S_{\alpha\beta}(\omega) F_{\alpha\beta,kl}^{(2)}(\omega)$

with $$F_{\alpha\beta,kl}^{(2)}(\omega)$$ the second order filter function.

calculate_noise_operators_from_atomic(phases: numpy.ndarray, noise_operators_atomic: numpy.ndarray, propagators: numpy.ndarray, show_progressbar: bool = False)numpy.ndarray[source]

Calculate the interaction picutre noise operators from atomic segments.

Parameters
phases: array_like, shape (n_dt, n_omega)

The phase factors for $$g\in\{0, 1, \dots, G-1\}$$.

noise_operators_atomic: array_like, shape (n_dt, n_nops, d, d, n_omega)

The noise operators in the interaction picture of the g-th pulse, i.e. for $$g\in\{1, 2, \dots, G\}$$.

propagators: array_like, shape (n_dt, d, d)

The cumulative propagators of the pulses $$g\in\{0, 1, \dots, G-1\}$$.

show_progressbar: bool, optional

Show a progress bar for the calculation.

Returns
noise_operators: ndarray, shape (n_omega, n_nops, d, d)

The interaction picture noise operators $$\tilde{B}_\alpha(\omega)$$.

calculate_noise_operators_from_scratch

Compute the operators from scratch.

calculate_control_matrix_from_atomic

Same calculation in Liouville space.

Notes

The noise operators are calculated by evaluating the sum

$\tilde{B}_\alpha(\omega) = \sum_{g=1}^G e^{i\omega t_{g-1}} Q_{g-1}^\dagger\tilde{B}_\alpha^{(g)}(\omega) Q_{g-1}.$

The control matrix then corresponds to the coefficients of expansion in an operator basis $$\{C_k\}_k$$:

$\tilde{\mathcal{B}}_{k\alpha}(\omega) = \mathrm{tr}(\tilde{B}_\alpha(\omega) C_k).$

Due to differences in implementation (for performance reasons), the axes of the result are transposed compared to the control matrix:

>>> ctrlmat = calculate_control_matrix_from_atomic(...)
>>> ctrlmat.shape
(n_nops, d**2, n_omega)
>>> noiseops = calculate_noise_operators_from_atomic(...)
>>> noiseops.shape
(n_omega, n_nops, d, d)
>>> ctrlmat_from_noiseops = ff.basis.expand(noiseops, basis)
>>> np.allclose(ctrlmat, ctrlmat_from_noiseops.transpose(1, 2, 0))
True

calculate_noise_operators_from_scratch(eigvals: numpy.ndarray, eigvecs: numpy.ndarray, propagators: numpy.ndarray, omega: Sequence[float], n_opers: Sequence[Union[numpy.ndarray, qutip.qobj.Qobj]], n_coeffs: Sequence[Sequence[float]], dt: Sequence[float], t: Optional[Sequence[float]] = None, show_progressbar: bool = False, cache_intermediates: bool = False)Union[numpy.ndarray, Tuple[numpy.ndarray, Dict[str, numpy.ndarray]]][source]

Calculate the noise operators in interaction picture from scratch.

Parameters
eigvals: array_like, shape (n_dt, d)

Eigenvalue vectors for each time pulse segment g with the first axis counting the pulse segment, i.e. eigvals == array([D_0, D_1, ...]).

eigvecs: array_like, shape (n_dt, d, d)

Eigenvector matrices for each time pulse segment g with the first axis counting the pulse segment, i.e. eigvecs == array([V_0, V_1, ...]).

propagators: array_like, shape (n_dt+1, d, d)

The propagators $$Q_g = P_g P_{g-1}\cdots P_0$$ as a (d, d) array with d the dimension of the Hilbert space.

omega: array_like, shape (n_omega,)

Frequencies at which the pulse control matrix is to be evaluated.

n_opers: array_like, shape (n_nops, d, d)

Noise operators $$B_\alpha$$.

n_coeffs: array_like, shape (n_nops, n_dt)

The sensitivities of the system to the noise operators given by n_opers at the given time step.

dt: array_like, shape (n_dt)

Sequence duration, i.e. for the $$g$$-th pulse $$t_g - t_{g-1}$$.

t: array_like, shape (n_dt+1), optional

The absolute times of the different segments. Can also be computed from dt.

show_progressbar: bool, optional

Show a progress bar for the calculation.

cache_intermediates: bool, optional

Keep and return intermediate terms of the calculation that can be reused in other computations (second order or gradients). Otherwise the sum is performed in-place. The default is False.

Returns
noise_operators: ndarray, shape (n_omega, n_nops, d, d)

The interaction picture noise operators $$\tilde{B}_\alpha(\omega)$$.

intermediates: dict[str, ndarray]

Intermediate results of the calculation. Only if cache_intermediates is True.

calculate_noise_operators_from_atomic

Compute the operators from atomic segments.

calculate_control_matrix_from_scratch

Same calculation in Liouville space.

Notes

The interaction picture noise operators are calculated according to

$\tilde{B}_\alpha(\omega) = \sum_{g=1}^G e^{i\omega t_{g-1}} s_\alpha^{(g)} P^{(g)\dagger}\left[ \bar{B}^{(g)}_\alpha \circ I^{(g)}(\omega) \right] P^{(g)}$

where

$\begin{split}I^{(g)}_{nm}(\omega) &= \int_0^{t_g - t_{g-1}}\mathrm{d}t\, e^{i(\omega+\omega_n-\omega_m)t} \\ &= \frac{e^{i(\omega+\omega_n-\omega_m) (t_g - t_{g-1})} - 1} {i(\omega+\omega_n-\omega_m)}, \\ \bar{B}_\alpha^{(g)} &= V^{(g)\dagger} B_\alpha V^{(g)}, \\ P^{(g)} &= V^{(g)\dagger} Q_{g-1},\end{split}$

and $$V^{(g)}$$ is the matrix of eigenvectors that diagonalizes $$\tilde{\mathcal{H}}_n^{(g)}$$, $$B_\alpha$$ the $$\alpha$$-th noise operator, and $$s_\alpha^{(g)}$$ the noise sensitivity during interval $$g$$.

The control matrix then corresponds to the coefficients of expansion in an operator basis $$\{C_k\}_k$$:

$\tilde{\mathcal{B}}_{k\alpha}(\omega) = \mathrm{tr}(\tilde{B}_\alpha(\omega) C_k).$

Due to differences in implementation (for performance reasons), the axes of the result are transposed compared to the control matrix:

>>> ctrlmat = calculate_control_matrix_from_scratch(...)
>>> ctrlmat.shape
(n_nops, d**2, n_omega)
>>> noiseops = calculate_noise_operators_from_scratch(...)
>>> noiseops.shape
(n_omega, n_nops, d, d)
>>> ctrlmat_from_noiseops = ff.basis.expand(noiseops, basis)
>>> np.allclose(ctrlmat, ctrlmat_from_noiseops.transpose(1, 2, 0))
True

calculate_pulse_correlation_filter_function(control_matrix: numpy.ndarray, which: str = 'fidelity')numpy.ndarray[source]

Compute pulse correlation filter function from control matrix.

Parameters
control_matrix: array_like, shape (n_pulses, n_nops, d**2, n_omega)

The control matrix.

whichstr, optional

Which filter function to return. Either ‘fidelity’ (default) or ‘generalized’ (see Notes).

Returns
filter_function_pc: ndarray, shape (n_pls, n_pls, n_nops, n_nops, [d**2, d**2,] n_omega)

The pulse correlation filter functions for each pulse and noise operator correlations. The first two axes hold the pulse correlations, the second two the noise correlations.

whichstr, optional

Which filter function to return. Either ‘fidelity’ (default) or ‘generalized’ (see Notes).

calculate_control_matrix_from_scratch

Control matrix from scratch.

calculate_control_matrix_from_atomic

Control matrix from concatenation.

calculate_filter_function

Regular filter function.

Notes

The generalized pulse correlation filter function is given by

$F_{\alpha\beta,kl}^{(gg')}(\omega) = \bigl[ \mathcal{Q}^{(g'-1)\dagger} \tilde{\mathcal{B}}^{(g')\dagger}(\omega) \bigr]_{k\alpha} \bigl[ \tilde{\mathcal{B}}^{(g)}(\omega)\mathcal{Q}^{(g-1)} \bigr]_{\beta l} e^{i\omega(t_{g-1} - t_{g'-1})},$

with $$\tilde{\mathcal{B}}^{(g)}$$ the control matrix of the $$g$$-th pulse. The fidelity pulse correlation function is obtained by tracing out the basis indices,

$F_{\alpha\beta}^{(gg')}(\omega) = \sum_{k} F_{\alpha\beta,kk}^{(gg')}(\omega)$
calculate_second_order_filter_function(eigvals: numpy.ndarray, eigvecs: numpy.ndarray, propagators: numpy.ndarray, omega: Sequence[float], basis: filter_functions.basis.Basis, n_opers: Sequence[Union[numpy.ndarray, qutip.qobj.Qobj]], n_coeffs: Sequence[Sequence[float]], dt: Sequence[float], intermediates: Optional[Dict[str, numpy.ndarray]] = None, show_progressbar: bool = False)numpy.ndarray[source]

Calculate the second order filter function for frequency shifts.

Parameters
eigvals: array_like, shape (n_dt, d)

Eigenvalue vectors for each time pulse segment l with the first axis counting the pulse segment, i.e. eigvals == array([D_0, D_1, ...]).

eigvecs: array_like, shape (n_dt, d, d)

Eigenvector matrices for each time pulse segment l with the first axis counting the pulse segment, i.e. eigvecs == array([V_0, V_1, ...]).

propagators: array_like, shape (n_dt+1, d, d)

The propagators $$Q_l = P_l P_{l-1}\cdots P_0$$ as a (d, d) array with d the dimension of the Hilbert space.

omega: array_like, shape (n_omega,)

Frequencies at which the pulse control matrix is to be evaluated.

basis: Basis, shape (d**2, d, d)

The basis elements in which the pulse control matrix will be expanded.

n_opers: array_like, shape (n_nops, d, d)

Noise operators $$B_\alpha$$.

n_coeffs: array_like, shape (n_nops, n_dt)

The sensitivities of the system to the noise operators given by n_opers at the given time step.

dt: array_like, shape (n_dt)

Sequence duration, i.e. for the $$l$$-th pulse $$t_l - t_{l-1}$$.

intermediates: Dict[str, ndarray], optional

Intermediate terms of the calculation of the control matrix that can be reused here. If None (default), they are computed from scratch.

show_progressbar: bool, optional

Show a progress bar for the calculation.

Returns
second_order_filter_function: ndarray, shape (n_nops, n_nops, d**2, d**2, n_omega)

The second order filter function.

calculate_frequency_shifts

Integrate over filter function times spectrum.

calculate_decay_amplitudes

First order (dissipative) terms.

infidelity

Compute the infidelity directly.

pulse_sequence.concatenate

Concatenate PulseSequence objects.

calculate_pulse_correlation_filter_function

Notes

The second order filter function is given by

$F_{\alpha\beta, kl}^{(2)} = \sum_{g=1}^G\left[ \mathcal{G}_{\alpha k}^{(g)\ast}(\omega) \sum_{g'=1}^{g-1}\mathcal{G}_{\beta l}^{(g')}(\omega) + s_\alpha^{(g)}\bar{B}_{\alpha,ij}^{(g)}\bar{C}_{k,ji}^{(g)} I_{ijmn}^{(g)}(\omega)\bar{C}_{l,nm}^{(g) \bar{B}_{\beta,mn}^{(g)}s_\beta^{(g)}} \right]$

with

$\begin{split}\mathcal{G}^{(g)}(\omega) &= e^{i\omega t_{g-1}}\mathcal{B}^{(g)}(\omega) \mathcal{Q}^{(g-1)}, \\ I_{ijmn}^{(g)}(\omega) &= \int_{t_{g-1}}^{t_g}\mathrm{d}{t} e^{i\Omega_{ij}^{(g)}(t - t_{g-1}) - i\omega t} \int_{t_{g-1}}^{t}\mathrm{d}{t'} e^{i\Omega_{mn}^{(g)}(t' - t_{g-1}) + i\omega t'}.\end{split}$
diagonalize(hamiltonian: numpy.ndarray, dt: Sequence[float])Tuple[numpy.ndarray][source]

Diagonalize a Hamiltonian.

Diagonalize the Hamiltonian which is piecewise constant during the times given by dt and return eigenvalues, eigenvectors, and the cumulative propagators $$Q_l$$. Note that we calculate in units where $$\hbar\equiv 1$$ so that

$U(t, t_0) = \mathcal{T}\exp\left( -i\int_{t_0}^t\mathrm{d}t'\mathcal{H}(t') \right).$
Parameters
hamiltonian: array_like, shape (n_dt, d, d)

Hamiltonian of shape (n_dt, d, d) with d the dimensionality of the system

dt: array_like

The time differences

Returns
eigvals: ndarray

Array of eigenvalues of shape (n_dt, d)

eigvecs: ndarray

Array of eigenvectors of shape (n_dt, d, d)

propagators: ndarray

Array of cumulative propagators of shape (n_dt+1, d, d)

error_transfer_matrix(pulse: Optional[PulseSequence] = None, spectrum: Optional[numpy.ndarray] = None, omega: Optional[Sequence[float]] = None, n_oper_identifiers: Optional[Sequence[str]] = None, second_order: bool = False, cumulant_function: Optional[numpy.ndarray] = None, show_progressbar: bool = False, memory_parsimonious: bool = False, cache_intermediates: bool = False)numpy.ndarray[source]

Compute the error transfer matrix up to unitary rotations.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to compute the error transfer matrix.

spectrum: array_like, shape ([[n_nops,] n_nops,] n_omega)

The two-sided noise power spectral density in units of inverse frequencies as an array of shape (n_omega,), (n_nops, n_omega), or (n_nops, n_nops, n_omega). In the first case, the same spectrum is taken for all noise operators, in the second, it is assumed that there are no correlations between different noise sources and thus there is one spectrum for each noise operator. In the third and most general case, there may be a spectrum for each pair of noise operators corresponding to the correlations between them. n_nops is the number of noise operators considered and should be equal to len(n_oper_identifiers).

omega: array_like, shape (n_omega,)

The frequencies at which to calculate the filter functions.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which to evaluate the error transfer matrix. The default is all. Note that, since in general contributions from different noise operators won’t commute, not selecting all noise operators results in neglecting terms of order $$\xi^4$$.

second_order: bool, optional

Also take into account the frequency shifts $$\Delta$$ that correspond to second order Magnus expansion and constitute unitary terms. Default False.

cumulant_function: ndarray, shape ([[n_pls, n_pls,] n_nops,] n_nops, d**2, d**2)

A precomputed cumulant function. If given, pulse, spectrum, omega are not required.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

memory_parsimonious: bool, optional

Trade memory footprint for performance. See calculate_decay_amplitudes(). The default is False.

cache_intermediates: bool, optional

Keep and return intermediate terms of the calculation of the control matrix (if it is not already cached) that can be reused for second order or gradients. Can consume large amount of memory, but speed up the calculation.

Returns
error_transfer_matrix: ndarray, shape (d**2, d**2)

The error transfer matrix. The individual noise operator contributions are summed up before exponentiating as they might not commute.

calculate_cumulant_function

Calculate the cumulant function $$\mathcal{K}$$

calculate_decay_amplitudes

Calculate the $$\Gamma_{\alpha\beta,kl}$$

infidelity

Calculate only infidelity of a pulse.

Notes

The error transfer matrix is given by

$\tilde{\mathcal{U}} = \exp\mathcal{K}(\tau)$

with $$\mathcal{K}(\tau)$$ the cumulant function (see calculate_cumulant_function()). For Gaussian noise this expression is exact when taking into account the decay amplitudes $$\Gamma$$ and frequency shifts $$\Delta$$. As the latter effects coherent errors it can be neglected if we assume that the experimenter has calibrated their pulse.

For non-Gaussian noise the expression above is perturbative and includes noise up to order $$\xi^2$$ and hence $$\tilde{\mathcal{U}} = \mathbb{1} + \mathcal{K}(\tau) + \mathcal{O}(\xi^4)$$ (although it is evaluated as a matrix exponential in any case).

Given the above expression of the error transfer matrix, the entanglement fidelity is given by

$\mathcal{F}_\mathrm{e} = \frac{1}{d^2}\mathrm{tr}\,\tilde{\mathcal{U}}.$

## 2.5. filter_functions.plotting module¶

This module provides various plotting functions.

### 2.5.1. Functions¶

plot_bloch_vector_evolution()

Plot the evolution of the Bloch vector on a QuTiP-generated Bloch sphere

plot_filter_function()

Plot the filter function of a given PulseSequence

plot_infidelity_convergence()

Helper function called by infidelity() to plot the convergence of the infidelity

plot_pulse_correlation_filter_function()

Plot the pulse correlation filter function of a given PulseSequence

plot_pulse_train()

Plot the pulse train of a given PulseSequence

plot_cumulant_function()

Plot the cumulant function of a PulseSequence for a given spectrum as an image.

plot_bloch_vector_evolution(pulse: PulseSequence, psi0: Optional[Union[numpy.ndarray, qutip.qobj.Qobj]] = None, b: Optional[qutip.bloch.Bloch] = None, n_samples: Optional[int] = None, cmap: Optional[Union[matplotlib.colors.Colormap, str]] = None, show: bool = True, return_Bloch: bool = False, **bloch_kwargs)Union[None, qutip.bloch.Bloch][source]

Plot the evolution of the Bloch vector under the given pulse sequence.

Parameters
pulse: PulseSequence

The PulseSequence instance whose control Hamiltonian determines the time evolution of the Bloch vector.

psi0: Qobj or array_like, optional

The initial state before the pulse is applied. Defaults to $$|0\rangle$$.

b: qutip.Bloch, optional

If given, the QuTiP Bloch instance on which to plot the time evolution.

n_samples: int, optional

The number of time points to be sampled.

cmap: matplotlib colormap, optional

The colormap for the trajectory.

show: bool, optional**

Whether to show the sphere (by calling b.make_sphere()).

return_Bloch: bool, optional

Whether to return the qutip.bloch.Bloch instance

bloch_kwargs: dict, optional

A dictionary with keyword arguments to be fed into the qutip.Bloch constructor (if b not given).

Returns
b: qutip.Bloch

The qutip.Bloch instance

Raises
ValueError

If the pulse is for more than one qubit

qutip.bloch.Bloch

Qutip’s Bloch sphere implementation.

plot_cumulant_function(pulse: Optional[PulseSequence] = None, spectrum: Optional[numpy.ndarray] = None, omega: Optional[Sequence[float]] = None, cumulant_function: Optional[numpy.ndarray] = None, n_oper_identifiers: Optional[Sequence[int]] = None, second_order: bool = False, colorscale: str = 'linear', linthresh: Optional[float] = None, basis_labels: Optional[Sequence[str]] = None, basis_labelsize: Optional[int] = None, cbar_label: str = 'Cumulant Function', cbar_labelsize: Optional[int] = None, fig: Optional[matplotlib.figure.Figure] = None, grid: = None, cmap: Optional[Union[matplotlib.colors.Colormap, str]] = None, grid_kw: Optional[dict] = None, cbar_kw: Optional[dict] = None, imshow_kw: Optional[dict] = None, **figure_kw)Tuple[matplotlib.figure.Figure, mpl_toolkits.axes_grid1.axes_grid.ImageGrid][source]

Plot the cumulant function for a given noise spectrum as an image.

The cumulant function generates the error transfer matrix $$\tilde{\mathcal{U}}$$ exactly for Gaussian noise and to second order for non-Gaussian noise.

The function may be called with either a PulseSequence, a spectrum, and a list of frequencies in which case the cumulant function is calculated for those parameters, or with a precomputed cumulant function.

As of now, only auto-correlated spectra are implemented.

Parameters
pulse: ‘PulseSequence’

The pulse sequence.

spectrum: ndarray

The two-sided noise spectrum.

omega: array_like

The frequencies for which to evaluate the error transfer matrix.

cumulant_function: ndarray, shape (n_nops, d**2, d**2)

A precomputed cumulant function. If given, pulse, spectrum, omega are not required.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which the cumulant function should be plotted. All identifiers can be accessed via pulse.n_oper_identifiers. Defaults to all.

second_order: bool, optional

Also take into account the frequency shifts $$\Delta$$ that correspond to second order Magnus expansion and constitute unitary terms. Default False.

colorscale: str, optional

The scale of the color code (‘linear’ or ‘log’ (default))

linthresh: float, optional

The threshold below which the colorscale will be linear (only for ‘log’) colorscale

basis_labels: array_like (str), optional

Custom labels for the elements of the cumulant function (the basis elements). Grabbed from the basis by default.

basis_labelsize: int, optional

The size in points for the basis labels.

cbar_label: str, optional

The label for the colorbar. Default: ‘Cumulant Function’.

cbar_labelsize: int, optional

The size in points for the colorbar label.

fig: matplotlib figure, optional

A matplotlib figure instance to plot in

grid: matplotlib ImageGrid, optional

An ImageGrid instance to use for plotting.

cmap: matplotlib colormap, optional

The colormap for the matrix plot.

grid_kw: dict, optional

Dictionary with keyword arguments passed to the ImageGrid constructor.

cbar_kw: dict, optional

Dictionary with keyword arguments passed to the colorbar constructor.

imshow_kw: dict, optional

Dictionary with keyword arguments passed to imshow.

figure_kw: optional

Keyword argument dictionaries that are fed into the matplotlib.pyplot.figure() function if no fig instance is specified.

Returns
fig: matplotlib figure

The matplotlib figure instance used for plotting.

grid: matplotlib ImageGrid

The ImageGrid instance used for plotting.

plot_filter_function(pulse: PulseSequence, omega: Optional[Sequence[float]] = None, n_oper_identifiers: Optional[Sequence[int]] = None, fig: Optional[matplotlib.figure.Figure] = None, axes: Optional[matplotlib.axes._axes.Axes] = None, xscale: str = 'log', yscale: str = 'linear', omega_in_units_of_tau: bool = True, cycler: Optional[cycler.Cycler] = None, plot_kw: dict = {}, subplot_kw: Optional[dict] = None, gridspec_kw: Optional[dict] = None, **figure_kw)Tuple[matplotlib.figure.Figure, matplotlib.axes._axes.Axes, matplotlib.legend.Legend][source]

Plot the fidelity filter function(s) of the given PulseSequence for positive frequencies. As of now only the diagonal elements of $$F_{\alpha\beta}$$ are implemented, i.e. the filter functions corresponding to uncorrelated noise sources.

Parameters
pulse: PulseSequence

The pulse sequence whose filter function to plot.

omega: array_like, optional

The frequencies at which to evaluate the filter function. If not given, the pulse sequence’s omega attribute is used (if set) or sensible values are chosen automatically (if None)

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which the filter function should be plotted. All identifiers can be accessed via pulse.n_oper_identifiers. Defaults to all.

fig: matplotlib figure, optional

A matplotlib figure instance to plot in

axes: matplotlib axes, optional

A matplotlib axes instance to use for plotting.

xscale: str, optional

x-axis scaling. One of (‘linear’, ‘log’).

yscale: str, optional

y-axis scaling. One of (‘linear’, ‘log’).

omega_in_units_of_tau: bool, optional

Plot $$\omega\tau$$ or just $$\omega$$ on x-axis.

cycler: cycler.Cycler, optional

A Cycler instance used to set the style cycle if multiple lines are to be drawn

plot_kw: dict, optional

Dictionary with keyword arguments passed to the plot function

subplot_kw: dict, optional

Dictionary with keyword arguments passed to the subplots constructor

gridspec_kw: dict, optional

Dictionary with keyword arguments passed to the gridspec constructor

figure_kw: optional

Keyword argument dictionaries that are fed into the matplotlib.pyplot.subplots() function if no fig instance is specified.

Returns
fig: matplotlib figure

The matplotlib figure instance used for plotting.

axes: matplotlib axes

The matplotlib axes instance used for plotting.

legend: matplotlib legend

The matplotlib legend instance in the plot.

Raises
ValueError

If an invalid number of n_oper_labels were given

plot_infidelity_convergence(n_samples: Sequence[int], infids: Sequence[float], axes: Optional[matplotlib.axes._axes.Axes] = None)Tuple[matplotlib.figure.Figure, matplotlib.axes._axes.Axes][source]

Plot the convergence of the infidelity integral. The function arguments are those returned by infidelity() with the test_convergence flag set to True.

Parameters
n_samples: array_like

Array with the number of samples at which the integral was evaluated

infids: array_like, shape (n_samples, [n_oper_inds, optional])

Array with the calculated infidelities for each noise operator on the second axis or the second axis already traced out.

axes: sequence of two matplotlib axes, optional

Two axes that the result is plotted in.

Returns
fig: matplotlib figure

The matplotlib figure instance used for plotting.

axes: matplotlib axes

The matplotlib axes instances used for plotting.

plot_pulse_correlation_filter_function(pulse: PulseSequence, n_oper_identifiers: Optional[Sequence[int]] = None, fig: Optional[matplotlib.figure.Figure] = None, xscale: str = 'log', yscale: str = 'linear', omega_in_units_of_tau: bool = True, cycler: Optional[cycler.Cycler] = None, plot_kw: dict = {}, subplot_kw: Optional[dict] = None, gridspec_kw: Optional[dict] = None, **figure_kw)Tuple[matplotlib.figure.Figure, matplotlib.axes._axes.Axes, matplotlib.legend.Legend][source]

Plot the fidelity pulse correlation filter functions of the given PulseSequence if they were computed during concatenation for positive frequencies.

Returns a figure with n by n subplots where n is the number of pulses that were concatenated. As of now only the diagonal elements of $$F_{\alpha\beta}$$ are implemented, i.e. the filter functions corresponding to uncorrelated noise sources.

Parameters
pulse: PulseSequence

The pulse sequence whose filter function to plot.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which the filter function should be plotted. All identifiers can be accessed via pulse.n_oper_identifiers. Defaults to all.

fig: matplotlib figure, optional

A matplotlib figure instance to plot in

xscale: str, optional

x-axis scaling. One of (‘linear’, ‘log’).

yscale: str, optional

y-axis scaling. One of (‘linear’, ‘log’).

omega_in_units_of_tau: bool, optional

Plot $$\omega\tau$$ or just $$\omega$$ on x-axis.

cycler: cycler.Cycler, optional

A Cycler instance used to set the style cycle if multiple lines are to be drawn in one subplot. Used for all subplots.

plot_kw: dict, optional

Dictionary with keyword arguments passed to the plot function

subplot_kw: dict, optional

Dictionary with keyword arguments passed to the subplots constructor

gridspec_kw: dict, optional

Dictionary with keyword arguments passed to the gridspec constructor

figure_kw: optional

Keyword argument dictionaries that are fed into the matplotlib.pyplot.subplots() function if no fig instance if specified.

Returns
fig: matplotlib figure

The matplotlib figure instance used for plotting.

axes: matplotlib axes

The matplotlib axes instances used for plotting.

legend: matplotlib legend

The matplotlib legend instance in the plot.

Raises
CalculationError

If the pulse correlation filter function was not computed during concatenation.

plot_pulse_train(pulse: PulseSequence, c_oper_identifiers: Optional[Sequence[int]] = None, fig: Optional[matplotlib.figure.Figure] = None, axes: Optional[matplotlib.axes._axes.Axes] = None, cycler: Optional[cycler.Cycler] = None, plot_kw: Optional[dict] = {}, subplot_kw: Optional[dict] = None, gridspec_kw: Optional[dict] = None, **figure_kw)Tuple[matplotlib.figure.Figure, matplotlib.axes._axes.Axes, matplotlib.legend.Legend][source]

Plot the pulsetrain of the PulseSequence pulse.

Parameters
pulse: PulseSequence

The pulse sequence whose pulse train to plot.

c_oper_identifiers: array_like, optional

The identifiers of the control operators for which the pulse train should be plotted. All identifiers can be accessed via pulse.c_oper_identifiers. Defaults to all.

fig: matplotlib figure, optional

A matplotlib figure instance to plot in

axes: matplotlib axes, optional

A matplotlib axes instance to use for plotting.

cycler: cycler.Cycler, optional

A Cycler instance used to set the style cycle if multiple lines are to be drawn

plot_kw: dict, optional

Dictionary with keyword arguments passed to the plot function

subplot_kw: dict, optional

Dictionary with keyword arguments passed to the subplots constructor

gridspec_kw: dict, optional

Dictionary with keyword arguments passed to the gridspec constructor

figure_kw: optional

Keyword argument dictionaries that are fed into the matplotlib.pyplot.subplots() function if no fig instance is specified.

Returns
fig: matplotlib figure

The matplotlib figure instance used for plotting.

axes: matplotlib axes

The matplotlib axes instance used for plotting.

legend: matplotlib legend

The matplotlib legend instance in the plot.

Raises
ValueError

If an invalid number of c_oper_labels were given

## 2.6. filter_functions.pulse_sequence module¶

This module defines the PulseSequence class, the package’s central object.

### 2.6.1. Classes¶

PulseSequence

The pulse sequence defined by a Hamiltonian

### 2.6.2. Functions¶

concatenate()

Function to concatenate different PulseSequence instances and efficiently compute their joint filter function

concatenate_periodic()

Function to more efficiently concatenate many versions of the same PulseSequence instances and compute their joint filter function

extend()

Function to map several PulseSequence instances to different qubits, efficiently scaling up cached attributes.

class PulseSequence(*args, **kwargs)[source]

Bases: object

A class for pulse sequences and their filter functions.

The Hamiltonian is separated into a control and a noise part with

$\begin{split}\mathcal{H}_c &= \sum_i a_i(t) A_i \\ \mathcal{H}_n &= \sum_j s_j(t) b_j(t) B_j\end{split}$

where $$A_i$$ and $$B_j$$ are hermitian operators and $$b_j(t)$$ are classically fluctuating noise variables captured in a power spectral density and not needed at instantiation of a PulseSequence.

Parameters
H_c: list of lists

A nested list of n_cops nested lists as taken by QuTiP functions (see for example qutip.propagator.propagator()) describing the control part of the Hamiltonian. The i-th entry of the list should be a list consisting of the i-th operator $$A_i$$ making up the control Hamiltonian and a list or array $$a_i(t)$$ describing the magnitude of that operator during the time intervals dt. Optionally, the list may also include operator identifiers. That is, H_c should look something like this:

H = [[c_oper1, c_coeff1, c_oper_identifier1],
[c_oper2, c_coeff2, c_oper_identifier2], ...]


The operators may be given either as NumPy arrays or QuTiP Qobjs and each coefficient array should have the same number of elements as dt, and should be given in units of $$\hbar$$. If not every sublist (read: operator) was given a identifier, they are automatically filled up with ‘A_i’ where i is the position of the operator.

H_n: list of lists

A nested list of n_nops nested lists as taken by QuTiP functions (see for example qutip.propagator.propagator()) describing the noise part of the Hamiltonian. The j-th entry of the list should be a list consisting of the j-th operator $$B_j$$ making up the noise Hamiltonian and a list or array describing the sensitivity $$s_j(t)$$ of the system to the noise operator during the time intervals dt. Optionally, the list may also include operator identifiers. That is, H_n should look something like this:

H = [[n_oper1, n_coeff1, n_oper_identifier1],
[n_oper2, n_coeff2, n_oper_identifier2], ...]


The operators may be given either as NumPy arrays or QuTiP Qobjs and each coefficient array should have the same number of elements as dt, and should be given in units of $$\hbar$$. If not every sublist (read: operator) was given a identifier, they are automatically filled up with ‘B_i’ where i is the position of the operator.

dt: array_like, shape (n_dt,)

The segment durations of the Hamiltonian (i.e. durations of constant control). Internally, the control operation is taken to start at $$t_0\equiv 0$$, i.e. the edges of the constant control segments are at times t = [0, *np.cumsum(dt)].

basis: Basis, shape (d**2, d, d), optional

The operator basis in which to calculate. If a Generalized Gell-Mann basis (see ggm()) is chosen, some calculations will be faster for large dimensions due to a simpler basis expansion. However, when extending the pulse sequence to larger qubit registers, cached filter functions cannot be retained since the GGM basis does not factor into tensor products. In this case a Pauli basis is preferable.

Notes

Due to the heavy use of NumPy’s einsum() function, results have a floating point error of ~1e-13.

Examples

A rotation by $$\pi$$ around the axis between x and y preceeded and followed by a period of free evolution with the system subject to dephasing noise.

>>> import qutip as qt; import numpy as np
>>> H_c = [[qt.sigmax(), [0, np.pi, 0]],
[qt.sigmay(), [0, np.pi, 0]]]
>>> # Equivalent pulse:
>>> # H_c = [[qt.sigmax() + qt.sigmay(), [0, np.pi, 0]]]
>>> # The noise sensitivity is constant
>>> H_n = [[qt.sigmaz()/np.sqrt(2), [1, 1, 1], 'Z']]
>>> dt = [1, 1, 1]
>>> # Free evolution between t=0 and t=1, rotation between t=1 and t=2,
>>> # and free evolution again from t=2 to t=3.
>>> pulse = PulseSequence(H_c, H_n, dt)
>>> pulse.c_oper_identifiers
['A_0', 'A_1']
>>> pulse.n_oper_identifiers
['Z']
>>> omega = np.logspace(-1, 2, 500)
>>> F = pulse.get_filter_function(omega)    # shape (1, 500)
>>> # Plot the resulting filter function:
>>> from filter_functions import plotting
>>> fig, ax, leg = plotting.plot_filter_function(pulse)

Attributes
c_opers: ndarray, shape (n_cops, d, d)

Control operators

n_opers: ndarray, shape (n_nops, d, d)

Noise operators

c_oper_identifers: sequence of str

Identifiers for the control operators of the system

n_oper_identifers: sequence of str

Identifiers for the noise operators of the system

c_coeffs: ndarray, shape (n_cops, n_dt)

Control parameters in units of $$\hbar$$

n_coeffs: ndarray, shape (n_nops, n_dt)

Noise sensitivities in units of $$\hbar$$

dt: ndarray, shape (n_dt,)

Time steps

t: ndarray, shape (n_dt + 1,)

Absolute times taken to start at $$t_0\equiv 0$$

tau: float

Total duration. Equal to t[-1].

d: int

Dimension of the Hamiltonian

basis: Basis, shape (d**2, d, d)

The operator basis used for calculation

nbytes: int

An estimate of the memory consumed by the PulseSequence instance and its attributes

If the Hamiltonian was diagonalized, the eigenvalues and -vectors as
well as the cumulative propagators are cached:
eigvals: ndarray, shape (n_dt, d)

Eigenvalues $$D^{(g)}$$

eigvecs: ndarray, shape (n_dt, d, d)

Eigenvectors $$V^{(g)}$$

propagators: ndarray, shape (n_dt+1, d, d)

Cumulative propagators $$Q_g$$

total_propagator: ndarray, shape (d, d)

The total propagator $$Q$$ of the pulse alone. That is, $$|\psi(\tau)\rangle = propagators|\psi(0)\rangle$$.

total_propagator_liouville: array_like, shape (d**2, d**2)

The transfer matrix for the total propagator of the pulse. Given by liouville_representation(pulse.total_propagator, pulse.basis).

Furthermore, when the filter function is calculated, the frequencies
are cached as well as other relevant quantities.

Methods

 get_filter_function_derivative(omega[, …]) Calculate the pulse sequence’s filter function derivative.
 cleanup(method=’conservative’) Delete cached attributes is_cached(attr) Checks if a given attribute of the PulseSequence is cached diagonalize() Diagonalize the Hamiltonian of the pulse sequence, computing eigenvalues and -vectors as well as cumulative propagators get_control_matrix(omega, show_progressbar=False) Calculate the control matrix for frequencies omega get_filter_function(omega, which=’fidelity’, show_progressbar=False) Calculate the filter function for frequencies omega get_pulse_correlation_filter_function(which=’fidelity’) Get the pulse correlation filter function (only possible if computed during concatenation) propagator_at_arb_t(t) Calculate the propagator at arbitrary times

Initialize a PulseSequence instance.

cache_control_matrix(omega: Sequence[float], control_matrix: Optional[numpy.ndarray] = None, show_progressbar: bool = False, cache_intermediates: bool = False)None[source]

Cache the control matrix and the frequencies it was calculated for.

Parameters
omega: array_like, shape (n_omega,)

The frequencies for which to cache the filter function.

control_matrix: array_like, shape ([n_nops,] n_nops, d**2, n_omega), optional

The control matrix for the frequencies omega. If None, it is computed.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations. Only applies if control_matrix is not supplied.

cache_filter_function(omega: Sequence[float], control_matrix: Optional[numpy.ndarray] = None, filter_function: Optional[numpy.ndarray] = None, which: str = 'fidelity', order: int = 1, show_progressbar: bool = False, cache_intermediates: bool = False)None[source]

Cache the filter function. If control_matrix.ndim == 4, it is taken to be the ‘pulse correlation control matrix’ and summed along the first axis. In that case, also the pulse correlation filter function is calculated and cached. Total phase factors and transfer matrices of the the cumulative propagator are also cached so they can be reused during concatenation.

Parameters
omega: array_like, shape (n_omega,)

The frequencies for which to cache the filter function.

control_matrix: array_like, shape ([n_nops,] n_nops, d**2, n_omega), optional

The control matrix for the frequencies omega. If None, it is computed and the filter function derived from it.

filter_function: array_like, shape (n_nops, n_nops, [d**2, d**2,] n_omega), optional

The filter function for the frequencies omega. If None, it is computed from R in the case order == 1 and from scratch else.

which: str, optional

Which filter function to cache. Either ‘fidelity’ (default) or ‘generalized’.

order: int, optional

First or second order filter function.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations.

cache_total_phases(omega: Sequence[float], total_phases: Optional[numpy.ndarray] = None)None[source]

Cache the total phase factors for this pulse and omega.

Parameters
omega: array_like, shape (n_omega,)

The frequencies for which to cache the phase factors.

total_phases: array_like, shape (n_omega,), optional

The total phase factors for the frequencies omega. If None, they are computed.

cleanup(method: str = 'conservative')None[source]

Delete cached byproducts of the calculation of the filter function that are not necessarily needed anymore in order to free up memory.

Parameters
method: optional

If set to ‘conservative’ (the default), only the following attributes are deleted:

• _eigvals

• _eigvecs

• _propagators

If set to ‘greedy’, all of the above as well as the following attributes are deleted:

• _total_propagator

• _total_propagator_liouville

• _total_phases

• _control_matrix

• _control_matrix_pc

If set to ‘all’, all of the above as well as the following attributes are deleted:

• omega

• _filter_function

• _filter_function_gen

• _filter_function_pc

• _filter_function_pc_gen

• _filter_function_2

• _intermediates[‘control_matrix_step’]

If set to ‘frequency dependent’ only attributes that are functions of frequency are initalized to None.

Note that if this PulseSequence is concatenated with another one, some of the attributes might need to be calculated again, resulting in slower execution of the concatenation.

diagonalize()None[source]

Diagonalize the Hamiltonian defining the pulse sequence.

property duration: float

The duration of the pulse. Alias of tau.

property eigvals: numpy.ndarray

Get the eigenvalues of the pulse’s Hamiltonian.

property eigvecs: numpy.ndarray

Get the eigenvectors of the pulse’s Hamiltonian.

get_control_matrix(omega: Sequence[float], show_progressbar: bool = False, cache_intermediates: bool = False)numpy.ndarray[source]

Get the control matrix for the frequencies omega. If it has been cached for the same frequencies, the cached version is returned, otherwise it is calculated from scratch.

Parameters
omega: array_like, shape (n_omega,)

The frequencies at which to evaluate the control matrix.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations.

Returns
control_matrix: ndarray, shape (n_nops, d**2, n_omega)

The control matrix for the noise operators.

get_filter_function(omega: Sequence[float], which: str = 'fidelity', order: int = 1, show_progressbar: bool = False, cache_intermediates: bool = False)numpy.ndarray[source]

Get the first or second order filter function.

The filter function is cached so it doesn’t need to be calculated twice for the same frequencies.

Parameters
omega: array_like, shape (n_omega,)

The frequencies at which to evaluate the filter function.

which: str, optional

Which filter function to return. Either ‘fidelity’ (default) or ‘generalized’ (see Notes). Only if order == 1.

order: int, optional

First or second order filter function.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations.

Returns
filter_function: ndarray, shape (n_nops, n_nops, [d**2, d**2,] n_omega)

The filter function for each combination of noise operators as a function of omega.

Notes

The first-order generalized filter function is given by

$F_{\alpha\beta,kl}(\omega) = \tilde{\mathcal{B}}_{\alpha k}^\ast(\omega) \tilde{\mathcal{B}}_{\beta l}(\omega),$

where $$\alpha,\beta$$ are indices counting the noise operators $$B_\alpha$$ and $$k,l$$ indices counting the basis elements $$C_k$$.

The fidelity filter function is obtained by tracing over the basis indices:

$F_{\alpha\beta}(\omega) = \sum_{k} F_{\alpha\beta,kk}(\omega).$
get_filter_function_derivative(omega: Sequence[float], control_identifiers: Optional[Sequence[str]] = None, n_coeffs_deriv: Optional[Sequence[Sequence[float]]] = None)numpy.ndarray[source]

Calculate the pulse sequence’s filter function derivative.

Parameters
omega: array_like, shape (n_omega,)

Frequencies at which the pulse control matrix is to be evaluated.

control_identifiers: Sequence[str]

Sequence of strings with the control identifiern to distinguish between control and drift Hamiltonian. The default is None.

n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)

The derivatives of the noise susceptibilities by the control amplitudes. Defaults to None.

Returns
filter_function_deriv: ndarray, shape (n_nops, n_t, n_ctrl, n_omega)

The regular filter functions’ derivatives for variation in each control contribution.

get_pulse_correlation_control_matrix()numpy.ndarray[source]

Get the pulse correlation control matrix if it was cached.

get_pulse_correlation_filter_function(which: str = 'fidelity')numpy.ndarray[source]

Get the pulse correlation filter function given by

$F_{\alpha\beta}^{(gg')}(\omega) = e^{i\omega(t_{g-1} - t_{g'-1})} \tilde{\mathcal{B}}^{(g)}(\omega)\mathcal{Q}^{(g-1)} \mathcal{Q}^{(g'-1)\dagger} \tilde{\mathcal{B}}^{(g')\dagger}(\omega),$

where $$g,g'$$ index the pulse in the sequence and $$\alpha,\beta$$ index the noise operators, if it was computed during concatenation. Since the calculation requires the individual pulse’s control matrices and phase factors, which are not retained after concatenation, the pulse correlation filter function cannot be computed afterwards.

Note that the frequencies for which the filter function was calculated are not stored.

Returns
filter_function_pc: ndarray, shape (n_pls, n_pls, n_nops, n_nops, n_omega)

The pulse correlation filter function for each noise operator as a function of omega. The first two axes correspond to the pulses in the sequence, i.e. if the concatenated pulse sequence is $$C\circ B\circ A$$, the first two axes are arranged like

$\begin{split}F_{\alpha\beta}^{(gg')} &= \begin{pmatrix} F_{\alpha\beta}^{(AA)} & F_{\alpha\beta}^{(AB)} & F_{\alpha\beta}^{(AC)} \\ F_{\alpha\beta}^{(BA)} & F_{\alpha\beta}^{(BB)} & F_{\alpha\beta}^{(BC)} \\ F_{\alpha\beta}^{(CA)} & F_{\alpha\beta}^{(CB)} & F_{\alpha\beta}^{(CC)} \end{pmatrix}\end{split}$

for $$g,g'\in\{A, B, C\}$$.

get_total_phases(omega: Sequence[float])numpy.ndarray[source]

Get the (cached) total phase factors for this pulse and omega.

is_cached(attr: str)bool[source]

Returns True if the attribute is cached

property nbytes: int

Return an estimate of the amount of memory consumed by this object (or, more precisely, the array attributes of this object).

property omega: numpy.ndarray

Cached frequencies

propagator_at_arb_t(t: Sequence[float])numpy.ndarray[source]

Calculate the cumulative propagator Q(t) at times t by making use of the fact that we assume piecewise-constant control.

property propagators: numpy.ndarray

Get the eigenvectors of the pulse’s Hamiltonian.

property t: numpy.ndarray

The times of the pulse.

property tau: Union[float, int]

The duration of the pulse.

property total_propagator: numpy.ndarray

Get total propagator of the pulse.

property total_propagator_liouville: numpy.ndarray

Get the transfer matrix for the total propagator of the pulse.

concatenate(pulses: , calc_pulse_correlation_FF: bool = False, calc_filter_function: Optional[bool] = None, which: str = 'fidelity', omega: Optional[Sequence[float]] = None, show_progressbar: bool = False)filter_functions.pulse_sequence.PulseSequence[source]

Concatenate an arbitrary number of pulses. Note that pulses are concatenated left-to-right, that is,

$\mathtt{concatenate((A, B))} \equiv B \circ A$

so that $$A$$ is executed before $$B$$ when applying the concatenated pulse.

Parameters
pulses: sequence of PulseSequences

The PulseSequence instances to be concatenated. If any of the instances have a cached filter function, the filter function for the composite pulse will also be calculated in order to make use of the speedup gained from concatenating the filter functions. If omega is given, calculation of the composite filter function is forced.

calc_pulse_correlation_FF: bool, optional

Switch to control whether the pulse correlation filter function (see PulseSequence.get_pulse_correlation_filter_function()) is calculated. If omega is not given, the cached frequencies of all pulses need to be equal.

calc_filter_function: bool, optional

Switch to force the calculation of the filter function to be carried out or not. Overrides the automatic behavior of calculating it if at least one pulse has a cached control matrix. If True and no pulse has a cached control matrix, a list of frequencies must be supplied as omega.

which: str, optional

Which filter function to compute. Either ‘fidelity’ (default) or ‘generalized’ (see PulseSequence.get_filter_function() and PulseSequence.get_pulse_correlation_filter_function()).

omega: array_like, optional

Frequencies at which to evaluate the (pulse correlation) filter functions. If None, an attempt is made to use cached frequencies.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

Returns
pulse: PulseSequence

The concatenated pulse.

concatenate_periodic(pulse: filter_functions.pulse_sequence.PulseSequence, repeats: int)filter_functions.pulse_sequence.PulseSequence[source]

Concatenate a pulse sequence pulse whose Hamiltonian is periodic repeats times. Although performing the same task, this function is much faster for concatenating many identical pulses with filter functions than concatenate().

Note that for large dimensions, the calculation of the control matrix using this function might be very memory intensive.

Parameters
pulse: PulseSequence

The PulseSequence instance to be repeated. If it has a cached filter function, the filter function for the new pulse will also be computed.

repeats: int

The number of repetitions

Returns
newpulse: PulseSequence

The concatenated PulseSequence

concatenate

Concatenate arbitrary PulseSequences.

Notes

The total control matrix is given by

$\begin{split}\tilde{\mathcal{B}}(\omega) &= \tilde{\mathcal{B}}^{(1)}(\omega) \sum_{g=0}^{G-1} \left(e^{i\omega T}\right)^g \\ &= \tilde{\mathcal{B}}^{(1)}(\omega)\bigl( \mathbb{I} - e^{i\omega T} \mathcal{Q}^{(1)}\bigr)^{-1}\bigl( \mathbb{I} - \bigl(e^{i\omega T} \mathcal{Q}^{(1)}\bigr)^G\bigr).\end{split}$

with $$T$$ the period of the control Hamiltonian and $$G$$ the number of periods. The last equality is valid only if $$\mathbb{I} - e^{i\omega T}\mathcal{Q}^{(1)}$$ is invertible.

extend(pulse_to_qubit_mapping: Sequence[Sequence[Optional[Union[filter_functions.pulse_sequence.PulseSequence, Sequence[int], int, Mapping[str, str]]]]], N: Optional[int] = None, d_per_qubit: int = 2, additional_noise_Hamiltonian: Optional[Sequence[Sequence[Union[numpy.ndarray, qutip.qobj.Qobj, Sequence[float]]]]] = None, cache_diagonalization: Optional[bool] = None, cache_filter_function: Optional[bool] = None, omega: Optional[Sequence[float]] = None, show_progressbar: bool = False)filter_functions.pulse_sequence.PulseSequence[source]

Map one or more pulse sequences to different qubits.

Parameters
pulse_to_qubit_mapping: sequence of mapping tuples

A sequence of tuples with the first entry a PulseSequence instance and the second an int or tuple of ints indicating the qubits that the PulseSequence should be mapped to. A mapping of operator identifiers may optionally be given as a third element of each tuple. By default, the index of the qubit the operator is mapped to is appended to its identifier.

Pulse sequences defined for multiple qubits may also be extended to non-neighboring qubits. Note that for multi-qubit pulses the order of the qubits is respected, i.e. mapping a pulse to (1, 0) is different from mapping it to (0, 1).

N: int

The total number of qubits the new PulseSequence should be defined for. By default, this is inferred from pulse_to_qubit_mapping.

d_per_qubit: int

The size of the Hilbert space a single qubit requires.

Additional noise operators and corresponding sensitivities for the new pulse sequence.

cache_diagonalization: bool

Force diagonalizing the new pulse sequence. By default, diagonalization is cached if all pulses in pulse_to_qubit_mapping have been diagonalized since it is much cheaper to get the relevant quantities as tensor products from the mapped pulses instead of diagonalizing the new pulse.

cache_filter_function: bool

Force computing the filter functions for the new pulse sequence. Noise operators of individual pulses will be extended to the new Hilbert space. By default, this is done if all pulses in pulse_to_qubit_mapping have their filter functions cached.

Note that extending the filter functions is only possible if they the mapped pulses are using a separable basis like the Pauli basis.

omega: array_like

Frequencies for which to compute the filter functions if cache_filter_function == True. Defaults to None, in which case the cached frequencies of the individual pulses need to be the same.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

Returns
newpulse: PulseSequence

The new pulse sequence on the larger qubit register. The noise operators (and possibly filter functions) are stored in the following order: first those of the multi-qubit pulses in the order they appeared in pulse_to_qubit_mapping, then those of the single-qubit pulses, and lastly any additional ones that may be given by additional_noise_Hamiltonian.

remap

Map PulseSequence to a different qubit.

concatenate

Concatenate PulseSequences (in time).

concatenate_periodic

Periodically concatenate a PulseSequence.

Examples

>>> import filter_functions as ff
>>> I, X, Y, Z = ff.util.paulis
>>> X_pulse = ff.PulseSequence([[X, [np.pi/2], 'X']],
...                            [[X, [1], 'X'], [Z, [1], 'Z']],
...                            [1], basis=ff.Basis.pauli(1))
>>> XX_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)])
>>> XX_pulse.d
4
>>> XIX_pulse_1 = ff.extend([(X_pulse, 0), (X_pulse, 2)])
>>> XIX_pulse_1.d
8
>>> XXI_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)], N=3)
>>> XXI_pulse.d
8


Filter functions are automatically cached if they are for mapped pulses:

>>> omega = ff.util.get_sample_frequencies(X_pulse)
>>> X_pulse.cache_filter_function(omega)
>>> XX_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)])
>>> XX_pulse.is_cached('filter_function')
True


This behavior can also be overriden manually:

>>> XX_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)],
...                      cache_filter_function=False)
>>> XX_pulse.is_cached('filter_function')
False


Mapping pulses to non-neighboring qubits is also possible:

>>> Y_pulse = ff.PulseSequence([[Y, [np.pi/2], 'Y']],
...                            [[Y, [1], 'Y'], [Z, [1], 'Z']],
...                            [1], basis=ff.Basis.pauli(1))
>>> XXY_pulse = ff.extend([(XX_pulse, (0, 1)), (Y_pulse, 2)])
>>> XYX_pulse = ff.extend([(XX_pulse, (0, 2)), (Y_pulse, 1)])


Additionally, pulses can have the order of the qubits they are defined for permuted (see remap()):

>>> Z_pulse = ff.PulseSequence([[Z, [np.pi/2], 'Z']], [[Z, [1], 'Z']],
...                            [1], basis=ff.Basis.pauli(1))
>>> XY_pulse = ff.extend([(X_pulse, 0), (Y_pulse, 1)])
>>> YZX_pulse = ff.extend([(XY_pulse, (2, 0)), (Z_pulse, 1)])


Control and noise operator identifiers can be mapped according to a specified mapping:

>>> YX_pulse = ff.extend([(X_pulse, 1, {'X': 'IX', 'Z': 'IZ'}),
...                       (Y_pulse, 0, {'Y': 'YI', 'Z': 'ZI'})])
>>> YX_pulse.c_oper_identifiers
array(['IX', 'YI'], dtype='<U2')
>>> YX_pulse.n_oper_identifiers
array(['IX', 'IZ', 'YI', 'ZI'], dtype='<U2')


>>> H_n = [[ff.util.tensor(Z, Z, Z), [1], 'ZZZ']]
>>> XYX_pulse = ff.extend([(XX_pulse, (0, 2)), (Y_pulse, 1)],
>>> 'ZZZ' in XYX_pulse.n_oper_identifiers
True

remap(pulse: filter_functions.pulse_sequence.PulseSequence, order: Sequence[int], d_per_qubit: int = 2, oper_identifier_mapping: Optional[Mapping[str, str]] = None)filter_functions.pulse_sequence.PulseSequence[source]

Remap a PulseSequence by changing the order of qubits in the register. Cached attributes are automatically attempted to be retained.

Caution

This function simply permutes the order of the tensor product elements of control and noise operators. Thus, the resultant pulse will have its filter functions defined for different noise operators than the original one.

Parameters
pulse: PulseSequence

The pulse whose qubit order should be permuted.

order: sequence of ints

A list of permutation indices. E.g., if pulse is defined for two qubits, order == [1, 0] will reverse the order of qubits.

d_per_qubit: int (default: 2)

The size of the Hilbert space a single qubit inhabitates.

oper_identifier_mapping: dict_like

A mapping that maps operator identifiers from the old pulse to the remapped pulse. The default is the identity mapping.

Returns
remapped_pulse: PulseSequence

A new PulseSequence instance with the order of the qubits permuted according to order.

extend

Map PulseSequences to composite Hilbert spaces.

util.tensor_transpose

Transpose the order of a tensor product.

Examples

>>> X, Y = util.paulis[1:3]
>>> XY, YX = util.tensor(X, Y), util.tensor(Y, X)
>>> pulse = PulseSequence([[XY, [np.pi/2], 'XY']], [[YX, [1], 'YX']], [1],
...                       Basis.pauli(2))
>>> mapping = {'XY': 'YX', 'YX': 'XY'}
>>> remapped_pulse = remap(pulse, (1, 0), oper_identifier_mapping=mapping)
>>> target_pulse = PulseSequence([[YX, [np.pi/2], 'YX']],
...                              [[XY, [1], 'XY']], [1], Basis.pauli(2))
>>> remapped_pulse == target_pulse
True


Caching of attributes is automatically handled >>> remapped_pulse.is_cached(‘filter_function’) False >>> pulse.cache_filter_function(util.get_sample_frequencies(pulse)) >>> remapped_pulse = remap(pulse, (1, 0)) >>> remapped_pulse.is_cached(‘filter_function’) True

This module implements functions to calculate filter function and infidelity derivatives. Currently only auto-correlated noise (i.e. no cross-correlations) is implemented.

Throughout this documentation the following notation will be used:
• n_dt denotes the number of time steps,

• n_cops the number of all control operators,

• n_ctrl the number of accessible control operators (if identifiers are provided, otherwise n_ctrl=n_cops),

• n_nops the number of noise operators,

• n_omega the number of frequency samples, and

• d the dimension of the system.

### 2.7.1. Functions¶

calculate_derivative_of_control_matrix_from_scratch()

Calculate the derivative of the control matrix from scratch.

calculate_canonical_filter_function_derivative()

Compute the filter function derivative from the control matrix.

infidelity_derivative()

Calculate the infidelity derivative.

calculate_derivative_of_control_matrix_from_scratch(omega: Sequence[float], propagators: numpy.ndarray, eigvals: numpy.ndarray, eigvecs: numpy.ndarray, basis: filter_functions.basis.Basis, t: Sequence[float], dt: Sequence[float], n_opers: Sequence[Union[numpy.ndarray, qutip.qobj.Qobj]], n_coeffs: Sequence[Sequence[float]], c_opers: Sequence[Union[numpy.ndarray, qutip.qobj.Qobj]], all_identifiers: Sequence[str], control_identifiers: Optional[Sequence[str]] = None, n_coeffs_deriv: Optional[Sequence[Sequence[float]]] = None, intermediates: Optional[Dict[str, numpy.ndarray]] = None)numpy.ndarray[source]

Calculate the derivative of the control matrix from scratch.

Parameters
omega: array_like, shape (n_omega,)

Frequencies, at which the pulse control matrix is to be evaluated.

propagators: array_like, shape (n_dt+1, d, d)

The propagators $$Q_g = P_g P_{g-1}\cdots P_0$$ as a (d, d) array with d the dimension of the Hilbert space.

eigvals: array_like, shape (n_dt, d)

Eigenvalue vectors for each time pulse segment g with the first axis counting the pulse segment, i.e. HD == array([D_0, D_1, ...]).

eigvecs: array_like, shape (n_dt, d, d)

Eigenvector matrices for each time pulse segment g with the first axis counting the pulse segment, i.e. HV == array([V_0, V_1, ...]).

basis: Basis, shape (d**2, d, d)

The basis elements, in which the pulse control matrix will be expanded.

t: array_like, shape (n_dt+1), optional

The absolute times of the different segments.

dt: array_like, shape (n_dt)

Sequence duration, i.e. for the $$g$$-th pulse $$t_g - t_{g-1}$$.

n_opers: array_like, shape (n_nops, d, d)

Noise operators $$B_\alpha$$.

n_coeffs: array_like, shape (n_nops, n_dt)

The sensitivities of the system to the noise operators given by n_opers at the given time step.

c_opers: array_like, shape (n_cops, d, d)

Control operators $$H_k$$.

all_identifiers: array_like, shape (n_cops)

Identifiers of all control operators.

control_identifiers: Sequence[str], shape (n_ctrl), Optional

Sequence of strings with the control identifiers to distinguish between accessible control and drift Hamiltonian. The default is None.

n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)

The derivatives of the noise susceptibilities by the control amplitudes. Defaults to None.

intermediates: Dict[str, ndarray], optional

Optional dictionary containing intermediate results of the calculation of the control matrix.

Returns
ctrlmat_deriv: ndarray, shape (n_ctrl, n_omega, n_dt, n_nops, d**2)

Partial derivatives of the total control matrix with respect to each control direction $$\frac{\partial\mathcal{B}_{\alpha k}(\omega)}{\partial u_h(t_{g'})}$$.

Raises
ValueError

If the given identifiers control_identifier are not subset of the total identifiers all_identifiers of all control operators.

_liouville_derivative
_control_matrix_at_timestep_derivative

Notes

The derivative of the control matrix is calculated according to

$\frac{\partial\mathcal{B}_{\alpha k}(\omega)}{\partial u_h(t_{g'})} = \sum_{g=1}^G \mathrm{e}^{\mathrm{i}\omega t_{g-1}}\cdot\left(\sum_j \left[\frac{\partial\mathcal{B}_{\alpha j}^{(g)}(\omega)} {\partial u_h(t_{g'})} \cdot \mathcal{Q}_{jk}^{(g-1)} + \mathcal{B}_{\alpha j}^{(g)}(\omega) \cdot\frac{\partial \mathcal{Q}_{jk}^{(g-1)}} {\partial u_h(t_{g'})} \right] \right)$
calculate_filter_function_derivative(ctrlmat: numpy.ndarray, ctrlmat_deriv: numpy.ndarray)numpy.ndarray[source]

Compute the filter function derivative from the control matrix.

Parameters
ctrlmat: array_like, shape (n_nops, d**2, n_omega)

The control matrix.

ctrlmat_deriv: array_like, shape (n_nops, d**2, n_t, n_ctrl, n_omega)

The derivative of the control matrix.

Returns
filter_function_deriv: ndarray, shape (n_nops, n_dt, n_ctrl, n_omega)

The regular filter functions’ derivatives for variation in each control contribution $$\frac{\partial F_\alpha(\omega)}{\partial u_h(t_{g'})}$$.

Notes

The filter function derivative is calculated according to

$\frac{\partial F_\alpha(\omega)}{\partial u_h(t_{g'})} = 2 \mathrm{Re}\left(\sum_k \mathcal{B}_{\alpha k}^\ast(\omega) \frac{\partial\mathcal{B}_{\alpha k}(\omega)} {\partial u_h(t_{g'})} \right)$
infidelity_derivative(pulse: PulseSequence, spectrum: Sequence[float], omega: Sequence[float], control_identifiers: Optional[Sequence[str]] = None, n_coeffs_deriv: Optional[Sequence[Sequence[float]]] = None)numpy.ndarray[source]

Calculate the entanglement infidelity derivative of the PulseSequence pulse.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to calculate the infidelity.

spectrum: array_like, shape ([[n_nops,] n_nops,] omega)

The two-sided noise power spectral density in units of inverse frequencies as an array of shape (n_omega,), (n_nops, n_omega), or (n_nops, n_nops, n_omega). In the first case, the same spectrum is taken for all noise operators, in the second, it is assumed that there are no correlations between different noise sources and thus there is one spectrum for each noise operator. In the third and most general case, there may be a spectrum for each pair of noise operators corresponding to the correlations between them. n_nops is the number of noise operators considered and should be equal to len(n_oper_identifiers).

omega: array_like, shape (n_omega,)

The frequencies at which the integration is to be carried out.

control_identifiers: Sequence[str], shape (n_ctrl,)

Sequence of strings with the control identifiern to distinguish between accessible control and drift Hamiltonian.

n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)

The derivatives of the noise susceptibilities by the control amplitudes. Defaults to None.

Returns
infid_deriv: ndarray, shape (n_nops, n_dt, n_ctrl)

Array with the derivative of the infidelity for each noise source taken for each control direction at each time step $$\frac{\partial I_e}{\partial u_h(t_{g'})}$$.

Raises
ValueError

If the provided noise spectral density does not fit expected shape.

Notes

The infidelity’s derivative is given by

$\frac{\partial I_e}{\partial u_h(t_{g'})} = \frac{1}{d} \int_{-\infty}^\infty \frac{d\omega}{2\pi} S_\alpha(\omega) \frac{\partial F_\alpha(\omega)} {\partial u_h(t_{g'})}$

with $$S_{\alpha}(\omega)$$ the noise spectral density and $$F_{\alpha}(\omega)$$ the canonical filter function for noise source $$\alpha$$.

To convert to the average gate infidelity, use the following relation given by Horodecki et al. [Hor99] and Nielsen [Nie02]:

$\big\langle\mathcal{I}_\mathrm{avg}\big\rangle = \frac{d}{d+1} \big\langle\mathcal{I}_\mathrm{e}\big\rangle.$

References

Hor99

Horodecki, M., Horodecki, P., & Horodecki, R. (1999). General teleportation channel, singlet fraction, and quasidistillation. Physical Review A - Atomic, Molecular, and Optical Physics, 60(3), 1888–1898. https://doi.org/10.1103/PhysRevA.60.1888

Nie02

Nielsen, M. A. (2002). A simple formula for the average gate fidelity of a quantum dynamical operation. Physics Letters, Section A: General, Atomic and Solid State Physics, 303(4), 249–252. https://doi.org/10.1016/S0375-9601(02)01272-0

## 2.8. filter_functions.superoperator module¶

This module provides some functions related to superoperators and quantum maps.

### 2.8.1. Functions¶

liouville_representation()

Calculate the Liouville representation of a unitary with respect to a basis

liouville_to_choi()

Convert from Liouville to Choi matrix representation.

liouville_is_CP()

Check if superoperator in Liouville representation is completely positive.

liouville_is_cCP()

Check if superoperator in Liouville representation is conditional CP.

liouville_is_CP(superoperator: numpy.ndarray, basis: filter_functions.basis.Basis, return_eig: Optional[bool] = False, atol: Optional[float] = None)Union[bool, Tuple[bool, Tuple[numpy.ndarray, numpy.ndarray]]][source]

Test if a Liouville superoperator is completely positive (CP).

Parameters
superoperator: ndarray, shape (…, d**2, d**2)

The superoperator in Liouville representation to be checked for CPness.

basis: Basis, shape (d**2, d, d)

The operator basis defining the Liouville representation.

return_eig: bool, optional

Return the tuple of eigenvalues and eigenvectors of the Choi matrix. The default is False.

atol: float, optional

Absolute tolerance for the complete positivity.

Returns
CP: bool, (shape (…,))

The (array, if broadcasted) of bools indicating if superoperator is CP.

(D, V): Tuple[ndarray, ndarray]

The eigenvalues and eigenvectors of the Choi matrix (only if return_eig is True).

liouville_representation

Calculate Liouville representation of a unitary.

Liouville_to_choi

Convert from Liouville to Choi matrix representation.

liouville_is_cCP

Test if a superoperator is conditional CP.

Notes

A superoperator $$\mathcal{S}$$ is completely positive (CP) if and only if its Choi matrix representation is positive semidefinite:

$\mathcal{S}\text{ is CP }\Leftrightarrow \mathrm{choi}(\mathcal{S})\geq 0.$
liouville_is_cCP(superoperator: numpy.ndarray, basis: filter_functions.basis.Basis, return_eig: Optional[bool] = False, atol: Optional[float] = None)Union[bool, Tuple[bool, Tuple[numpy.ndarray, numpy.ndarray]]][source]

Test if a Liouville superoperator is conditional completely positive.

Parameters
superoperator: ndarray, shape (…, d**2, d**2)

The superoperator in Liouville representation to be checked for cCPness

basis: Basis, shape (d**2, d, d)

The operator basis defining the Liouville representation.

return_eig: bool, optional

Return the tuple of eigenvalues and eigenvectors of the Choi matrix projected on the complement of the maximally entangled state. The default is False.

atol: float, optional

Absolute tolerance for the complete positivity.

Returns
cCP: bool, (shape (…,))

The (array, if broadcasted) of bools indicating if superoperator is cCP

(D, V): Tuple[ndarray, ndarray]

The eigenvalues and eigenvectors of the projected Choi matrix (only if return_eig is True).

liouville_representation

Calculate Liouville representation of a unitary.

Liouville_to_choi

Convert from Liouville to Choi matrix representation.

liouville_is_CP

Test if a superoperator is CP.

Notes

A superoperator $$\mathcal{S}$$ is conditional completely positive (cCP) if and only if its Choi matrix projected on the complement of the maximally entangled state is positive semidefinite:

$\mathcal{S}\text{ is cCP }\Leftrightarrow Q\mathrm{choi}(\mathcal{S})Q\geq 0$

with $$Q = \mathbb{I} - |\Omega\rangle\langle\Omega|$$.

liouville_representation(U: numpy.ndarray, basis: filter_functions.basis.Basis)numpy.ndarray[source]

Get the Liouville representaion of the unitary U with respect to the basis.

Parameters
U: ndarray, shape (…, d, d)

The unitary.

basis: Basis, shape (d**2, d, d)

The basis used for the representation, e.g. a Pauli basis.

Returns
R: ndarray, shape (…, d**2, d**2)

The Liouville representation of U.

Notes

The Liouville representation of a unitary quantum operation $$\mathcal{U}:\rho\rightarrow U\rho U^\dagger$$ is given by

$\mathcal{U}_{ij} = \mathrm{tr}(C_i U C_j U^\dagger)$

with $$C_i$$ elements of the basis spanning $$\mathbb{C}^{d\times d}$$ with $$d$$ the dimension of the Hilbert space.

liouville_to_choi(superoperator: numpy.ndarray, basis: filter_functions.basis.Basis)numpy.ndarray[source]

Convert from Liouville to Choi matrix representation.

Parameters
superoperator: ndarray, shape (…, d**2, d**2)

The Liouville representation of a superoperator.

basis: Basis, shape (d**2, d, d)

The operator basis defining the Liouville representation.

Returns
choi: ndarray, shape (…, d**2, d**2)

The Choi matrix representation of the superoperator.

liouville_representation

Calculate Liouville representation of a unitary.

liouville_is_CP

Test if a superoperator is completely positive (CP).

liouville_is_cCP

Test if a superoperator is conditional CP.

Notes

The Choi matrix is given by

$\begin{split}\mathrm{choi}(\mathcal{S}) &= (\mathbb{I}\otimes\mathcal{S}) (|\Omega\rangle\langle\Omega|) \\ &= \sum_{ij} E_{ij}\otimes\mathcal{S}(E_{ij}) \\ &= \sum_{ij}\mathcal{S}_{ij} C_j^T\otimes C_i\end{split}$

where $$|\Omega\rangle$$ is a maximally entangled state, $$E_{ij} = |i\rangle\langle j|$$, and $$C_i$$ are the basis elements that define the Liouville representation $$\mathcal{S}_{ij}$$ [Mer13].

References

Mer13

Merkel, S. T. et al. Self-consistent quantum process tomography. Physical Review A - Atomic, Molecular, and Optical Physics, 87, 062119 (2013). https://doi.org/10.1103/PhysRevA.87.062119

## 2.9. filter_functions.types module¶

Defines custom types for the package.

## 2.10. filter_functions.util module¶

This module provides various helper functions.

### 2.10.1. Functions¶

abs2()

Absolute value squared

get_indices_from_identifiers()

The the indices of a subset of identifiers within a list of identifiers.

tensor()

Fast, flexible tensor product of an arbitrary number of inputs using einsum()

tensor_insert()

For an array that is known to be a tensor product, insert arrays at a given position in the product chain

tensor_merge()

For two arrays that are tensor products of known dimensions, merge them at arbitary positions in the product chain

tensor_transpose()

For a tensor product, transpose the order of the constituents in the product chain

mdot()

Multiple matrix product

remove_float_errors()

Set entries whose absolute value is below a certain threshold to zero

oper_equiv()

Determine if two vectors or operators are equal up to a global phase

dot_HS()

Hilbert-Schmidt inner product

get_sample_frequencies()

Get frequencies with typical infrared and ultraviolet cutoffs for a PulseSequence

progressbar()

A progress bar for loops. Uses tqdm if available and a simple custom one if not.

hash_array_along_axis()

Return a list of hashes along a given axis

all_array_equal()

Check if all arrays in an iterable are equal

### 2.10.2. Exceptions¶

CalculationError

Exception raised if trying to fetch the pulse correlation function when it was not computed during concatenation

abs2(x: numpy.ndarray)numpy.ndarray[source]

Fast function to calculate the absolute value squared,

$|\cdot|^2 := \Re(\cdot)^2 + \Im(\cdot)^2$

Equivalent to:

np.abs(x)**2

all_array_equal(it: Iterable)bool[source]

Return True if all array elements of it are equal by hashing the bytes representation of each array. Note that this is not thread-proof.

dot_HS(U: Union[numpy.ndarray, qutip.qobj.Qobj], V: Union[numpy.ndarray, qutip.qobj.Qobj], eps: Optional[float] = None)float[source]

Return the Hilbert-Schmidt inner product of U and V,

$\langle U, V\rangle_\mathrm{HS} := \mathrm{tr}(U^\dagger V).$
Parameters
U, V: qutip.Qobj or ndarray

Objects to compute the inner product of.

eps: float

The floating point precision. The result is rounded to abs(int(np.log10(eps))) decimals if eps > 0.

Returns
result: float, complex

The result rounded to precision eps.

Examples

>>> U, V = paulis[1:3]
>>> dot_HS(U, V)
0.0
>>> dot_HS(U, U)
2.0

get_sample_frequencies(pulse: PulseSequence, n_samples: int = 300, spacing: str = 'log', include_quasistatic: bool = False)numpy.ndarray[source]

Get n_samples sample frequencies spaced ‘linear’ or ‘log’.

The ultraviolet cutoff is taken to be two orders of magnitude larger than the timescale of the pulse tau. In the case of log spacing, the values are clipped in the infrared at two orders of magnitude below the timescale of the pulse.

Parameters
pulse: PulseSequence

The pulse to get frequencies for.

n_samples: int, optional

The number of frequency samples. Default is 300.

spacing: str, optional

The spacing of the frequencies. Either ‘log’ or ‘linear’, default is ‘log’.

include_quasistatic: bool, optional

Include zero frequency. Default is False.

Returns
omega: ndarray

The frequencies.

hash_array_along_axis(arr: numpy.ndarray, axis: int = 0)List[int][source]

Return the hashes of arr along the first axis

mdot(arr: Sequence, axis: int = 0)numpy.ndarray[source]

Multiple matrix products along axis

oper_equiv(psi: Union[numpy.ndarray, qutip.qobj.Qobj], phi: Union[numpy.ndarray, qutip.qobj.Qobj], eps: Optional[float] = None, normalized: bool = False)Tuple[bool, float][source]

Checks whether psi and phi are equal up to a global phase, i.e.

$|\psi\rangle = e^{i\chi}|\phi\rangle \Leftrightarrow \langle \phi|\psi\rangle = e^{i\chi},$

and returns the phase. If the first return value is false, the second is meaningless in this context. psi and phi can also be operators.

Parameters
psi, phi: qutip.Qobj or array_like

Vectors or operators to be compared

eps: float

The tolerance below which the two objects are treated as equal, i.e., the function returns True if abs(1 - modulus) <= eps.

normalized: bool

Flag indicating if psi and phi are normalized with respect to the Hilbert-Schmidt inner product dot_HS().

Examples

>>> psi = paulis[1]
>>> phi = paulis[1]*np.exp(1j*1.2345)
>>> oper_equiv(psi, phi)
(True, 1.2345)

progressbar(iterable: Iterable, *args, **kwargs)[source]

Progress bar for loops. Uses tqdm.

Usage:

for i in progressbar(range(10)):
do_something()

remove_float_errors(arr: numpy.ndarray, eps_scale: Optional[float] = None)numpy.ndarray[source]

Clean up arr by removing floating point numbers smaller than the dtype’s precision multiplied by eps_scale. Treats real and imaginary parts separately.

Obviously only works for arrays with norm ~1.

tensor(*args, rank: int = 2, optimize: Union[bool, str] = False)numpy.ndarray[source]

Fast, flexible tensor product using einsum. The product is taken over the last rank axes and broadcast over the remaining axes which thus need to follow numpy broadcasting rules. Note that vectors are treated as rank 2 tensors with shape (1, x) or (x, 1).

For example, the following shapes are compatible:

• rank == 2 (e.g. matrices or vectors):

(a, b, c, d, d), (a, b, c, e, e) -> (a, b, c, d*e, d*e)
(a, b, c), (a, d, e) -> (a, b*d, c*e)
(a, b), (c, d, e) -> (c, a*d, b*e)
(1, a), (b, 1, c) -> (b, 1, a*c)

• rank == 1:

(a, b), (a, c) -> (a, b*c)
(a, b, 1), (a, c) -> (a, b, c)

Parameters
args: array_like

The elements of the tensor product

rank: int, optional (default: 2)

The rank of the tensors. E.g., for a Kronecker product between two matrices rank == 2. The remaining axes are broadcast over.

optimize: bool|str, optional (default: False)

Optimize the tensor contraction order. Passed through to numpy.einsum().

numpy.kron

NumPy tensor product.

tensor_insert

Insert array at given position in tensor product chain.

tensor_merge

Merge tensor product chains.

tensor_transpose

Transpose the order of a tensor product chain.

Examples

>>> Z = np.diag([1, -1])
>>> np.array_equal(tensor(Z, Z), np.kron(Z, Z))
True

>>> A, B = np.arange(2), np.arange(2, 5)
>>> tensor(A, B, rank=1)
array([[0, 0, 0, 2, 3, 4]])

>>> args = np.random.randn(4, 10, 3, 2)
>>> result = tensor(*args, rank=1)
>>> result.shape == (10, 3, 2**4)
True
>>> result = tensor(*args, rank=2)
>>> result.shape == (10, 3**4, 2**4)
True

>>> A, B = np.random.randn(1, 3), np.random.randn(3, 4)
>>> result = tensor(A, B)
>>> result.shape == (1*3, 3*4)
True

>>> A, B = np.random.randn(3, 1, 2), np.random.randn(2, 2, 2)
>>> try:
...     result = tensor(A, B, rank=2)
... except ValueError as err:  # cannot broadcast over axis 0
...     print(err)
Incompatible shapes (3, 1, 2) and (2, 2, 2) for tensor product of rank 2.
>>> result = tensor(A, B, rank=3)
>>> result.shape == (3*2, 1*2, 2*2)
True

tensor_insert(arr: numpy.ndarray, *args, pos: Union[int, Sequence[int]], arr_dims: Sequence[Sequence[int]], rank: int = 2, optimize: Union[bool, str] = False)numpy.ndarray[source]

For a tensor product arr, insert args into the product chain at pos. E.g, if $$\verb|arr|\equiv A\otimes B\otimes C$$ and $$\verb|pos|\equiv 2$$, the result will be the tensor product

$A\otimes B\otimes\left[\bigotimes_{X\in\verb|args|}X\right] \otimes C.$

This function works in a similar way to numpy.insert() and the following would be functionally equivalent in the case that the constituent tensors of the product arr are known:

>>> tensor_insert(tensor(*arrs, rank=rank), *args, pos=pos, arr_dims=...,
...               rank=rank)

>>> tensor(*np.insert(arrs, pos, args, axis=0), rank=rank)

Parameters
arr: ndarray

The tensor product in whose chain the other args should be inserted

*args: ndarray

The tensors to be inserted in the product chain

pos: int|sequence of ints

The position(s) at which the args are inserted in the product chain. If an int and len(args) > 1, it is repeated so that all args are inserted in a row. If a sequence, it should indicate the indices in the original tensor product chain that led to arr before which args should be inserted.

arr_dims: array_like, shape (rank, n_const)

The last rank dimensions of the n_const constituent tensors of the tensor product arr as a list of lists with the list at position i containing the i-th relevant dimension of all args. Since the remaing axes are broadcast over, their shape is irrelevant.

For example, if arr = tensor(a, b, c, rank=2) and a,b,c have shapes (2, 3, 4), (5, 2, 2, 1), (2, 2), arr_dims = [[3, 2, 2], [4, 1, 2]].

rank: int, optional (default: 2)

The rank of the tensors. E.g., for a Kronecker product between two vectors, rank == 1, and between two matrices rank == 2. The remaining axes are broadcast over.

optimize: bool|str, optional (default: False)

Optimize the tensor contraction order. Passed through to numpy.einsum().

numpy.insert

NumPy array insertion with similar syntax.

numpy.kron

NumPy tensor product.

tensor_insert

Insert array at given position in tensor product chain.

tensor_merge

Merge tensor product chains.

tensor_transpose

Transpose the order of a tensor product chain.

Examples

>>> I, X, Y, Z = paulis
>>> arr = tensor(X, I)
>>> r = tensor_insert(arr, Y, Z, arr_dims=[[2, 2], [2, 2]], pos=0)
>>> np.allclose(r, tensor(Y, Z, X, I))
True
>>> r = tensor_insert(arr, Y, Z, arr_dims=[[2, 2], [2, 2]], pos=1)
>>> np.allclose(r, tensor(X, Y, Z, I))
True
>>> r = tensor_insert(arr, Y, Z, arr_dims=[[2, 2], [2, 2]], pos=2)
>>> np.allclose(r, tensor(X, I, Y, Z))
True


Other ranks and different dimensions:

>>> from numpy.random import randn
>>> A, B, C = randn(2, 3, 1, 2), randn(2, 2, 2, 2), randn(3, 2, 1)
>>> arr = tensor(A, C, rank=3)
>>> r = tensor_insert(arr, B, pos=1, rank=3,
...                   arr_dims=[[3, 3], [1, 2], [2, 1]])
>>> np.allclose(r, tensor(A, B, C, rank=3))
True

>>> arrs, args = randn(2, 2, 2), randn(2, 2, 2)
>>> arr_dims = [[2, 2], [2, 2]]
>>> r = tensor_insert(tensor(*arrs), *args, pos=(0, 1), arr_dims=arr_dims)
>>> np.allclose(r, tensor(args[0], arrs[0], args[1], arrs[1]))
True
>>> r = tensor_insert(tensor(*arrs), *args, pos=(0, 0), arr_dims=arr_dims)
>>> np.allclose(r, tensor(*args, *arrs))
True
>>> r = tensor_insert(tensor(*arrs), *args, pos=(1, 2), arr_dims=arr_dims)
>>> np.allclose(r, tensor(*np.insert(arrs, (1, 2), args, axis=0)))
True

tensor_merge(arr: numpy.ndarray, ins: numpy.ndarray, pos: Sequence[int], arr_dims: Sequence[Sequence[int]], ins_dims: Sequence[Sequence[int]], rank: int = 2, optimize: Union[bool, str] = False)numpy.ndarray[source]

For two tensor products arr and ins, merge ins into the product chain at indices pos. E.g, if $$\verb|arr|\equiv A\otimes B\otimes C$$, $$\verb|ins|\equiv D\otimes E$$, and $$\verb|pos|\equiv [1, 2]$$, the result will be the tensor product

$A\otimes D\otimes B\otimes E\otimes C.$

This function works in a similar way to numpy.insert() and tensor_insert().

Parameters
arr: ndarray

The tensor product in whose chain the other args should be inserted

ins: ndarray

The tensor product to be inserted in the product chain

pos: sequence of ints

The positions at which the constituent tensors of ins are inserted in the product chain. Should indicate the indices in the original tensor product chain that led to arr before which the constituents of ins should be inserted.

arr_dims: array_like, shape (rank, n_const)

The last rank dimensions of the n_const constituent tensors of the tensor product arr as a list of lists with the list at position i containing the i-th relevant dimension of all args. Since the remaing axes are broadcast over, their shape is irrelevant.

For example, if arr = tensor(a, b, c, rank=2) and a,b,c have shapes (2, 3, 4), (5, 2, 2, 1), (2, 2), arr_dims = [[3, 2, 2], [4, 1, 2]].

ins_dims: array_like, shape (rank, n_const)

The last rank dimensions of the n_const constituent tensors of the tensor product ins as a list of lists with the list at position i containing the i-th relevant dimension of ins. Since the remaing axes are broadcast over, their shape is irrelevant.

rank: int, optional (default: 2)

The rank of the tensors. E.g., for a Kronecker product between two vectors, rank == 1, and between two matrices rank == 2. The remaining axes are broadcast over.

optimize: bool|str, optional (default: False)

Optimize the tensor contraction order. Passed through to numpy.einsum().

numpy.insert

NumPy array insertion with similar syntax.

numpy.kron

NumPy tensor product.

tensor

tensor_insert

Insert array at given position in tensor product chain.

tensor_transpose

Transpose the order of a tensor product chain.

Examples

>>> I, X, Y, Z = paulis
>>> arr = tensor(X, Y, Z)
>>> ins = tensor(I, I)
>>> r1 = tensor_merge(arr, ins, pos=[1, 2], arr_dims=[[2]*3, [2]*3],
...                   ins_dims=[[2]*2, [2]*2])
>>> np.allclose(r1, tensor(X, I, Y, I, Z))
True
>>> r2 = tensor_merge(ins, arr, pos=[0, 1, 2], arr_dims=[[2]*2, [2]*2],
...                   ins_dims=[[2]*3, [2]*3])
>>> np.allclose(r1, r2)
True


tensor_insert() can provide the same functionality in some cases:

>>> arr = tensor(Y, Z)
>>> ins = tensor(I, X)
>>> r1 = tensor_merge(arr, ins, pos=[0, 0], arr_dims=[[2]*2, [2]*2],
...                   ins_dims=[[2]*2, [2]*2])
>>> r2 = tensor_insert(arr, I, X, pos=[0, 0], arr_dims=[[2]*2, [2]*2])
>>> np.allclose(r1, r2)
True


Also tensors of rank other than 2 and numpy broadcasting are supported:

>>> arr = np.random.randn(2, 10, 3, 4)
>>> ins = np.random.randn(2, 10, 3, 2)
>>> r = tensor_merge(tensor(*arr, rank=1), tensor(*ins, rank=1), [0, 1],
...                  arr_dims=[[4, 4]], ins_dims=[[2, 2]], rank=1)
>>> np.allclose(r, tensor(ins[0], arr[0], ins[1], arr[1], rank=1))
True

tensor_transpose(arr: numpy.ndarray, order: Sequence[int], arr_dims: Sequence[Sequence[int]], rank: int = 2)numpy.ndarray[source]

Transpose the order of a tensor product chain.

Parameters
arr: ndarray

The tensor product whose chain should be reordered.

order: sequence of ints

The transposition order. If arr == tensor(A, B) and order == (1, 0), the result will be tensor(B, A).

arr_dims: array_like, shape (rank, n_const)

The last rank dimensions of the n_const constituent tensors of the tensor product arr as a list of lists with the list at position i containing the i-th relevant dimension of all args. Since the remaing axes are broadcast over, their shape is irrelevant.

For example, if arr = tensor(a, b, c, rank=2) and a,b,c have shapes (2, 3, 4), (5, 2, 2, 1), (2, 2), arr_dims = [[3, 2, 2], [4, 1, 2]].

rank: int, optional (default: 2)

The rank of the tensors. E.g., for a Kronecker product between two vectors, rank == 1, and between two matrices rank == 2. The remaining axes are broadcast over.

Returns
transposed_arr: ndarray

The tensor product arr with its order transposed according to order

numpy.insert

NumPy array insertion with similar syntax.

numpy.kron

NumPy tensor product.

tensor

tensor_insert

Insert array at given position in tensor product chain.

tensor_merge

Merge tensor product chains.

Examples

>>> I, X, Y, Z = paulis
>>> arr = tensor(X, Y, Z)
>>> transposed = tensor_transpose(arr, [1, 2, 0], arr_dims=[[2, 2, 2]]*2)
>>> np.allclose(transposed, tensor(Y, Z, X))
True


## 2.11. Module contents¶

Package for efficient calculation of generalized filter functions

class Basis(basis_array: Sequence, traceless: Optional[bool] = None, btype: Optional[str] = None, labels: Optional[Sequence[str]] = None)[source]

Bases: numpy.ndarray

Class for operator bases. There are several ways to instantiate a Basis object:

• by just calling this constructor with a (possibly incomplete) array of basis matrices. No checks regarding orthonormality or hermiticity are performed.

• by calling one of the classes alternative constructors (classmethods):

These bases guarantee the following properties:

• hermitian

• orthonormal

• [traceless] (can be controlled by a flag)

Since Basis is a subclass of NumPy’s ndarray, it inherits all of its attributes, e.g. shape. The following attributes behave slightly differently to a ndarray, however

• A == B is True if all elements evaluate almost equal, i.e. equivalent to np.allclose(A, B).

• basis.T transposes the last two axes of basis. For a full basis, this corresponds to transposing each element individually. For a basis element, it corresponds to normal transposition.

Parameters
basis_array: array_like, shape (n, d, d)

An array or list of square matrices that are elements of an operator basis spanning $$\mathbb{C}^{d\times d}$$. n should be smaller than or equal to d**2.

traceless: bool, optional (default: auto)

Controls whether a traceless basis is forced. Here, traceless means that the first element of the basis is the identity and the remaining elements are matrices of trace zero. If an element of basis_array is neither traceless nor the identity and traceless == True, an exception will be raised. Defaults to True if basis_array is traceless and False if not.

btype: str, optional (default: ’custom’)

A string describing the basis type. For example, a basis created by the factory method pauli() has btype ‘pauli’.

labels: sequence of str, optional

A list of labels for the individual basis elements. Defaults to ‘C_0’, ‘C_1’, …

Attributes
Other than the attributes inherited from ndarray, a Basis
instance has the following attributes:
btype: str

Basis type.

labels: sequence of str

The labels for the basis elements.

d: int

Dimension of the space spanned by the basis.

H: Basis

Hermitian conjugate.

isherm: bool

If the basis is hermitian.

isorthonorm: bool

If the basis is orthonormal.

istraceless: bool

If the basis is traceless except for an identity element

iscomplete: bool

If the basis is complete, ie spans the full space.

sparse: COO, shape (n, d, d)

Representation in the COO format supplied by the sparse package.

four_element_traces: COO, shape (n, n, n, n)

Traces over all possible combinations of four elements of self. This is required for the calculation of the error transfer matrix and thus cached in the Basis instance.

Most of the attributes above are properties which are lazily
evaluated and cached.

Methods

 Other than the methods inherited from ndarray, a Basis instance has the following methods: normalize(b) Normalizes the basis (used internally when creating a basis from elements) tidyup(eps_scale=None) Cleans up floating point errors in-place to make zeros actual zeros. eps_scale is an optional argument multiplied to the data type’s eps to get the absolute tolerance.

Constructor.

property H: filter_functions.basis.Basis

Return the basis hermitian conjugated element-wise.

property T: filter_functions.basis.Basis

Return the basis transposed element-wise.

property four_element_traces: sparse._coo.core.COO

Return all traces of the form $$\mathrm{tr}(C_i C_j C_k C_l)$$ as a sparse COO array for $$i,j,k,l > 0$$ (i.e. excluding the identity).

classmethod from_partial(partial_basis_array: Sequence, traceless: Optional[bool] = None, btype: Optional[str] = None, labels: Optional[Sequence[str]] = None)filter_functions.basis.Basis[source]

Generate complete and orthonormal basis from a partial set.

The basis is completed using singular value decomposition to determine the null space of the expansion coefficients of the partial basis with respect to another complete basis.

Parameters
partial_basis_array: array_like

A sequence of basis elements.

traceless: bool, optional

If a traceless basis should be generated (i.e. the first element is the identity and all the others have trace zero).

btype: str, optional

A custom identifier.

labels: Sequence[str], optional

A list of custom labels for each element. If len(labels) == len(partial_basis_array), the newly created elements get labels ‘C_i’.

Returns
basis: Basis, shape (d**2, d, d)

The orthonormal basis.

Raises
ValueError

If the given elements are not orthonormal.

ValueError

If the given elements are not traceless but traceless==True.

ValueError

If not len(partial_basis_array) or d**2 labels were given.

classmethod ggm(d: int)filter_functions.basis.Basis[source]

Returns a generalized Gell-Mann basis in $$d$$ dimensions [Bert08] where the elements $$\Lambda_i$$ are normalized with respect to the Hilbert-Schmidt inner product,

$\begin{split}\langle\Lambda_i,\Lambda_j\rangle &= \mathrm{Tr}\,\Lambda_i^\dagger\Lambda_j \\ &= \delta_{ij}.\end{split}$
Parameters
d: int

The dimensionality of the space spanned by the basis

Returns
basis: Basis

The Basis object representing the GGM.

References

Bert08

Bertlmann, R. A., & Krammer, P. (2008). Bloch vectors for qudits. Journal of Physics A: Mathematical and Theoretical, 41(23). https://doi.org/10.1088/1751-8113/41/23/235303

property iscomplete: bool

Returns True if basis is complete.

property isherm: bool

Returns True if all basis elements are hermitian.

property isorthonorm: bool

Returns True if basis is orthonormal.

property istraceless: bool

Returns True if basis is traceless except for possibly the identity.

normalize(copy: bool = False)Union[None, filter_functions.basis.Basis][source]

Normalize the basis.

classmethod pauli(n: int)filter_functions.basis.Basis[source]

Returns a Pauli basis for $$n$$ qubits, i.e. the basis spans the space $$\mathbb{C}^{d\times d}$$ with $$d = 2^n$$:

$\mathcal{P} = \{I, X, Y, Z\}^{\otimes n}.$

The elements $$\sigma_i$$ are normalized with respect to the Hilbert-Schmidt inner product,

$\begin{split}\langle\sigma_i,\sigma_j\rangle &= \mathrm{Tr}\,\sigma_i^\dagger\sigma_j \\ &= \delta_{ij}.\end{split}$
Parameters
n: int

The number of qubits.

Returns
basis: Basis

The Basis object representing the Pauli basis.

property sparse: sparse._coo.core.COO

Return the basis as a sparse COO array

tidyup(eps_scale: Optional[float] = None)None[source]

Wraps util.remove_float_errors.

class PulseSequence(*args, **kwargs)[source]

Bases: object

A class for pulse sequences and their filter functions.

The Hamiltonian is separated into a control and a noise part with

$\begin{split}\mathcal{H}_c &= \sum_i a_i(t) A_i \\ \mathcal{H}_n &= \sum_j s_j(t) b_j(t) B_j\end{split}$

where $$A_i$$ and $$B_j$$ are hermitian operators and $$b_j(t)$$ are classically fluctuating noise variables captured in a power spectral density and not needed at instantiation of a PulseSequence.

Parameters
H_c: list of lists

A nested list of n_cops nested lists as taken by QuTiP functions (see for example qutip.propagator.propagator()) describing the control part of the Hamiltonian. The i-th entry of the list should be a list consisting of the i-th operator $$A_i$$ making up the control Hamiltonian and a list or array $$a_i(t)$$ describing the magnitude of that operator during the time intervals dt. Optionally, the list may also include operator identifiers. That is, H_c should look something like this:

H = [[c_oper1, c_coeff1, c_oper_identifier1],
[c_oper2, c_coeff2, c_oper_identifier2], ...]


The operators may be given either as NumPy arrays or QuTiP Qobjs and each coefficient array should have the same number of elements as dt, and should be given in units of $$\hbar$$. If not every sublist (read: operator) was given a identifier, they are automatically filled up with ‘A_i’ where i is the position of the operator.

H_n: list of lists

A nested list of n_nops nested lists as taken by QuTiP functions (see for example qutip.propagator.propagator()) describing the noise part of the Hamiltonian. The j-th entry of the list should be a list consisting of the j-th operator $$B_j$$ making up the noise Hamiltonian and a list or array describing the sensitivity $$s_j(t)$$ of the system to the noise operator during the time intervals dt. Optionally, the list may also include operator identifiers. That is, H_n should look something like this:

H = [[n_oper1, n_coeff1, n_oper_identifier1],
[n_oper2, n_coeff2, n_oper_identifier2], ...]


The operators may be given either as NumPy arrays or QuTiP Qobjs and each coefficient array should have the same number of elements as dt, and should be given in units of $$\hbar$$. If not every sublist (read: operator) was given a identifier, they are automatically filled up with ‘B_i’ where i is the position of the operator.

dt: array_like, shape (n_dt,)

The segment durations of the Hamiltonian (i.e. durations of constant control). Internally, the control operation is taken to start at $$t_0\equiv 0$$, i.e. the edges of the constant control segments are at times t = [0, *np.cumsum(dt)].

basis: Basis, shape (d**2, d, d), optional

The operator basis in which to calculate. If a Generalized Gell-Mann basis (see ggm()) is chosen, some calculations will be faster for large dimensions due to a simpler basis expansion. However, when extending the pulse sequence to larger qubit registers, cached filter functions cannot be retained since the GGM basis does not factor into tensor products. In this case a Pauli basis is preferable.

Notes

Due to the heavy use of NumPy’s einsum() function, results have a floating point error of ~1e-13.

Examples

A rotation by $$\pi$$ around the axis between x and y preceeded and followed by a period of free evolution with the system subject to dephasing noise.

>>> import qutip as qt; import numpy as np
>>> H_c = [[qt.sigmax(), [0, np.pi, 0]],
[qt.sigmay(), [0, np.pi, 0]]]
>>> # Equivalent pulse:
>>> # H_c = [[qt.sigmax() + qt.sigmay(), [0, np.pi, 0]]]
>>> # The noise sensitivity is constant
>>> H_n = [[qt.sigmaz()/np.sqrt(2), [1, 1, 1], 'Z']]
>>> dt = [1, 1, 1]
>>> # Free evolution between t=0 and t=1, rotation between t=1 and t=2,
>>> # and free evolution again from t=2 to t=3.
>>> pulse = PulseSequence(H_c, H_n, dt)
>>> pulse.c_oper_identifiers
['A_0', 'A_1']
>>> pulse.n_oper_identifiers
['Z']
>>> omega = np.logspace(-1, 2, 500)
>>> F = pulse.get_filter_function(omega)    # shape (1, 500)
>>> # Plot the resulting filter function:
>>> from filter_functions import plotting
>>> fig, ax, leg = plotting.plot_filter_function(pulse)

Attributes
c_opers: ndarray, shape (n_cops, d, d)

Control operators

n_opers: ndarray, shape (n_nops, d, d)

Noise operators

c_oper_identifers: sequence of str

Identifiers for the control operators of the system

n_oper_identifers: sequence of str

Identifiers for the noise operators of the system

c_coeffs: ndarray, shape (n_cops, n_dt)

Control parameters in units of $$\hbar$$

n_coeffs: ndarray, shape (n_nops, n_dt)

Noise sensitivities in units of $$\hbar$$

dt: ndarray, shape (n_dt,)

Time steps

t: ndarray, shape (n_dt + 1,)

Absolute times taken to start at $$t_0\equiv 0$$

tau: float

Total duration. Equal to t[-1].

d: int

Dimension of the Hamiltonian

basis: Basis, shape (d**2, d, d)

The operator basis used for calculation

nbytes: int

An estimate of the memory consumed by the PulseSequence instance and its attributes

If the Hamiltonian was diagonalized, the eigenvalues and -vectors as
well as the cumulative propagators are cached:
eigvals: ndarray, shape (n_dt, d)

Eigenvalues $$D^{(g)}$$

eigvecs: ndarray, shape (n_dt, d, d)

Eigenvectors $$V^{(g)}$$

propagators: ndarray, shape (n_dt+1, d, d)

Cumulative propagators $$Q_g$$

total_propagator: ndarray, shape (d, d)

The total propagator $$Q$$ of the pulse alone. That is, $$|\psi(\tau)\rangle = propagators|\psi(0)\rangle$$.

total_propagator_liouville: array_like, shape (d**2, d**2)

The transfer matrix for the total propagator of the pulse. Given by liouville_representation(pulse.total_propagator, pulse.basis).

Furthermore, when the filter function is calculated, the frequencies
are cached as well as other relevant quantities.

Methods

 get_filter_function_derivative(omega[, …]) Calculate the pulse sequence’s filter function derivative.
 cleanup(method=’conservative’) Delete cached attributes is_cached(attr) Checks if a given attribute of the PulseSequence is cached diagonalize() Diagonalize the Hamiltonian of the pulse sequence, computing eigenvalues and -vectors as well as cumulative propagators get_control_matrix(omega, show_progressbar=False) Calculate the control matrix for frequencies omega get_filter_function(omega, which=’fidelity’, show_progressbar=False) Calculate the filter function for frequencies omega get_pulse_correlation_filter_function(which=’fidelity’) Get the pulse correlation filter function (only possible if computed during concatenation) propagator_at_arb_t(t) Calculate the propagator at arbitrary times

Initialize a PulseSequence instance.

cache_control_matrix(omega: Sequence[float], control_matrix: Optional[numpy.ndarray] = None, show_progressbar: bool = False, cache_intermediates: bool = False)None[source]

Cache the control matrix and the frequencies it was calculated for.

Parameters
omega: array_like, shape (n_omega,)

The frequencies for which to cache the filter function.

control_matrix: array_like, shape ([n_nops,] n_nops, d**2, n_omega), optional

The control matrix for the frequencies omega. If None, it is computed.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations. Only applies if control_matrix is not supplied.

cache_filter_function(omega: Sequence[float], control_matrix: Optional[numpy.ndarray] = None, filter_function: Optional[numpy.ndarray] = None, which: str = 'fidelity', order: int = 1, show_progressbar: bool = False, cache_intermediates: bool = False)None[source]

Cache the filter function. If control_matrix.ndim == 4, it is taken to be the ‘pulse correlation control matrix’ and summed along the first axis. In that case, also the pulse correlation filter function is calculated and cached. Total phase factors and transfer matrices of the the cumulative propagator are also cached so they can be reused during concatenation.

Parameters
omega: array_like, shape (n_omega,)

The frequencies for which to cache the filter function.

control_matrix: array_like, shape ([n_nops,] n_nops, d**2, n_omega), optional

The control matrix for the frequencies omega. If None, it is computed and the filter function derived from it.

filter_function: array_like, shape (n_nops, n_nops, [d**2, d**2,] n_omega), optional

The filter function for the frequencies omega. If None, it is computed from R in the case order == 1 and from scratch else.

which: str, optional

Which filter function to cache. Either ‘fidelity’ (default) or ‘generalized’.

order: int, optional

First or second order filter function.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations.

cache_total_phases(omega: Sequence[float], total_phases: Optional[numpy.ndarray] = None)None[source]

Cache the total phase factors for this pulse and omega.

Parameters
omega: array_like, shape (n_omega,)

The frequencies for which to cache the phase factors.

total_phases: array_like, shape (n_omega,), optional

The total phase factors for the frequencies omega. If None, they are computed.

cleanup(method: str = 'conservative')None[source]

Delete cached byproducts of the calculation of the filter function that are not necessarily needed anymore in order to free up memory.

Parameters
method: optional

If set to ‘conservative’ (the default), only the following attributes are deleted:

• _eigvals

• _eigvecs

• _propagators

If set to ‘greedy’, all of the above as well as the following attributes are deleted:

• _total_propagator

• _total_propagator_liouville

• _total_phases

• _control_matrix

• _control_matrix_pc

If set to ‘all’, all of the above as well as the following attributes are deleted:

• omega

• _filter_function

• _filter_function_gen

• _filter_function_pc

• _filter_function_pc_gen

• _filter_function_2

• _intermediates[‘control_matrix_step’]

If set to ‘frequency dependent’ only attributes that are functions of frequency are initalized to None.

Note that if this PulseSequence is concatenated with another one, some of the attributes might need to be calculated again, resulting in slower execution of the concatenation.

diagonalize()None[source]

Diagonalize the Hamiltonian defining the pulse sequence.

property duration: float

The duration of the pulse. Alias of tau.

property eigvals: numpy.ndarray

Get the eigenvalues of the pulse’s Hamiltonian.

property eigvecs: numpy.ndarray

Get the eigenvectors of the pulse’s Hamiltonian.

get_control_matrix(omega: Sequence[float], show_progressbar: bool = False, cache_intermediates: bool = False)numpy.ndarray[source]

Get the control matrix for the frequencies omega. If it has been cached for the same frequencies, the cached version is returned, otherwise it is calculated from scratch.

Parameters
omega: array_like, shape (n_omega,)

The frequencies at which to evaluate the control matrix.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations.

Returns
control_matrix: ndarray, shape (n_nops, d**2, n_omega)

The control matrix for the noise operators.

get_filter_function(omega: Sequence[float], which: str = 'fidelity', order: int = 1, show_progressbar: bool = False, cache_intermediates: bool = False)numpy.ndarray[source]

Get the first or second order filter function.

The filter function is cached so it doesn’t need to be calculated twice for the same frequencies.

Parameters
omega: array_like, shape (n_omega,)

The frequencies at which to evaluate the filter function.

which: str, optional

Which filter function to return. Either ‘fidelity’ (default) or ‘generalized’ (see Notes). Only if order == 1.

order: int, optional

First or second order filter function.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

cache_intermediates: bool, optional

Keep intermediate terms of the calculation that are also required by other computations.

Returns
filter_function: ndarray, shape (n_nops, n_nops, [d**2, d**2,] n_omega)

The filter function for each combination of noise operators as a function of omega.

Notes

The first-order generalized filter function is given by

$F_{\alpha\beta,kl}(\omega) = \tilde{\mathcal{B}}_{\alpha k}^\ast(\omega) \tilde{\mathcal{B}}_{\beta l}(\omega),$

where $$\alpha,\beta$$ are indices counting the noise operators $$B_\alpha$$ and $$k,l$$ indices counting the basis elements $$C_k$$.

The fidelity filter function is obtained by tracing over the basis indices:

$F_{\alpha\beta}(\omega) = \sum_{k} F_{\alpha\beta,kk}(\omega).$
get_filter_function_derivative(omega: Sequence[float], control_identifiers: Optional[Sequence[str]] = None, n_coeffs_deriv: Optional[Sequence[Sequence[float]]] = None)numpy.ndarray[source]

Calculate the pulse sequence’s filter function derivative.

Parameters
omega: array_like, shape (n_omega,)

Frequencies at which the pulse control matrix is to be evaluated.

control_identifiers: Sequence[str]

Sequence of strings with the control identifiern to distinguish between control and drift Hamiltonian. The default is None.

n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)

The derivatives of the noise susceptibilities by the control amplitudes. Defaults to None.

Returns
filter_function_deriv: ndarray, shape (n_nops, n_t, n_ctrl, n_omega)

The regular filter functions’ derivatives for variation in each control contribution.

get_pulse_correlation_control_matrix()numpy.ndarray[source]

Get the pulse correlation control matrix if it was cached.

get_pulse_correlation_filter_function(which: str = 'fidelity')numpy.ndarray[source]

Get the pulse correlation filter function given by

$F_{\alpha\beta}^{(gg')}(\omega) = e^{i\omega(t_{g-1} - t_{g'-1})} \tilde{\mathcal{B}}^{(g)}(\omega)\mathcal{Q}^{(g-1)} \mathcal{Q}^{(g'-1)\dagger} \tilde{\mathcal{B}}^{(g')\dagger}(\omega),$

where $$g,g'$$ index the pulse in the sequence and $$\alpha,\beta$$ index the noise operators, if it was computed during concatenation. Since the calculation requires the individual pulse’s control matrices and phase factors, which are not retained after concatenation, the pulse correlation filter function cannot be computed afterwards.

Note that the frequencies for which the filter function was calculated are not stored.

Returns
filter_function_pc: ndarray, shape (n_pls, n_pls, n_nops, n_nops, n_omega)

The pulse correlation filter function for each noise operator as a function of omega. The first two axes correspond to the pulses in the sequence, i.e. if the concatenated pulse sequence is $$C\circ B\circ A$$, the first two axes are arranged like

$\begin{split}F_{\alpha\beta}^{(gg')} &= \begin{pmatrix} F_{\alpha\beta}^{(AA)} & F_{\alpha\beta}^{(AB)} & F_{\alpha\beta}^{(AC)} \\ F_{\alpha\beta}^{(BA)} & F_{\alpha\beta}^{(BB)} & F_{\alpha\beta}^{(BC)} \\ F_{\alpha\beta}^{(CA)} & F_{\alpha\beta}^{(CB)} & F_{\alpha\beta}^{(CC)} \end{pmatrix}\end{split}$

for $$g,g'\in\{A, B, C\}$$.

get_total_phases(omega: Sequence[float])numpy.ndarray[source]

Get the (cached) total phase factors for this pulse and omega.

is_cached(attr: str)bool[source]

Returns True if the attribute is cached

property nbytes: int

Return an estimate of the amount of memory consumed by this object (or, more precisely, the array attributes of this object).

property omega: numpy.ndarray

Cached frequencies

propagator_at_arb_t(t: Sequence[float])numpy.ndarray[source]

Calculate the cumulative propagator Q(t) at times t by making use of the fact that we assume piecewise-constant control.

property propagators: numpy.ndarray

Get the eigenvectors of the pulse’s Hamiltonian.

property t: numpy.ndarray

The times of the pulse.

property tau: Union[float, int]

The duration of the pulse.

property total_propagator: numpy.ndarray

Get total propagator of the pulse.

property total_propagator_liouville: numpy.ndarray

Get the transfer matrix for the total propagator of the pulse.

concatenate(pulses: , calc_pulse_correlation_FF: bool = False, calc_filter_function: Optional[bool] = None, which: str = 'fidelity', omega: Optional[Sequence[float]] = None, show_progressbar: bool = False)filter_functions.pulse_sequence.PulseSequence[source]

Concatenate an arbitrary number of pulses. Note that pulses are concatenated left-to-right, that is,

$\mathtt{concatenate((A, B))} \equiv B \circ A$

so that $$A$$ is executed before $$B$$ when applying the concatenated pulse.

Parameters
pulses: sequence of PulseSequences

The PulseSequence instances to be concatenated. If any of the instances have a cached filter function, the filter function for the composite pulse will also be calculated in order to make use of the speedup gained from concatenating the filter functions. If omega is given, calculation of the composite filter function is forced.

calc_pulse_correlation_FF: bool, optional

Switch to control whether the pulse correlation filter function (see PulseSequence.get_pulse_correlation_filter_function()) is calculated. If omega is not given, the cached frequencies of all pulses need to be equal.

calc_filter_function: bool, optional

Switch to force the calculation of the filter function to be carried out or not. Overrides the automatic behavior of calculating it if at least one pulse has a cached control matrix. If True and no pulse has a cached control matrix, a list of frequencies must be supplied as omega.

which: str, optional

Which filter function to compute. Either ‘fidelity’ (default) or ‘generalized’ (see PulseSequence.get_filter_function() and PulseSequence.get_pulse_correlation_filter_function()).

omega: array_like, optional

Frequencies at which to evaluate the (pulse correlation) filter functions. If None, an attempt is made to use cached frequencies.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

Returns
pulse: PulseSequence

The concatenated pulse.

concatenate_periodic(pulse: filter_functions.pulse_sequence.PulseSequence, repeats: int)filter_functions.pulse_sequence.PulseSequence[source]

Concatenate a pulse sequence pulse whose Hamiltonian is periodic repeats times. Although performing the same task, this function is much faster for concatenating many identical pulses with filter functions than concatenate().

Note that for large dimensions, the calculation of the control matrix using this function might be very memory intensive.

Parameters
pulse: PulseSequence

The PulseSequence instance to be repeated. If it has a cached filter function, the filter function for the new pulse will also be computed.

repeats: int

The number of repetitions

Returns
newpulse: PulseSequence

The concatenated PulseSequence

concatenate

Concatenate arbitrary PulseSequences.

Notes

The total control matrix is given by

$\begin{split}\tilde{\mathcal{B}}(\omega) &= \tilde{\mathcal{B}}^{(1)}(\omega) \sum_{g=0}^{G-1} \left(e^{i\omega T}\right)^g \\ &= \tilde{\mathcal{B}}^{(1)}(\omega)\bigl( \mathbb{I} - e^{i\omega T} \mathcal{Q}^{(1)}\bigr)^{-1}\bigl( \mathbb{I} - \bigl(e^{i\omega T} \mathcal{Q}^{(1)}\bigr)^G\bigr).\end{split}$

with $$T$$ the period of the control Hamiltonian and $$G$$ the number of periods. The last equality is valid only if $$\mathbb{I} - e^{i\omega T}\mathcal{Q}^{(1)}$$ is invertible.

error_transfer_matrix(pulse: Optional[PulseSequence] = None, spectrum: Optional[numpy.ndarray] = None, omega: Optional[Sequence[float]] = None, n_oper_identifiers: Optional[Sequence[str]] = None, second_order: bool = False, cumulant_function: Optional[numpy.ndarray] = None, show_progressbar: bool = False, memory_parsimonious: bool = False, cache_intermediates: bool = False)numpy.ndarray[source]

Compute the error transfer matrix up to unitary rotations.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to compute the error transfer matrix.

spectrum: array_like, shape ([[n_nops,] n_nops,] n_omega)

The two-sided noise power spectral density in units of inverse frequencies as an array of shape (n_omega,), (n_nops, n_omega), or (n_nops, n_nops, n_omega). In the first case, the same spectrum is taken for all noise operators, in the second, it is assumed that there are no correlations between different noise sources and thus there is one spectrum for each noise operator. In the third and most general case, there may be a spectrum for each pair of noise operators corresponding to the correlations between them. n_nops is the number of noise operators considered and should be equal to len(n_oper_identifiers).

omega: array_like, shape (n_omega,)

The frequencies at which to calculate the filter functions.

n_oper_identifiers: array_like, optional

The identifiers of the noise operators for which to evaluate the error transfer matrix. The default is all. Note that, since in general contributions from different noise operators won’t commute, not selecting all noise operators results in neglecting terms of order $$\xi^4$$.

second_order: bool, optional

Also take into account the frequency shifts $$\Delta$$ that correspond to second order Magnus expansion and constitute unitary terms. Default False.

cumulant_function: ndarray, shape ([[n_pls, n_pls,] n_nops,] n_nops, d**2, d**2)

A precomputed cumulant function. If given, pulse, spectrum, omega are not required.

show_progressbar: bool, optional

Show a progress bar for the calculation of the control matrix.

memory_parsimonious: bool, optional

Trade memory footprint for performance. See calculate_decay_amplitudes(). The default is False.

cache_intermediates: bool, optional

Keep and return intermediate terms of the calculation of the control matrix (if it is not already cached) that can be reused for second order or gradients. Can consume large amount of memory, but speed up the calculation.

Returns
error_transfer_matrix: ndarray, shape (d**2, d**2)

The error transfer matrix. The individual noise operator contributions are summed up before exponentiating as they might not commute.

calculate_cumulant_function

Calculate the cumulant function $$\mathcal{K}$$

calculate_decay_amplitudes

Calculate the $$\Gamma_{\alpha\beta,kl}$$

infidelity

Calculate only infidelity of a pulse.

Notes

The error transfer matrix is given by

$\tilde{\mathcal{U}} = \exp\mathcal{K}(\tau)$

with $$\mathcal{K}(\tau)$$ the cumulant function (see calculate_cumulant_function()). For Gaussian noise this expression is exact when taking into account the decay amplitudes $$\Gamma$$ and frequency shifts $$\Delta$$. As the latter effects coherent errors it can be neglected if we assume that the experimenter has calibrated their pulse.

For non-Gaussian noise the expression above is perturbative and includes noise up to order $$\xi^2$$ and hence $$\tilde{\mathcal{U}} = \mathbb{1} + \mathcal{K}(\tau) + \mathcal{O}(\xi^4)$$ (although it is evaluated as a matrix exponential in any case).

Given the above expression of the error transfer matrix, the entanglement fidelity is given by

$\mathcal{F}_\mathrm{e} = \frac{1}{d^2}\mathrm{tr}\,\tilde{\mathcal{U}}.$
extend(pulse_to_qubit_mapping: Sequence[Sequence[Optional[Union[filter_functions.pulse_sequence.PulseSequence, Sequence[int], int, Mapping[str, str]]]]], N: Optional[int] = None, d_per_qubit: int = 2, additional_noise_Hamiltonian: Optional[Sequence[Sequence[Union[numpy.ndarray, qutip.qobj.Qobj, Sequence[float]]]]] = None, cache_diagonalization: Optional[bool] = None, cache_filter_function: Optional[bool] = None, omega: Optional[Sequence[float]] = None, show_progressbar: bool = False)filter_functions.pulse_sequence.PulseSequence[source]

Map one or more pulse sequences to different qubits.

Parameters
pulse_to_qubit_mapping: sequence of mapping tuples

A sequence of tuples with the first entry a PulseSequence instance and the second an int or tuple of ints indicating the qubits that the PulseSequence should be mapped to. A mapping of operator identifiers may optionally be given as a third element of each tuple. By default, the index of the qubit the operator is mapped to is appended to its identifier.

Pulse sequences defined for multiple qubits may also be extended to non-neighboring qubits. Note that for multi-qubit pulses the order of the qubits is respected, i.e. mapping a pulse to (1, 0) is different from mapping it to (0, 1).

N: int

The total number of qubits the new PulseSequence should be defined for. By default, this is inferred from pulse_to_qubit_mapping.

d_per_qubit: int

The size of the Hilbert space a single qubit requires.

Additional noise operators and corresponding sensitivities for the new pulse sequence.

cache_diagonalization: bool

Force diagonalizing the new pulse sequence. By default, diagonalization is cached if all pulses in pulse_to_qubit_mapping have been diagonalized since it is much cheaper to get the relevant quantities as tensor products from the mapped pulses instead of diagonalizing the new pulse.

cache_filter_function: bool

Force computing the filter functions for the new pulse sequence. Noise operators of individual pulses will be extended to the new Hilbert space. By default, this is done if all pulses in pulse_to_qubit_mapping have their filter functions cached.

Note that extending the filter functions is only possible if they the mapped pulses are using a separable basis like the Pauli basis.

omega: array_like

Frequencies for which to compute the filter functions if cache_filter_function == True. Defaults to None, in which case the cached frequencies of the individual pulses need to be the same.

show_progressbar: bool

Show a progress bar for the calculation of the control matrix.

Returns
newpulse: PulseSequence

The new pulse sequence on the larger qubit register. The noise operators (and possibly filter functions) are stored in the following order: first those of the multi-qubit pulses in the order they appeared in pulse_to_qubit_mapping, then those of the single-qubit pulses, and lastly any additional ones that may be given by additional_noise_Hamiltonian.

remap

Map PulseSequence to a different qubit.

concatenate

Concatenate PulseSequences (in time).

concatenate_periodic

Periodically concatenate a PulseSequence.

Examples

>>> import filter_functions as ff
>>> I, X, Y, Z = ff.util.paulis
>>> X_pulse = ff.PulseSequence([[X, [np.pi/2], 'X']],
...                            [[X, [1], 'X'], [Z, [1], 'Z']],
...                            [1], basis=ff.Basis.pauli(1))
>>> XX_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)])
>>> XX_pulse.d
4
>>> XIX_pulse_1 = ff.extend([(X_pulse, 0), (X_pulse, 2)])
>>> XIX_pulse_1.d
8
>>> XXI_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)], N=3)
>>> XXI_pulse.d
8


Filter functions are automatically cached if they are for mapped pulses:

>>> omega = ff.util.get_sample_frequencies(X_pulse)
>>> X_pulse.cache_filter_function(omega)
>>> XX_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)])
>>> XX_pulse.is_cached('filter_function')
True


This behavior can also be overriden manually:

>>> XX_pulse = ff.extend([(X_pulse, 0), (X_pulse, 1)],
...                      cache_filter_function=False)
>>> XX_pulse.is_cached('filter_function')
False


Mapping pulses to non-neighboring qubits is also possible:

>>> Y_pulse = ff.PulseSequence([[Y, [np.pi/2], 'Y']],
...                            [[Y, [1], 'Y'], [Z, [1], 'Z']],
...                            [1], basis=ff.Basis.pauli(1))
>>> XXY_pulse = ff.extend([(XX_pulse, (0, 1)), (Y_pulse, 2)])
>>> XYX_pulse = ff.extend([(XX_pulse, (0, 2)), (Y_pulse, 1)])


Additionally, pulses can have the order of the qubits they are defined for permuted (see remap()):

>>> Z_pulse = ff.PulseSequence([[Z, [np.pi/2], 'Z']], [[Z, [1], 'Z']],
...                            [1], basis=ff.Basis.pauli(1))
>>> XY_pulse = ff.extend([(X_pulse, 0), (Y_pulse, 1)])
>>> YZX_pulse = ff.extend([(XY_pulse, (2, 0)), (Z_pulse, 1)])


Control and noise operator identifiers can be mapped according to a specified mapping:

>>> YX_pulse = ff.extend([(X_pulse, 1, {'X': 'IX', 'Z': 'IZ'}),
...                       (Y_pulse, 0, {'Y': 'YI', 'Z': 'ZI'})])
>>> YX_pulse.c_oper_identifiers
array(['IX', 'YI'], dtype='<U2')
>>> YX_pulse.n_oper_identifiers
array(['IX', 'IZ', 'YI', 'ZI'], dtype='<U2')


>>> H_n = [[ff.util.tensor(Z, Z, Z), [1], 'ZZZ']]
>>> XYX_pulse = ff.extend([(XX_pulse, (0, 2)), (Y_pulse, 1)],
>>> 'ZZZ' in XYX_pulse.n_oper_identifiers
True

infidelity_derivative(pulse: PulseSequence, spectrum: Sequence[float], omega: Sequence[float], control_identifiers: Optional[Sequence[str]] = None, n_coeffs_deriv: Optional[Sequence[Sequence[float]]] = None)numpy.ndarray[source]

Calculate the entanglement infidelity derivative of the PulseSequence pulse.

Parameters
pulse: PulseSequence

The PulseSequence instance for which to calculate the infidelity.

spectrum: array_like, shape ([[n_nops,] n_nops,] omega)

The two-sided noise power spectral density in units of inverse frequencies as an array of shape (n_omega,), (n_nops, n_omega), or (n_nops, n_nops, n_omega). In the first case, the same spectrum is taken for all noise operators, in the second, it is assumed that there are no correlations between different noise sources and thus there is one spectrum for each noise operator. In the third and most general case, there may be a spectrum for each pair of noise operators corresponding to the correlations between them. n_nops is the number of noise operators considered and should be equal to len(n_oper_identifiers).

omega: array_like, shape (n_omega,)

The frequencies at which the integration is to be carried out.

control_identifiers: Sequence[str], shape (n_ctrl,)

Sequence of strings with the control identifiern to distinguish between accessible control and drift Hamiltonian.

n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)

The derivatives of the noise susceptibilities by the control amplitudes. Defaults to None.

Returns
infid_deriv: ndarray, shape (n_nops, n_dt, n_ctrl)

Array with the derivative of the infidelity for each noise source taken for each control direction at each time step $$\frac{\partial I_e}{\partial u_h(t_{g'})}$$.

Raises
ValueError

If the provided noise spectral density does not fit expected shape.

Notes

The infidelity’s derivative is given by

$\frac{\partial I_e}{\partial u_h(t_{g'})} = \frac{1}{d} \int_{-\infty}^\infty \frac{d\omega}{2\pi} S_\alpha(\omega) \frac{\partial F_\alpha(\omega)} {\partial u_h(t_{g'})}$

with $$S_{\alpha}(\omega)$$ the noise spectral density and $$F_{\alpha}(\omega)$$ the canonical filter function for noise source $$\alpha$$.

To convert to the average gate infidelity, use the following relation given by Horodecki et al. [Hor99] and Nielsen [Nie02]:

$\big\langle\mathcal{I}_\mathrm{avg}\big\rangle = \frac{d}{d+1} \big\langle\mathcal{I}_\mathrm{e}\big\rangle.$

References

Hor99

Horodecki, M., Horodecki, P., & Horodecki, R. (1999). General teleportation channel, singlet fraction, and quasidistillation. Physical Review A - Atomic, Molecular, and Optical Physics, 60(3), 1888–1898. https://doi.org/10.1103/PhysRevA.60.1888

Nie02

Nielsen, M. A. (2002). A simple formula for the average gate fidelity of a quantum dynamical operation. Physics Letters, Section A: General, Atomic and Solid State Physics, 303(4), 249–252. https://doi.org/10.1016/S0375-9601(02)01272-0

liouville_representation(U: numpy.ndarray, basis: filter_functions.basis.Basis)numpy.ndarray[source]

Get the Liouville representaion of the unitary U with respect to the basis.

Parameters
U: ndarray, shape (…, d, d)

The unitary.

basis: Basis, shape (d**2, d, d)

The basis used for the representation, e.g. a Pauli basis.

Returns
R: ndarray, shape (…, d**2, d**2)

The Liouville representation of U.

Notes

The Liouville representation of a unitary quantum operation $$\mathcal{U}:\rho\rightarrow U\rho U^\dagger$$ is given by

$\mathcal{U}_{ij} = \mathrm{tr}(C_i U C_j U^\dagger)$

with $$C_i$$ elements of the basis spanning $$\mathbb{C}^{d\times d}$$ with $$d$$ the dimension of the Hilbert space.

remap(pulse: filter_functions.pulse_sequence.PulseSequence, order: Sequence[int], d_per_qubit: int = 2, oper_identifier_mapping: Optional[Mapping[str, str]] = None)filter_functions.pulse_sequence.PulseSequence[source]

Remap a PulseSequence by changing the order of qubits in the register. Cached attributes are automatically attempted to be retained.

Caution

This function simply permutes the order of the tensor product elements of control and noise operators. Thus, the resultant pulse will have its filter functions defined for different noise operators than the original one.

Parameters
pulse: PulseSequence

The pulse whose qubit order should be permuted.

order: sequence of ints

A list of permutation indices. E.g., if pulse is defined for two qubits, order == [1, 0] will reverse the order of qubits.

d_per_qubit: int (default: 2)

The size of the Hilbert space a single qubit inhabitates.

oper_identifier_mapping: dict_like

A mapping that maps operator identifiers from the old pulse to the remapped pulse. The default is the identity mapping.

Returns
remapped_pulse: PulseSequence

A new PulseSequence instance with the order of the qubits permuted according to order.

extend

Map PulseSequences to composite Hilbert spaces.

util.tensor_transpose

Transpose the order of a tensor product.

Examples

>>> X, Y = util.paulis[1:3]
>>> XY, YX = util.tensor(X, Y), util.tensor(Y, X)
>>> pulse = PulseSequence([[XY, [np.pi/2], 'XY']], [[YX, [1], 'YX']], [1],
...                       Basis.pauli(2))
>>> mapping = {'XY': 'YX', 'YX': 'XY'}
>>> remapped_pulse = remap(pulse, (1, 0), oper_identifier_mapping=mapping)
>>> target_pulse = PulseSequence([[YX, [np.pi/2], 'YX']],
...                              [[XY, [1], 'XY']], [1], Basis.pauli(2))
>>> remapped_pulse == target_pulse
True


Caching of attributes is automatically handled >>> remapped_pulse.is_cached(‘filter_function’) False >>> pulse.cache_filter_function(util.get_sample_frequencies(pulse)) >>> remapped_pulse = remap(pulse, (1, 0)) >>> remapped_pulse.is_cached(‘filter_function’) True