Running the Simulations

Setting up the Simulator Object

Documentation for the pycce.Simulator - main class for conducting CCE Simulations.

class Simulator(spin, position=None, alpha=None, beta=None, gyro=None, magnetic_field=None, D=None, E=0.0, r_dipole=None, order=None, bath=None, pulses=None, as_delay=False, n_clusters=None, **bath_kw)

The main class for CCE calculations.

The typical usage includes:

  1. Read array of the bath spins. This is done with Simulator.read_bath method which accepts either reading from .xyz file or from the BathArray instance with defined positions and names of the bath spins. In the process, the subset of the array within the distance of r_dipole from the central spin is taken and for this subset the Hyperfine couplings can be generated.

    If no hyperfine keyword is provided and there are some hyperfine couplings already, then no changes are done to the hyperfine tensors. If hyperfine='pd', the hyperfine couplings are computed assuming point dipole approximation. For all accepted arguments, see Simulator.read_bath.

  2. Generate set of clusters with Simulator.generate_clusters, determined by the maximum connectivity radius r_dipole and the maximum size of the cluster order.

  3. Compute the desired property with Simulator.compute method.

Note

Directly setting up the attribute values will rerun Simulator.read_bath and/or Simulator.generate_clusters to reflect updated value of the given attribute.

E.g. If Simulator.r_bath is set to some new value after initialization, then Simulator.read_bath and Simulator.generate_clusters are called with the increased bath.

First two steps are usually done during the initialization of the Simulator object by providing the necessary arguments.

Notes

Depending on the number of provided arguments, in the initialization process will call the following methods to setup the calculation engine.

  • If bath is provided, Simulator.read_bath is called with additional keywords in **bath_kw.

  • If both r_dipole and order are provided and bath is not None, the Simulator.generate_clusters is called.

See the corresponding method documentation for details.

Examples:

>>> atoms = random_bath('13C', 100, number=2000, seed=10)
>>> calc = Simulator(1, bath=atoms, r_bath=40, r_dipole=6,
>>>                  order=2, D=2.88 * 1e6,
>>>                  magnetic_field=500, pulses=1)
>>> print(calc)
Simulator for center array of size 1.
Parameters of cluster expansion:
r_bath: 40
r_dipole: 6
order: 2

Bath consists of 549 spins.

Clusters include:
549  clusters of order 1.
457  clusters of order 2.
Parameters:
  • spin (CenterArray or float or array with shape (n,)) –

    CenterArray containing properties of all central spins.

    OR

    Total spin of the central spin (Assumes one central spin).

    OR

    Array of total spins of the central spins (Assumes n central spins).

  • position (ndarray) – Cartesian coordinates ar array of coordinates in Angstrom of the central spin(s). Default (0., 0., 0.). If provided, overrides the position in CenterArray.

  • alpha (float or ndarray with shape (S, )) –

    \(\ket{0}\) state of the qubit in \(S_z\) basis or the index of eigenstate to be used as one.

    Default: Lowest energy eigenstate of the central spin Hamiltonian.

    If provided, overrides the alpha state in the CenterArray.

  • beta (float or ndarray with shape (S, )) –

    \(\ket{1}\) state of the qubit in \(S_z\) basis or the index of the eigenstate to be used as one.

    Default: Second lowest energy eigenstate of the central spin Hamiltonian.

    If provided, overrides the beta state in the CenterArray.

  • gyro (float or ndarray with shape (3, 3)) –

    Gyromagnetic ratio of the central spin(s) in rad / ms / G.

    OR

    Tensor describing central spin interactions with the magnetic field.

    Default -17608.597050 kHz * rad / G - gyromagnetic ratio of the free electron spin.

    If provided, overrides the gyro value in CenterArray.

  • D (float or ndarray with shape (3, 3)) –

    D (longitudinal splitting) parameter of central spin in ZFS tensor of central spin in kHz.

    OR

    Total ZFS tensor. Default 0.

    If provided, overrides the ZFS value in CenterArray.

  • E (float) – E (transverse splitting) parameter of central spin in ZFS tensor of central spin in kHz. Default 0. Ignored if D is None or tensor.

  • bath (ndarray or str) –

    First positional argument of the Simulator.read_bath method.

    Either:

    • Instance of BathArray class;

    • ndarray with dtype([('N', np.unicode_, 16), ('xyz', np.float64, (3,))]) containing names of bath spins (same ones as stored in self.ntype) and positions of the spins in angstroms;

    • the name of the .xyz text file containing 4 columns: name of the bath spin and xyz coordinates in A.

  • r_dipole (float) – Maximum connectivity distance between two bath spins.

  • order (int) – Maximum size of the cluster to be considered in CCE expansion.

  • n_clusters (dict) –

    Dictionary which contain maximum number of clusters of the given size. Has the form n_clusters = {order: number}, where order is the size of the cluster, number is the maximum number of clusters with this size.

    If provided, sort the clusters by the “strength” of cluster. Then the strongest number of clusters is taken.

    We define the strength of the cluster \(s\) as an inverse of the sum over inverse pairwise interaction strengths of the minimal cluster:

    \[s = (\sum_{i<j\in C} \frac{r^3}{\gamma_i\gamma_j})^{-1}\]

    Where \(\gamma_i\) is the gyromagnetic ration of a spin \(i\), \(r\) is the distance between two spins, and the summation of \(i, j\) goes only over the edges of the minimally connected cluster.

    We define minimally connected cluster as a cluster with lowest possible number of edges that still forms a connected graph. If multiple choices of the minimally connected cluster for the same cluster are possible, the one with the larger strength \(s\) is chosen.

  • pulses (list or int or Sequence) – Number of pulses in CPMG sequence or list with pulses.

  • **bath_kw – Additional keyword arguments for the Simulator.read_bath method.

center

Array of central spins.

Type:

CenterArray

clusters

Dictionary containing information about cluster structure of the bath.

Each keys n correspond to the size of the cluster. Each Simulator.clusters[n] contains ndarray of shape (m, n), where m is the number of clusters of given size, n is the size of the cluster. Each row of this array contains indexes of the bath spins included in the given cluster. Generated during .generate_clusters call.

Type:

dict

as_delay

True if time points are delay between pulses (for equispaced pulses), False if time points are total time. Ignored if pulses contains the time delays.

Type:

bool

interlaced

True if use hybrid CCE approach - for each cluster sample over states of the supercluster.

Type:

bool

seed

Seed for random number generator, used in random bath states sampling.

Type:

int

nbstates

Number or random bath states to sample over.

Type:

int

fixstates

If not None, shows which bath states to fix in random bath states.

Each key is the index of bath spin, value - fixed \(\hat S_z\) projection of the mixed state of nuclear spin.

Type:

dict

masked

True if mask numerically unstable points (with coherence > 1) in the averaging over bath states.

Note

It is up to user to check whether the possible instability is due to numerical error or unphysical assumptions of the calculations.

Type:

bool

second_order

True if add second order perturbation theory correction to the cluster Hamiltonian in conventional CCE. Relevant only for conventional CCE calculations.

Type:

bool

level_confidence

Maximum fidelity of the qubit state to be considered eigenstate of the central spin Hamiltonian when second_order set to True.

Type:

float

projected_bath_state

Array with z-projections of the bath spins states.

Type:

ndarray with shape (n,)

bath_state

Array of bath states.

Type:

bath_state (ndarray)

timespace

Time points at which compute the desired property.

Type:

timespace (ndarray with shape (n,))

fulldm

True if return full density matrix as the result in the gCCE calculations.

Type:

bool

property alpha

\(\ket{0}\) qubit state of the central spin in Sz basis OR index of the energy state to be considered as one.

Type:

ndarray or int

Type:

Returns .center.alpha property

property beta

\(\ket{1}\) qubit state of the central spin in Sz basis OR index of the energy state to be considered as one.

Type:

ndarray or int

Type:

Returns .center.beta property

property magnetic_field

Array containing external magnetic field as (Bx, By, Bz) or callable with signature magnetic_field(pos), where pos is an array with shape (3,) with the position of either bath or central spin. Default is (0, 0, 0).

Type:

ndarray

property order

Maximum size of the cluster.

Type:

int

property n_clusters

Dictionary which contain maximum number of clusters of the given size. Has the form n_clusters = {order: number}, where order is the size of the cluster, number is the maximum number of clusters with this size.

If provided, sorts the clusters by the strength of cluster interaction, equal to the inverse of a sum of inverse pairwise interaction in the minimally connected cluster. Then the strongest number of clusters is taken.

Type:

dict

property r_dipole

Maximum connectivity distance.

Type:

float

property pulses

List-like object, containing the sequence of the instantaneous ideal control pulses.

Each item is Pulse object, containing the following attributes:

  • which (array-like): Indexes of the central spins to be rotated by the pulse. Default is all.

  • x (float): Rotation angle of the central spin about x-axis in radians.

  • y (float): Rotation angle of the central spin about y-axis in radians.

  • z (float): Rotation angle of the central spin about z-axis in radians.

  • delay (float or ndarray): Delay before the pulse or array of delays with the same shape as time points.

Additionally, act as a container object for the pulses, applied to the bath.

The bath pulses can be accessed as items of the Pulse object, with name of the item corresponding to the name of the bath spin impacted, and the item corresponding to the BasePulse object with attributes:

  • x (float): Rotation angle of the central spin about x-axis in radians.

  • y (float): Rotation angle of the central spin about y-axis in radians.

  • z (float): Rotation angle of the central spin about z-axis in radians.

Examples

>>> p = Pulse('x', 'pi')
>>> print(p)
Pulse((x: 3.14, y: 0.00, z: 0.00))
>>> pb = Pulse('x', 'pi', bath_names=['13C', '14C'])
>>> print(pb)
Pulse((x: 3.14, y: 0.00, z: 0.00), {13C: (x: 3.14, y: 0.00, z: 0.00),
                                    14C: (x: 3.14, y: 0.00, z: 0.00)})
>>> print(pb['13C'])
(x: 3.14, y: 0.00, z: 0.00)

If delay is not provided in all pulses, assumes even delay of CPMG sequence. If only some delays are provided, assumes 0 delay in the pulses without delay provided.

For the full list of properties, see Pulse and Sequence documentations.

Type:

Sequence

property r_bath

Cutoff size of the spin bath. If len(r_bath) > 1, uses different cutoff sizes for each of the central spins. The total bath then is the sum of all bath spins, that are close to at least one of the central spins.

Type:

float or array-like

property external_bath

Array with spins read from DFT output (see pycce.io).

Type:

BathArray

property ext_r_bath

Maximum distance from the central spins of the bath spins for which to use the data from external_bath.

Type:

float

property error_range

Maximum distance between positions in bath and external bath to consider two positions the same (default 0.2).

Type:

float

property hyperfine

This argument tells the code how to generate hyperfine couplings. If (hyperfine = None and all A in provided bath are 0) or (hyperfine = 'pd'), use point dipole approximation. Otherwise can be an instance of Cube object, or callable with signature:

func(array)

where array is the BathArray object.

Type:

str, func, or Cube instance

property bath

Array of bath spins used in CCE simulations.

Type:

BathArray

set_zfs(D=None, E=0)

Set Zero Field Splitting of the central spin from longitudinal ZFS D and transverse ZFS E.

Parameters:
  • D (float or ndarray with shape (3, 3)) –

    D (longitudinal splitting) parameter of central spin in ZFS tensor of central spin in kHz.

    OR

    Total ZFS tensor. Default 0.

  • E (float) – E (transverse splitting) parameter of central spin in ZFS tensor of central spin in kHz. Default 0. Ignored if D is None or tensor.

set_magnetic_field(magnetic_field=None)

Set magnetic field from either value of the magnetic field along z-direction or full magnetic field vector.

Parameters:

magnetic_field (float or array-like) –

Magnetic field along z-axis.

OR

Array containing external magnetic field as (Bx, By, Bz). Default (0, 0, 0).

Reading the Bath

Documentation for the Simulator.read_bath and Simulator.generate_clusters method. These methods are called automatically on the initialization of the Simulator object if the necessary keywords are provided. Otherwise they can also be called by themselves to update the properties of the spin bath in Simulator object.

Simulator.read_bath(bath=None, r_bath=None, skiprows=1, external_bath=None, hyperfine=None, types=None, error_range=None, ext_r_bath=None, imap=None, func_kw=None)

Read spin bath from the file or from the BathArray.

Parameters:
  • bath (ndarray, BathArray or str) –

    Either:

    • Instance of BathArray class;

    • ndarray with dtype([('N', np.unicode_, 16), ('xyz', np.float64, (3,))]) containing names of bath spins (same ones as stored in self.ntype) and positions of the spins in angstroms;

    • the name of the xyz text file containing 4 cols: name of the bath spin and xyz coordinates in A.

  • r_bath (float or array-like) – Cutoff size of the spin bath. If len(r_bath) > 1, uses different cutoff sizes for each of the central spins. The total bath then is the sum of all bath spins, that are close to at least one of the central spins.

  • skiprows (int, optional) – If bath is name of the file, this argument gives number of rows to skip while reading the .xyz file (default 1).

  • external_bath (BathArray, optional) – BathArray containing spins read from DFT output (see pycce.io).

  • hyperfine (str, func, or Cube instance, optional) –

    This argument tells the code how to generate hyperfine couplings.

    If (hyperfine = None and all A in provided bath are 0) or (hyperfine = 'pd'), use point dipole approximation.

    Otherwise can be an instance of Cube object, or callable with signature:

    func(array, *args, **kwargs)
    

    where array is array of the bath spins,

  • func_kw (dict) – Additional keywords if for generating hyperfine couplings if hyperfine is callable.

  • types (SpinDict) –

    SpinDict or input to create one. Contains either SpinTypes of the bath spins or tuples which will initialize those.

    See pycce.bath.SpinDict documentation for details.

  • error_range (float, optional) – Maximum distance between positions in bath and external bath to consider two positions the same (default 0.2).

  • ext_r_bath (float, optional) – Maximum distance from the central spins of the bath spins for which to use the DFT positions.

  • imap (InteractionMap) –

    Instance of InteractionMap class, containing interaction tensors for bath spins. Each key of the InteractionMap is a tuple with indexes of two bath spins. The value is the 3x3 tensor describing the interaction between two spins in a format:

    \[I^iJI^j = I^i_{x}J_{xx}I^j_{x} + I^i_{x}J_{xy}I^j_{y} ...\]

Note

For each bath spin pair without interaction tensor in imap, coupling is approximated assuming magnetic point dipole–dipole interaction. If imap = None all interactions between bath spins are approximated in this way. Then interaction tensor between spins i and j is computed as:

\[\mathbf{J}_{ij} = -\gamma_{i} \gamma_{j} \frac{\hbar^2}{4\pi \mu_0} \left[ \frac{3 \vec{r}_{ij} \otimes \vec{r}_{ij} - |r_{ij}|^2 I}{|r_{ij}|^5} \right]\]

Where \(\gamma_{i}\) is gyromagnetic ratio of i spin, \(I\) is 3x3 identity matrix, and \(\vec{r}_{ij}\) is distance between two spins.

Returns:

The view of Simulator.bath attribute, generated by the method.

Return type:

BathArray

Simulator.generate_clusters(order=None, r_dipole=None, r_inner=0, strong=False, ignore=None, n_clusters=None)

Generate set of clusters used in CCE calculations.

The clusters are generated from the following procedure:

  • Each bath spin \(i\) forms a cluster of one.

  • Bath spins \(i\) and \(j\) form cluster of two if there is an edge between them (distance \(d_{ij} \le\) r_dipole).

  • Bath spins \(i\), \(j\), and \(j\) form a cluster of three if enough edges connect them (e.g., there are two edges \(ij\) and \(jk\))

  • And so on.

In general, we assume that spins \(\{i..n\}\) form clusters if they form a connected graph. Only clusters up to the size imposed by the order parameter (equal to CCE order) are included.

Parameters:
  • order (int) – Maximum size of the cluster.

  • r_dipole (float) – Maximum connectivity distance.

  • r_inner (float) – Minimum connectivity distance.

  • strong (bool) – True - generate only clusters with “strong” connectivity (all nodes should be interconnected). Default False.

  • ignore (list or str, optional) – If not None, includes the names of bath spins which are ignored in the cluster generation.

  • n_clusters (dict) –

    Dictionary which contain maximum number of clusters of the given size. Has the form n_clusters = {order: number}, where order is the size of the cluster, number is the maximum number of clusters with this size.

    If provided, sort the clusters by the strength of cluster interaction, Then the strongest number of clusters is taken.

    Strength of the cluster \(s\) is defined as an inverse of a sum of inverse pairwise interactions of the minimal cluster:

    \[s = (\sum_{i<j\in C} \frac{r^3}{\gamma_i\gamma_j})^{-1}\]

Returns:

View of Simulator.clusters. Simulator.clusters is a dictionary

with keys corresponding to size of the cluster.

I.e. Simulator.clusters[n] contains ndarray of shape (m, n), where m is the number of clusters of given size, n is the size of the cluster. Each row contains indexes of the bath spins included in the given cluster.

Return type:

dict

Calculate Properties with Simulator

Documentation for the Simulator.compute method - the interface to run calculations with PyCCE.

Simulator.compute(timespace, quantity='coherence', method='cce', **kwargs)

General function for computing properties with CCE.

The dynamics are simulated using the Hamiltonian:

\[\begin{split}&\hat H_S = \mathbf{SDS} + \mathbf{B\gamma}_{S}\mathbf{S} \\ &\hat H_{SB} = \sum_i \mathbf{S}\mathbf{A}_i\mathbf{I}_i \\ &\hat H_{B} = \sum_i{\mathbf{I}_i\mathbf{P}_i \mathbf{I}_i + \mathbf{B}\mathbf{\gamma}_i\mathbf{I}_i} + \sum_{i>j} \mathbf{I}_i\mathbf{J}_{ij}\mathbf{I}_j\end{split}\]

Here \(\hat H_S\) is the central spin Hamiltonian with Zero Field splitting tensor \(\mathbf{D}\) and gyromagnetic ratio tensor \(\mathbf{\gamma_S} = \mu_S \mathbf{g}\) are read from Simulator.zfs and Simulator.gyro respectively.

The \(\hat H_{SB}\) is the Hamiltonian describing interactions between central spin and the bath. The hyperfine coupling tensors \(\mathbf{A}_i\) are read from the BathArray stored in Simulator.bath['A']. They can be generated using point dipole approximation or provided by the user (see Simulator.read_bath for details).

The \(\hat H_{B}\) is the Hamiltonian describing interactions between the bath spins. The self interaction tensors \(\mathbf{P}_i\) are read from the BathArray stored in Simulator.bath['Q'] and have to be provided by the user.

The gyromagnetic ratios \(\mathbf{\gamma}_i\) are read from the BathArray.gyros attribuite, which is generated from the properties of the types of bath spins, stored in BathArray.types. They can either be provided by user or read from the pycce.common_isotopes object.

The interaction tensors \(\mathbf{J}_{ij}\) are assumed from point dipole approximation or can be provided in BathArray.imap attrubite.

Note

The compute method takes two keyword arguments to determine which quantity to compute and how:

  • method can take ‘cce’ or ‘gcce’ values, and determines which method to use - conventional or generalized CCE.

  • quantity can take ‘coherence’ or ‘noise’ values, and determines which quantity to compute - coherence function or autocorrelation function of the noise.

Each of the methods can be performed with monte carlo bath state sampling (if nbstates keyword is non zero) and with interlaced averaging (If interlaced keyword is set to True).

Examples

First set up Simulator object using random bath of 1000 13C nuclear spins.

>>> import pycce as pc
>>> import numpy as np
>>> atoms = pc.random_bath('13C', 100, number=2000, seed=10) # Random spin bath
>>> calc = pc.Simulator(1, bath=atoms, r_bath=40, r_dipole=6,
>>>                     order=2, D=2.88 * 1e6, # D of NV in GHz -> kHz
>>>                     magnetic_field=500, pulses=1)
>>> ts = np.linspace(0, 2, 101) # timesteps

We set magnetic field to 500 G along z-axis and chose 1 decoupling pulse (Hahn-echo) in this example. The zero field splitting is set to the one of NV center in diamond.

Run conventional CCE calculation at time points timespace to obtain coherence without second order effects:

>>> calc.compute(ts)

This will call Simulator.cce_coherence method with default keyword values.

Compute the coherence conventional CCE coherence with second order interactions between bath spins:

>>> calc.compute(ts, second_order=True)

Compute the coherence with conventional CCE with bath state sampling (over 10 states):

>>> calc.compute(ts, nbstates=10)

Compute the coherence with generalized CCE:

>>> calc.compute(ts, method='gcce')

This will call Simulator.gcce_dm method with default keyword values and obtain off diagonal element as \(\bra{0} \hat \rho_{C} \ket{1}\), where \(\hat \rho_{C}\) is the density matrix of the qubit.

Compute the coherence with generalized CCE with bath state sampling (over 10 states):

>>> calc.compute(ts, method='gcce', nbstates=10)
Parameters:
  • timespace (ndarray with shape (n,)) – Time points at which compute the desired property.

  • quantity (str) –

    Which quantity to compute. Case insensitive.

    Possible values:

    • ’coherence’: compute coherence function.

    • ’noise’: compute noise autocorrelation function.

  • method (str) –

    Which implementation of CCE to use. Case insensitive.

    Possible values:

    • ’cce’: conventional CCE, where interactions are mapped on 2 level pseudospin.

    • ’gcce’: Generalized CCE where central spin is included in each cluster.

  • magnetic_field (ndarray with shape (3,) or callable) –

    Magnetic field vector of form (Bx, By, Bz) or callable with signature magnetic_field(pos), where pos is an array with shape (3,) with the position of the spin.

    Default is None. Overrides Simulator.magnetic_field if provided.

  • pulses (list or int or Sequence) –

    Number of pulses in CPMG sequence.

    OR

    Sequence of the instantaneous ideal control pulses. It can be provided as an instance of Sequence class or a list with Pulse objects. (See documentation for pycce.Sequence).

    pulses can be provided as a list with tuples or dictionaries, each tuple or dictionary is used to initialize Pulse class instance.

    For example, for only central spin pulses the pulses argument can be provided as a list of tuples, containing:

    1. axis the rotation is about;

    2. angle of rotation;

    3. (optional) Time before the pulse. Can be as fixed, as well as varied. If varied, it should be provided as an array with the same length as timespace.

    E.g. for Hahn-Echo the pulses can be defined as [('x', np.pi)] or [('x', np.pi, timespace / 2)].

    Note

    If delay is not provided in all pulses, assumes even delay of CPMG sequence. If only some delays are provided, assumes delay = 0 in the pulses without delay.

    Then total experiment is assumed to be:

    tau – pulse – 2tau – pulse – … – 2tau – pulse – tau

    Where tau is the delay between pulses.

    The sum of delays at each time point should be less or equal to the total time of the experiment at the same time point, provided in timespace argument.

    Warning

    In conventional CCE calculations, only \(pi\) pulses on the central spin are allowed.

    In the calculations of noise autocorrelation this parameter is ignored.

    Default is None. Overrides``Simulator.pulses`` if provided.

  • i (int or ndarray with shape (2s+1, ) or callable) –

    Used in gCCE calculations. Along with j parameter indicates which density matrix element to compute with gCCE as:

    \[L=\bra{i}\hat \rho \ket{j}\]

    By default is equal to \(R\ket{0}\) state of the .center where \(R\) is a product of all rotations applied in the pulse sequence. Can be set as a vector in \(S_z\) basis, the index of the central spin Hamiltonian eigenstate, or as a callable with call signature i(dm), where dm is a density matrix of the central spin. If callable, j parameter is ignored.

  • j (int or ndarray with shape (2s+1, ) or callable) –

    Used in gCCE calculations. Along with i parameter indicates which density matrix element to compute.

    By default is equal to \(R\ket{1}\) state of the .center where \(R\) is a product of all rotations applied in the pulse sequence. Can be set as a vector in \(S_z\) basis, the index of the central spin Hamiltonian eigenstate, or as a callable with call signature j(dm), where dm is a density matrix of the central spin. If callable, i parameter is ignored.

  • as_delay (bool) –

    True if time points are delay between pulses (for equispaced pulses), False if time points are total time. Ignored in gCCE if pulses contains the time fractions. Conventional CCE calculations do not support custom time fractions.

    Default is False.

  • interlaced (bool) –

    True if use hybrid CCE approach - for each cluster sample over states of the supercluster.

    Default is False.

  • state (ndarray with shape (2s+1,)) –

    Initial state of the central spin, used in gCCE and noise autocorrelation calculations.

    Defaults to \(\frac{1}{N}(\ket{0} + \ket{1})\) if not set.

  • bath_state (array-like) –

    List of bath spin states. If len(shape) == 1, contains \(I_z\) projections of \(I_z\) eigenstates. Otherwise, contains array of initial density matrices of bath spins.

    Default is None. If not set, the code assumes completely random spin bath (density matrix of each nuclear spin is proportional to identity, \(\hat {\mathbb{1}}/N\)).

  • nbstates (int) –

    Number or random bath states to sample over.

    If provided, sampling of random states is carried and bath_states values are ignored.

    Default is 0.

  • seed (int) –

    Seed for random number generator, used in random bath states sampling.

    Default is None.

  • masked (bool) –

    True if mask numerically unstable points (with coherence > 1) in the averaging over bath states.

    Note

    It is up to user to check whether the possible instability is due to numerical error or unphysical assumptions of the calculations.

    Default is True for coherence calculations, False for noise calculations.

  • parallel_states (bool) –

    True if to use MPI to parallelize the calculations of density matrix equally over present mpi processes for random bath state sampling calculations.

    Compared to parallel keyword, when this argument is True each process is given a fraction of random bath states. This makes the implementation faster. Works best when the number of bath states is divisible by the number of processes, nbstates % size == 0.

    Default is False.

  • second_order (bool) –

    True if add second order perturbation theory correction to the cluster Hamiltonian in conventional CCE. Relevant only for conventional CCE calculations.

    If set to True sets the qubit states as eigenstates of central spin Hamiltonian from the following procedure. If qubit states are provided as vectors in \(S_z\) basis, for each qubit state compute the fidelity of the qubit state and all eigenstates of the central spin and chose the one with fidelity higher than level_confidence. If such state is not found, raises an error.

    Warning

    Second order corrections are not implemented as mean field.

    I.e., setting second_order=True and nbstates != 0 leads to the calculation, when mean field effect is accounted only from dipolar interactions within the bath.

    Default is False.

  • level_confidence (float) –

    Maximum fidelity of the qubit state to be considered eigenstate of the central spin Hamiltonian.

    Default is 0.95.

  • direct (bool) –

    True if use direct approach (requires way more memory but might be more numerically stable). False if use memory efficient approach.

    Default is False.

  • parallel (bool) –

    True if parallelize calculation of cluster contributions over different mpi processes.

    Default is False.

  • fulldm (bool) –

    Used in gCCE. True if return full density matrix of the central spin. False if only return the computed property.

    Default is False.

  • Returns – ndarray: Computed property.

Pulse sequences

Documentation of the Pulse and Sequence classes, used in definition of the complicated pulse sequences.

class BasePulse(x=None, y=None, z=None)

Base class for Pulse.

Parameters:
  • x (float) – Rotation angle about x-axis in radians.

  • y (float) – Rotation angle about y-axis in radians.

  • z (float) – Rotation angle about z-axis in radians.

set_angle(axis, angle)

Set rotation angle angle about axis axis.

Parameters:
  • axis (str) – Axis of the rotation.

  • angle (float) – Rotation angle in radians.

Returns:

check_flip()

Check if the rotation is about single cartesian axis by an angle \(\pi\).

property naxes

Number of axes the rotation is defined for.

Type:

int

property flip

True if the angle == pi.

Type:

bool

property x

Angle of rotation of the spin about x axis in rad.

Type:

float

property y

Angle of rotation of the spin about y axis in rad.

Type:

float

property z

Angle of rotation of the spin about z axis in rad.

Type:

float

generate_rotation(spinvec, spin_half=False)

Generate rotation matrix given spin vector.

Parameters:
  • spinvec (ndarray with shape (3, n, n)) – Spin vector.

  • spin_half (bool) – True if spin vector is for a spin-1/2. Default is False.

Returns:

Rotation operator.

Return type:

ndarray with shape (n, n)

class Pulse(axis=None, angle=None, delay=None, which=None, bath_names=None, bath_axes=None, bath_angles=None, **kwargs)

Class containing properties of each control pulse, applied to the system.

The properties of the pulse, applied on the central spin(s) can be accessed as attributes, while bath spin pulses can be acessed as elements of the Pulse instance.

Parameters:
  • axis (str) – Axis of rotation of the central spin. Can be ‘x’, ‘y’, or ‘z’. Default is None.

  • angle (float or str) – Angle of rotation of central spin. Can be provided in rad, or as a string, containing fraction of pi: 'pi', 'pi/2', '2*pi' etc. Default is None.

  • delay (float or ndarray) – Delay before the pulse or array of delays with the same shape as time points. Default is None.

  • which (array-like) – Indexes of the central spins to be rotated by the pulse. Default is all. Separated indexes are supported only if qubit states are provided separately for all center spins.

  • bath_names (str or array-like of str) – Name or array of names of bath spin types, impacted by the bath pulse. Default is None.

  • bath_axes (str or array-like of str) – Axis of rotation or array of axes of the bath spins. Default is None. If bath_names is provided, but bath_axes and bath_angles are not, assumes the same axis and angle as the one of the central spin

  • bath_angles (float or str or array-like) – Angle of rotation or array of axes of rotations of the bath spins.

  • x (float) – Rotation angle of the central spin about x-axis in radians.

  • y (float) – Rotation angle of the central spin about y-axis in radians.

  • z (float) – Rotation angle of the central spin about z-axis in radians.

Examples

>>> Pulse('x', 'pi')
Pulse((x: 3.14, y: 0.00, z: 0.00))
>>> Pulse('x', 'pi', bath_names=['13C', '14C'])
Pulse((x: 3.14, y: 0.00, z: 0.00), {13C: (x: 3.14, y: 0.00, z: 0.00), 14C: (x: 3.14, y: 0.00, z: 0.00)})
>>> import numpy as np
>>> p = Pulse('x', 'pi', delay=np.linspace(0, 1, 5), bath_names=['13C', '14C'],
>>>           bath_axes='x', bath_angles='pi/2')
>>> print(p)
Pulse((x: 3.14, y: 0.00, z: 0.00), {13C: (x: 1.57, y: 0.00, z: 0.00), 14C: (x: 1.57, y: 0.00, z: 0.00)},
t = [0.   0.25 0.5  0.75 1.  ])
>>> print(p['13C'])
(x: 1.57, y: 0.00, z: 0.00)
which

Indexes of the central spins to be rotated by the pulse.

Type:

iterable

bath_names

Array of names of bath spin types, impacted by the bath pulse.

Type:

ndarray

bath_axes

Array of axes of rotation of the bath spins.

Type:

ndarray

bath_angles

Array of angles of rotation of the bath spins.

Type:

ndarray

rotation

Matrix representation of the pulse for the given cluster. Generated by Run object.

Type:

ndarray

property delay

Delay or array of delays before the pulse.

Type:

ndarray or float

class Sequence(t=None)

List-like object, which contains the sequence of the pulses.

Each element is a Pulse instance, which can be generated from either the tuple with positional arguments or from the dictionary, or set manually.

If delay is not provided in all pulses in the sequence, assume equispaced pulse sequence:

t - pulse - 2t - pulse - 2t - … - pulse - t

If only some delays are provided, assumes 0 delay in the pulses without delay provided.

Examples

>>> import numpy as np
>>> Sequence([('x', np.pi, 0),
>>>           {'axis': 'y', 'angle': 'pi', 'delay': np.linspace(0, 1, 3), 'bath_names': '13C'},
>>>           Pulse('x', 'pi', 1)])
[Pulse((x, 3.14), t = 0), Pulse((y, 3.14), {13C: (y, 3.14)}, t = [0.  0.5 1. ]), Pulse((x, 3.14), t = 1)]
append(item)

S.append(value) – append value to the end of the sequence