pygsti.tools.symplectic
Symplectic representation utility functions
Module Contents
Functions
|
Creates the symplectic form for the number of qubits specified. |
|
Maps the input symplectic matrix between the 'standard' and 'directsum' symplectic form conventions. |
|
Checks whether a matrix is symplectic. |
Returns the inverse of a symplectic matrix over the integers mod 2. |
|
|
Returns the inverse of a Clifford gate in the symplectic representation. |
|
Checks if a symplectic matrix - phase vector pair (s,p) is the symplectic representation of a Clifford. |
|
Constructs a phase vector that, when paired with the provided symplectic matrix, defines a Clifford gate. |
|
Finds the Pauli layer that should be appended to a circuit to implement a given Clifford. |
|
Finds the Pauli layer that should be prepended to a circuit to implement a given Clifford. |
|
TODO: docstring |
|
TODO: docstring |
|
Multiplies two cliffords in the symplectic representation. |
|
Takes a kronecker product of symplectic representations. |
|
Contruct the (s,p) stabilizer representation for a computational basis state given by zvals. |
|
Applies a clifford in the symplectic representation to a stabilizer state in the standard stabilizer representation. |
|
Computes the probabilities of 0/1 (+/-) outcomes from measuring a Pauli operator on a stabilizer state. |
|
A helper routine used for manipulating stabilizer state representations. |
|
A helper routine used for manipulating stabilizer state representations. |
|
Compute the probability of a given outcome when measuring some or all of the qubits in a stabilizer state. |
|
Embeds the (s,p) Clifford symplectic representation into a larger symplectic representation. |
Creates a dictionary of the symplectic representations of 'standard' Clifford gates. |
|
|
Returns the symplectic representation of the composite Clifford implemented by the specified Clifford circuit. |
|
Constructs the symplectic representation of the n-qubit Clifford implemented by a single quantum circuit layer. |
Gives the group relationship between the 'I', 'H', 'P' 'HP', 'PH', and 'HPH' up-to-Paulis operators. |
|
|
Returns True if the unitary is a Clifford gate (w.r.t the standard basis), and False otherwise. |
|
Returns the symplectic representation of a one-qubit or two-qubit Clifford unitary. |
|
Returns a symplectic matrix of dimensions 2n x 2n sampled uniformly at random from the symplectic group S(n). |
|
Returns a Clifford, in the symplectic representation, sampled uniformly at random from the n-qubit Clifford group. |
|
Generates a uniformly random phase vector for a n-qubit Clifford. |
Get the bitstring corresponding to a Pauli. |
|
|
Applies a Clifford gate to the n-qubit Clifford gate specified by the 2n x 2n symplectic matrix. |
The number of Clifford gates in the n-qubit Clifford group. |
|
The number of elements in the symplectic group S(n) over the 2-element finite field. |
|
Returns the number of different cosets for the symplectic group S(n) over the 2-element finite field. |
|
|
Returns the symplectic inner product of two vectors in F_2^(2n). |
|
Applies transvection Z k to v. |
|
Converts integer i to an length n array of bits. |
|
Converts an n-bit string b to an integer between 0 and 2^`n` - 1. |
A utility function for selecting a random Clifford element. |
|
Returns the 2n x 2n symplectic matrix, over the finite field containing 0 and 1, with the "canonical" index i. |
|
|
Returns the "canonical" index of 2n x 2n symplectic matrix gn over the finite field containing 0 and 1. |
|
The index of a uniformly random 2n x 2n symplectic matrix over the finite field containing 0 and 1. |
- pygsti.tools.symplectic.symplectic_form(n, convention='standard')
Creates the symplectic form for the number of qubits specified.
There are two variants, of the sympletic form over the finite field of the integers modulo 2, used in pyGSTi. These corresponding to the ‘standard’ and ‘directsum’ conventions. In the case of ‘standard’, the symplectic form is the 2n x 2n matrix of ((0,1),(1,0)), where ‘1’ and ‘0’ are the identity and all-zeros matrices of size n x n. The ‘standard’ symplectic form is probably the most commonly used, and it is the definition used throughout most of the code, including the Clifford compilers. In the case of ‘directsum’, the symplectic form is the direct sum of n 2x2 bit-flip matrices. This is only used in pyGSTi for sampling from the symplectic group.
Parameters
- nint
The number of qubits the symplectic form should be constructed for. That is, the function creates a 2n x 2n matrix that is a sympletic form
- conventionstr, optional
Can be either ‘standard’ or ‘directsum’, which correspond to two different definitions for the symplectic form.
Returns
- numpy array
The specified symplectic form.
- pygsti.tools.symplectic.change_symplectic_form_convention(s, outconvention='standard')
Maps the input symplectic matrix between the ‘standard’ and ‘directsum’ symplectic form conventions.
That is, if the input is a symplectic matrix with respect to the ‘directsum’ convention and outconvention =’standard’ the output of this function is the equivalent symplectic matrix in the ‘standard’ symplectic form convention. Similarily, if the input is a symplectic matrix with respect to the ‘standard’ convention and outconvention = ‘directsum’ the output of this function is the equivalent symplectic matrix in the ‘directsum’ symplectic form convention.
Parameters
- snumpy.ndarray
The input symplectic matrix.
- outconventionstr, optional
Can be either ‘standard’ or ‘directsum’, which correspond to two different definitions for the symplectic form. This is the convention the input is being converted to (and so the input should be a symplectic matrix in the other convention).
Returns
- numpy array
The matrix s converted to outconvention.
- pygsti.tools.symplectic.check_symplectic(m, convention='standard')
Checks whether a matrix is symplectic.
Parameters
- mnumpy array
The matrix to check.
- conventionstr, optional
Can be either ‘standard’ or ‘directsum’, Specifies the convention of the symplectic form with respect to which the matrix should be sympletic.
Returns
- bool
A bool specifying whether the matrix is symplectic
- pygsti.tools.symplectic.inverse_symplectic(s)
Returns the inverse of a symplectic matrix over the integers mod 2.
Parameters
- snumpy array
The matrix to invert
Returns
- numpy array
The inverse of s, over the field of the integers mod 2.
- pygsti.tools.symplectic.inverse_clifford(s, p)
Returns the inverse of a Clifford gate in the symplectic representation.
This uses the formualas derived in Hostens and De Moor PRA 71, 042315 (2005).
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the Clifford
Returns
- sinversenumpy array
The symplectic matrix representing the inverse of the input Clifford.
- pinversenumpy array
The ‘phase vector’ representing the inverse of the input Clifford.
- pygsti.tools.symplectic.check_valid_clifford(s, p)
Checks if a symplectic matrix - phase vector pair (s,p) is the symplectic representation of a Clifford.
This uses the formualas derived in Hostens and De Moor PRA 71, 042315 (2005).
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the Clifford
Returns
- bool
True if (s,p) is the symplectic representation of some Clifford.
- pygsti.tools.symplectic.construct_valid_phase_vector(s, pseed)
Constructs a phase vector that, when paired with the provided symplectic matrix, defines a Clifford gate.
If the seed phase vector, when paired with s, represents some Clifford this seed is returned. Otherwise 1 mod 4 is added to the required elements of the pseed in order to make it at valid phase vector (which is one of many possible phase vectors that, together with s, define a valid Clifford).
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford
- pseednumpy array
The seed ‘phase vector’ over the integers mod 4.
Returns
- numpy array
Some p such that (s,p) is the symplectic representation of some Clifford.
- pygsti.tools.symplectic.find_postmultipled_pauli(s, p_implemented, p_target, qubit_labels=None)
Finds the Pauli layer that should be appended to a circuit to implement a given Clifford.
If some circuit implements the clifford described by the symplectic matrix s and the vector p_implemented, this function returns the Pauli layer that should be appended to this circuit to implement the clifford described by s and the vector p_target.
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford implemented by the circuit
- p_implementednumpy array
The ‘phase vector’ over the integers mod 4 representing the Clifford implemented by the circuit
- p_targetnumpy array
The ‘phase vector’ over the integers mod 4 that, together with s represents the Clifford that you want to implement. Together with s, this vector must define a valid Clifford.
- qubit_labelslist, optional
A list of qubit labels, that are strings or ints. The length of this list should be equal to the number of qubits the Clifford acts on. The ith element of the list is the label corresponding to the qubit at the ith index of s and the two phase vectors. If None, defaults to the integers from 0 to number of qubits - 1.
Returns
- list
A list that defines a Pauli layer, with the ith element containig one of the 4 tuples (P,qubit_labels[i]) with P = ‘I’, ‘Z’, ‘Y’ and ‘Z’
- pygsti.tools.symplectic.find_premultipled_pauli(s, p_implemented, p_target, qubit_labels=None)
Finds the Pauli layer that should be prepended to a circuit to implement a given Clifford.
If some circuit implements the clifford described by the symplectic matrix s and the vector p_implemented, this function returns the Pauli layer that should be prefixed to this circuit to implement the clifford described by s and the vector p_target.
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford implemented by the circuit
- p_implementednumpy array
The ‘phase vector’ over the integers mod 4 representing the Clifford implemented by the circuit
- p_targetnumpy array
The ‘phase vector’ over the integers mod 4 that, together with s represents the Clifford that you want to implement. Together with s, this vector must define a valid Clifford.
- qubit_labelslist, optional
A list of qubit labels, that are strings or ints. The length of this list should be equal to the number of qubits the Clifford acts on. The ith element of the list is the label corresponding to the qubit at the ith index of s and the two phase vectors. If None, defaults to the integers from 0 to number of qubits - 1.
Returns
- list
A list that defines a Pauli layer, with the ith element containig one of the 4 tuples (‘I’,i), (‘X’,i), (‘Y’,i), (‘Z’,i).
- pygsti.tools.symplectic.find_pauli_layer(pvec, qubit_labels, pauli_labels=None)
TODO: docstring pauli_labels defaults to [‘I’, ‘X’, ‘Y’, ‘Z’].
- pygsti.tools.symplectic.find_pauli_number(pvec)
TODO: docstring
- pygsti.tools.symplectic.compose_cliffords(s1, p1, s2, p2, do_checks=True)
Multiplies two cliffords in the symplectic representation.
The output corresponds to the symplectic representation of C2 times C1 (i.e., C1 acts first) where s1 (s2) and p1 (p2) are the symplectic matrix and phase vector, respectively, for Clifford C1 (C2). This uses the formualas derived in Hostens and De Moor PRA 71, 042315 (2005).
Parameters
- s1numpy array
The symplectic matrix over the integers mod 2 representing the first Clifford
- p1numpy array
The ‘phase vector’ over the integers mod 4 representing the first Clifford
- s2numpy array
The symplectic matrix over the integers mod 2 representing the second Clifford
- p2numpy array
The ‘phase vector’ over the integers mod 4 representing the second Clifford
- do_checksbool
If True (default), check inputs and output are valid cliffords. If False, these checks are skipped (for speed)
Returns
- snumpy array
The symplectic matrix over the integers mod 2 representing the composite Clifford
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the compsite Clifford
- pygsti.tools.symplectic.symplectic_kronecker(sp_factors)
Takes a kronecker product of symplectic representations.
Construct a single (s,p) symplectic (or stabilizer) representation that corresponds to the tensor (kronecker) product of the objects represented by each (s,p) element of sp_factors.
This is performed by inserting each factor’s s and p elements into the appropriate places of the final (large) s and p arrays. This operation works for combining Clifford operations AND also stabilizer states.
Parameters
- sp_factorsiterable
A list of (s,p) symplectic (or stabilizer) representation factors.
Returns
- snumpy.ndarray
An array of shape (2n,2n) where n is the total number of qubits (the sum of the number of qubits in each sp_factors element).
- pnumpy.ndarray
A 1D array of length 2n.
- pygsti.tools.symplectic.prep_stabilizer_state(nqubits, zvals=None)
Contruct the (s,p) stabilizer representation for a computational basis state given by zvals.
Parameters
- nqubitsint
Number of qubits
- zvalsiterable, optional
An iterable over anything that can be cast as True/False to indicate the 0/1 value of each qubit in the Z basis. If None, the all-zeros state is created. If None, then all zeros is assumed.
Returns
- s,pnumpy.ndarray
The stabilizer “matrix” and phase vector corresponding to the desired state. s has shape (2n,2n) (it includes antistabilizers) and p has shape 2n, where n equals nqubits.
- pygsti.tools.symplectic.apply_clifford_to_stabilizer_state(s, p, state_s, state_p)
Applies a clifford in the symplectic representation to a stabilizer state in the standard stabilizer representation.
The output corresponds to the stabilizer representation of the output state.
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the Clifford
- state_snumpy array
The matrix over the integers mod 2 representing the stabilizer state
- state_pnumpy array
The ‘phase vector’ over the integers mod 4 representing the stabilizer state
Returns
- out_snumpy array
The symplectic matrix over the integers mod 2 representing the output state
- out_pnumpy array
The ‘phase vector’ over the integers mod 4 representing the output state
- pygsti.tools.symplectic.pauli_z_measurement(state_s, state_p, qubit_index)
Computes the probabilities of 0/1 (+/-) outcomes from measuring a Pauli operator on a stabilizer state.
Parameters
- state_snumpy array
The matrix over the integers mod 2 representing the stabilizer state
- state_pnumpy array
The ‘phase vector’ over the integers mod 4 representing the stabilizer state
- qubit_indexint
The index of the qubit being measured
Returns
- p0, p1float
Probabilities of 0 (+ eigenvalue) and 1 (- eigenvalue) outcomes.
- state_s_0, state_s_1numpy array
Matrix over the integers mod 2 representing the output stabilizer states.
- state_p_0, state_p_1numpy array
Phase vectors over the integers mod 4 representing the output stabilizer states.
- pygsti.tools.symplectic.colsum(i, j, s, p, n)
A helper routine used for manipulating stabilizer state representations.
Updates the i-th stabilizer generator (column of s and element of p) with the group-action product of the j-th and the i-th generators, i.e.
generator[i] -> generator[j] + generator[i]
Parameters
- iint
Destination generator index.
- jint
Sournce generator index.
- snumpy array
The matrix over the integers mod 2 representing the stabilizer state
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the stabilizer state
- nint
The number of qubits. s must be shape (2n,2n) and p must be length 2n.
Returns
None
- pygsti.tools.symplectic.colsum_acc(acc_s, acc_p, j, s, p, n)
A helper routine used for manipulating stabilizer state representations.
Similar to
colsum()
except a separate “accumulator” column is used instead of the i-th column of s and element of p. I.e., this performs:acc[0] -> generator[j] + acc[0]
Parameters
- acc_snumpy array
The matrix over the integers mod 2 representing the “accumulator” stabilizer state
- acc_pnumpy array
The ‘phase vector’ over the integers mod 4 representing the “accumulator” stabilizer state
- jint
Index of the stabilizer generator being accumulated (see above).
- snumpy array
The matrix over the integers mod 2 representing the stabilizer state
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the stabilizer state
- nint
The number of qubits. s must be shape (2n,2n) and p must be length 2n.
Returns
None
- pygsti.tools.symplectic.stabilizer_measurement_prob(state_sp_tuple, moutcomes, qubit_filter=None, return_state=False)
Compute the probability of a given outcome when measuring some or all of the qubits in a stabilizer state.
Returns this probability, optionally along with the updated (post-measurement) stabilizer state.
Parameters
- state_sp_tupletuple
A (s,p) tuple giving the stabilizer state to measure.
- moutcomesarray-like
The z-values identifying which measurement outcome (a computational basis state) to compute the probability for.
- qubit_filteriterable, optional
If not None, a list of qubit indices which are measured. len(qubit_filter) should always equal len(moutcomes). If None, then assume all qubits are measured (len(moutcomes) == num_qubits).
- return_statebool, optional
Whether the post-measurement (w/outcome moutcomes) state is also returned.
Returns
- pfloat
The probability of the given measurement outcome.
- state_s,state_pnumpy.ndarray
Only returned when return_state=True. The post-measurement stabilizer state representation (an updated version of state_sp_tuple).
- pygsti.tools.symplectic.embed_clifford(s, p, qubit_inds, n)
Embeds the (s,p) Clifford symplectic representation into a larger symplectic representation.
The action of (s,p) takes place on the qubit indices specified by qubit_inds.
Parameters
- snumpy array
The symplectic matrix over the integers mod 2 representing the Clifford
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the Clifford
- qubit_indslist
A list or array of integers specifying which qubits s and p act on.
- nint
The total number of qubits
Returns
- snumpy array
The symplectic matrix over the integers mod 2 representing the embedded Clifford
- pnumpy array
The ‘phase vector’ over the integers mod 4 representing the embedded Clifford
- pygsti.tools.symplectic.compute_internal_gate_symplectic_representations(gllist=None)
Creates a dictionary of the symplectic representations of ‘standard’ Clifford gates.
Returns a dictionary containing the symplectic matrices and phase vectors that represent the specified ‘standard’ Clifford gates, or the representations of all the standard gates if no list of operation labels is supplied. These ‘standard’ Clifford gates are those gates that are already known to the code gates (e.g., the label ‘CNOT’ has a specfic meaning in the code), and are recorded as unitaries in “internalgates.py”.
Parameters
- gllistlist, optional
If not None, a list of strings corresponding to operation labels for any of the standard gates that have fixed meaning for the code (e.g., ‘CNOT’ corresponds to the CNOT gate with the first qubit the target). For example, this list could be gllist = [‘CNOT’,’H’,’P’,’I’,’X’].
Returns
- srep_dictdict
dictionary of (smatrix,svector) tuples, where smatrix and svector are numpy arrays containing the symplectic matrix and phase vector representing the operation label given by the key.
- pygsti.tools.symplectic.symplectic_rep_of_clifford_circuit(circuit, srep_dict=None, pspec=None)
Returns the symplectic representation of the composite Clifford implemented by the specified Clifford circuit.
This uses the formualas derived in Hostens and De Moor PRA 71, 042315 (2005).
Parameters
- circuitCircuit
The Clifford circuit to calculate the global action of, input as a Circuit object.
- srep_dictdict, optional
If not None, a dictionary providing the (symplectic matrix, phase vector) tuples associated with each operation label. If the circuit layer contains only ‘standard’ gates which have a hard-coded symplectic representation this may be None. Alternatively, if pspec is specifed and it contains the gates in circuit in a Clifford model, it also does not need to be specified (and it is ignored if it is specified). Otherwise it must be specified.
- pspecQubitProcessorSpec, optional
A QubitProcessorSpec that contains a Clifford model that defines the symplectic action of all of the gates in circuit. If this is not None it over-rides srep_dict. Both pspec and srep_dict can only be None if the circuit contains only gates with names that are hard-coded into pyGSTi.
Returns
- snumpy array
The symplectic matrix representing the Clifford implement by the input circuit
- pdictionary of numpy arrays
The phase vector representing the Clifford implement by the input circuit
- pygsti.tools.symplectic.symplectic_rep_of_clifford_layer(layer, n=None, q_labels=None, srep_dict=None, add_internal_sreps=True)
Constructs the symplectic representation of the n-qubit Clifford implemented by a single quantum circuit layer.
(Gates in a “single layer” must act on disjoint sets of qubits, but not all qubits need to be acted upon in the layer.)
Parameters
- layerLabel
A layer label, often a compound label with components. Specifies The Clifford gate(s) to calculate the global action of.
- nint, optional
The total number of qubits. Must be specified if q_labels is None.
- q_labelslist, optional
A list of all the qubit labels. If the layer is over qubits that are not labelled by integers 0 to n-1 then it is necessary to specify this list. Note that this should contain all the qubit labels for the circuit that this is a layer from, and they should be ordered as in that circuit, otherwise the symplectic rep returned might not be of the correct dimension or of the correct order.
- srep_dictdict, optional
If not None, a dictionary providing the (symplectic matrix, phase vector) tuples associated with each operation label. If the circuit layer contains only ‘standard’ gates which have a hard-coded symplectic representation this may be None. Otherwise it must be specified. If the layer contains some standard gates it is not necesary to specify the symplectic represenation for those gates.
- add_internal_srepsbool, optional
If True, the symplectic reps for internal gates are calculated and added to srep_dict. For speed, calculate these reps once, store them in srep_dict, and set this to False.
Returns
- snumpy array
The symplectic matrix representing the Clifford implement by specified circuit layer
- pnumpy array
The phase vector representing the Clifford implement by specified circuit layer
- pygsti.tools.symplectic.one_q_clifford_symplectic_group_relations()
Gives the group relationship between the ‘I’, ‘H’, ‘P’ ‘HP’, ‘PH’, and ‘HPH’ up-to-Paulis operators.
The returned dictionary contains keys (A,B) for all A and B in the above list. The value for key (A,B) is C if BA = C x some Pauli operator. E,g, (‘P’,’P’) = ‘I’.
This dictionary is important for Compiling multi-qubit Clifford gates without unneccessary 1-qubit gate over-heads. But note that this dictionary should not be used for compressing circuits containing these gates when the exact action of the circuit is of importance (not only the up-to-Paulis action of the circuit).
Returns
dict
- pygsti.tools.symplectic.unitary_is_clifford(unitary)
Returns True if the unitary is a Clifford gate (w.r.t the standard basis), and False otherwise.
Parameters
- unitarynumpy.ndarray
A unitary matrix to test.
Returns
bool
- pygsti.tools.symplectic.unitary_to_symplectic(u, flagnonclifford=True)
Returns the symplectic representation of a one-qubit or two-qubit Clifford unitary.
The Clifford is input as a complex matrix in the standard computational basis.
Parameters
- unumpy array
The unitary matrix to construct the symplectic representation for. This must be a one-qubit or two-qubit gate (so, it is a 2 x 2 or 4 x 4 matrix), and it must be provided in the standard computational basis. It also must be a Clifford gate in the standard sense.
- flagnoncliffordbool, opt
If True, a ValueError is raised when the input unitary is not a Clifford gate. If False, when the unitary is not a Clifford the returned s and p are None.
Returns
- snumpy array or None
The symplectic matrix representing the unitary, or None if the input unitary is not a Clifford and flagnonclifford is False
- pnumpy array or None
The phase vector representing the unitary, or None if the input unitary is not a Clifford and flagnonclifford is False
- pygsti.tools.symplectic.random_symplectic_matrix(n, convention='standard', rand_state=None)
Returns a symplectic matrix of dimensions 2n x 2n sampled uniformly at random from the symplectic group S(n).
This uses the method of Robert Koenig and John A. Smolin, presented in “How to efficiently select an arbitrary Clifford group element”.
Parameters
- nint
The size of the symplectic group to sample from.
- conventionstr, optional
Can be either ‘standard’ or ‘directsum’, which correspond to two different definitions for the symplectic form. In the case of ‘standard’, the symplectic form is the 2n x 2n matrix of ((0,1),(1,0)), where ‘1’ and ‘0’ are the identity and all-zeros matrices of size n x n. The ‘standard’ symplectic form is the convention used throughout most of the code. In the case of ‘directsum’, the symplectic form is the direct sum of n 2x2 bit-flip matrices.
- rand_state: RandomState, optional
A np.random.RandomState object for seeding RNG
Returns
- snumpy array
A uniformly sampled random symplectic matrix.
- pygsti.tools.symplectic.random_clifford(n, rand_state=None)
Returns a Clifford, in the symplectic representation, sampled uniformly at random from the n-qubit Clifford group.
The core of this function uses the method of Robert Koenig and John A. Smolin, presented in “How to efficiently select an arbitrary Clifford group element”, for sampling a uniformly random symplectic matrix.
Parameters
- nint
The number of qubits the Clifford group is over.
- rand_state: RandomState, optional
A np.random.RandomState object for seeding RNG
Returns
- snumpy array
The symplectic matrix representating the uniformly sampled random Clifford.
- pnumpy array
The phase vector representating the uniformly sampled random Clifford.
- pygsti.tools.symplectic.random_phase_vector(s, n, rand_state=None)
Generates a uniformly random phase vector for a n-qubit Clifford.
(This vector, together with the provided symplectic matrix, define a valid Clifford operation.) In combination with a uniformly random s the returned p defines a uniformly random Clifford gate.
Parameters
- snumpy array
The symplectic matrix to construct a random phase vector
- nint
The number of qubits the Clifford group is over.
- rand_state: RandomState, optional
A np.random.RandomState object for seeding RNG
Returns
- pnumpy array
A phase vector sampled uniformly at random from all those phase vectors that, as a pair with s, define a valid n-qubit Clifford.
- pygsti.tools.symplectic.bitstring_for_pauli(p)
Get the bitstring corresponding to a Pauli.
The state, represented by a bitstring, that the Pauli operator represented by the phase-vector p creates when acting on the standard input state.
Parameters
- pnumpy.ndarray
Phase vector of a symplectic representation, encoding a Pauli operation.
Returns
- list
A list of 0 or 1 elements.
- pygsti.tools.symplectic.apply_internal_gate_to_symplectic(s, gate_name, qindex_list, optype='row')
Applies a Clifford gate to the n-qubit Clifford gate specified by the 2n x 2n symplectic matrix.
The Clifford gate is specified by the internally hard-coded name gate_name. This gate is applied to the qubits with indices in qindex_list, where these indices are w.r.t to indeices of s. This gate is applied from the left (right) of s if optype is ‘row’ (‘column’), and has a row-action (column-action) on s. E.g., the Hadmard (‘H’) on qubit with index i swaps the ith row (or column) with the (i+n)th row (or column) of s; CNOT adds rows, etc.
Note that this function updates s, and returns None.
Parameters
- snp.array
A even-dimension square array over [0,1] that is the symplectic representation of some (normally multi-qubit) Clifford gate.
- gate_namestr
The gate name. Should be one of the gate-names of the hard-coded gates used internally in pyGSTi that is also a Clifford gate. Currently not all of those gates are supported, and gate_name must be one of: ‘H’, ‘P’, ‘CNOT’, ‘SWAP’.
- qindex_listlist or tuple
The qubit indices that gate_name acts on (can be either length 1 or 2 depending on whether the gate acts on 1 or 2 qubits).
- optype{‘row’, ‘column’}, optional
Whether the symplectic operator type uses rows or columns: TODO: docstring - better explanation.
Returns
None
- pygsti.tools.symplectic.compute_num_cliffords(n)
The number of Clifford gates in the n-qubit Clifford group.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- nint
The number of qubits the Clifford group is over.
Returns
- long integer
The cardinality of the n-qubit Clifford group.
- pygsti.tools.symplectic.compute_num_symplectics(n)
The number of elements in the symplectic group S(n) over the 2-element finite field.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- nint
S(n) group parameter.
Returns
int
- pygsti.tools.symplectic.compute_num_cosets(n)
Returns the number of different cosets for the symplectic group S(n) over the 2-element finite field.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- nint
S(n) group parameter.
Returns
int
- pygsti.tools.symplectic.symplectic_innerproduct(v, w)
Returns the symplectic inner product of two vectors in F_2^(2n).
Here F_2 is the finite field containing 0 and 1, and 2n is the length of the vectors. Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- vnumpy.ndarray
A length-2n vector.
- wnumpy.ndarray
A length-2n vector.
Returns
int
- pygsti.tools.symplectic.symplectic_transvection(k, v)
Applies transvection Z k to v.
Code from “How to efficiently select an arbitrary Clifford group element by Robert Koenig and John A. Smolin.
Parameters
- knumpy.ndarray
A length-2n vector.
- vnumpy.ndarray
A length-2n vector.
Returns
numpy.ndarray
- pygsti.tools.symplectic.int_to_bitstring(i, n)
Converts integer i to an length n array of bits.
Code from “How to efficiently select an arbitrary Clifford group element by Robert Koenig and John A. Smolin.
Parameters
- iint
Any integer.
- nint
Number of bits
Returns
- numpy.ndarray
Integer array of 0s and 1s.
- pygsti.tools.symplectic.bitstring_to_int(b, n)
Converts an n-bit string b to an integer between 0 and 2^`n` - 1.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- blist, tuple, or array
Sequence of bits (a bitstring).
- nint
Number of bits.
Returns
int
- pygsti.tools.symplectic.find_symplectic_transvection(x, y)
A utility function for selecting a random Clifford element.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- xnumpy.ndarray
A length-2n vector.
- ynumpy.ndarray
A length-2n vector.
Returns
numpy.ndarray
- pygsti.tools.symplectic.compute_symplectic_matrix(i, n)
Returns the 2n x 2n symplectic matrix, over the finite field containing 0 and 1, with the “canonical” index i.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- iint
Canonical index.
- nint
Number of qubits.
Returns
numpy.ndarray
- pygsti.tools.symplectic.compute_symplectic_label(gn, n=None)
Returns the “canonical” index of 2n x 2n symplectic matrix gn over the finite field containing 0 and 1.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- gnnumpy.ndarray
symplectic matrix
- nint, optional
Number of qubits (if None, use gn.shape[0] // 2).
Returns
- int
The canonical index of gn.
- pygsti.tools.symplectic.random_symplectic_index(n, rand_state=None)
The index of a uniformly random 2n x 2n symplectic matrix over the finite field containing 0 and 1.
Code from “How to efficiently select an arbitrary Clifford group element” by Robert Koenig and John A. Smolin.
Parameters
- nint
Number of qubits (half dimension of symplectic matrix).
- rand_state: RandomState, optional
A np.random.RandomState object for seeding RNG
Returns
numpy.ndarray