pygsti.modelmembers.operations.linearop

The LinearOperator class and supporting functionality.

Module Contents

Classes

LinearOperator

Base class for all operation representations

Functions

finite_difference_deriv_wrt_params(operation, wrt_filter)

Computes a finite-difference Jacobian for a LinearOperator object.

finite_difference_hessian_wrt_params(operation, ...[, eps])

Computes a finite-difference Hessian for a LinearOperator object.

class pygsti.modelmembers.operations.linearop.LinearOperator(rep, evotype)

Bases: pygsti.modelmembers.modelmember.ModelMember

Base class for all operation representations

Parameters

repobject

A representation object containing the core data for this operator.

evotypeEvotype

The evolution type.

Attributes

sizeint

Return the number of independent elements in this operation (when viewed as a dense array)

dirtybool

Whether this object has been modified in a way that could have affected its parameters. A parent OpModel uses this information to know when it needs to refresh it’s model-wide parameter vector.

Initialize a new LinearOperator

property dim

Return the dimension of this operation (when viewed as a dense matrix)

Returns

int

property hilbert_schmidt_size

Return the number of independent elements in this operation as a dense Hilbert-Schmidt superoperator.

Returns

int

property chp_str

A string suitable for printing to a CHP input file after probabilistically selecting operation.

Returns
sstr

String of CHP code

set_dense(m)

Set the dense-matrix value of this operation.

Attempts to modify operation parameters so that the specified raw operation matrix becomes mx. Will raise ValueError if this operation is not possible.

Parameters
marray_like or LinearOperator

An array of shape (dim, dim) or LinearOperator representing the operation action.

Returns

None

set_time(t)

Sets the current time for a time-dependent operator.

For time-independent operators (the default), this function does nothing.

Parameters
tfloat

The current time.

Returns

None

abstract to_dense(on_space='minimal')

Return this operation as a dense matrix.

Parameters
on_space{‘minimal’, ‘Hilbert’, ‘HilbertSchmidt’}

The space that the returned dense operation acts upon. For unitary matrices and bra/ket vectors, use ‘Hilbert’. For superoperator matrices and super-bra/super-ket vectors use ‘HilbertSchmidt’. ‘minimal’ means that ‘Hilbert’ is used if possible given this operator’s evolution type, and otherwise ‘HilbertSchmidt’ is used.

Returns

numpy.ndarray

acton(state, on_space='minimal')

Act with this operator upon state

Parameters
stateState

The state to act on

Returns
State

The output state

abstract to_sparse(on_space='minimal')

Return this operation as a sparse matrix.

Returns

scipy.sparse.csr_matrix

abstract taylor_order_terms(order, max_polynomial_vars=100, return_coeff_polys=False)

Get the order-th order Taylor-expansion terms of this operation.

This function either constructs or returns a cached list of the terms at the given order. Each term is “rank-1”, meaning that its action on a density matrix rho can be written:

rho -> A rho B

The coefficients of these terms are typically polynomials of the operation’s parameters, where the polynomial’s variable indices index the global parameters of the operation’s parent (usually a Model), not the operation’s local parameter array (i.e. that returned from to_vector).

Parameters
orderint

The order of terms to get.

max_polynomial_varsint, optional

maximum number of variables the created polynomials can have.

return_coeff_polysbool

Whether a parallel list of locally-indexed (using variable indices corresponding to this object’s parameters rather than its parent’s) polynomial coefficients should be returned as well.

Returns
termslist

A list of RankOneTerm objects.

coefficientslist

Only present when return_coeff_polys == True. A list of compact polynomial objects, meaning that each element is a (vtape,ctape) 2-tuple formed by concatenating together the output of Polynomial.compact().

highmagnitude_terms(min_term_mag, force_firstorder=True, max_taylor_order=3, max_polynomial_vars=100)

Get terms with magnitude above min_term_mag.

Get the terms (from a Taylor expansion of this operator) that have magnitude above min_term_mag (the magnitude of a term is taken to be the absolute value of its coefficient), considering only those terms up to some maximum Taylor expansion order, max_taylor_order.

Note that this function also sets the magnitudes of the returned terms (by calling term.set_magnitude(…)) based on the current values of this operator’s parameters. This is an essential step to using these terms in pruned-path-integral calculations later on.

Parameters
min_term_magfloat

the threshold for term magnitudes: only terms with magnitudes above this value are returned.

force_firstorderbool, optional

if True, then always return all the first-order Taylor-series terms, even if they have magnitudes smaller than min_term_mag. This behavior is needed for using GST with pruned-term calculations, as we may begin with a guess model that has no error (all terms have zero magnitude!) and still need to compute a meaningful jacobian at this point.

max_taylor_orderint, optional

the maximum Taylor-order to consider when checking whether term- magnitudes exceed min_term_mag.

max_polynomial_varsint, optional

maximum number of variables the created polynomials can have.

Returns
highmag_termslist

A list of the high-magnitude terms that were found. These terms are sorted in descending order by term-magnitude.

first_order_indiceslist

A list of the indices into highmag_terms that mark which of these terms are first-order Taylor terms (useful when we’re forcing these terms to always be present).

taylor_order_terms_above_mag(order, max_polynomial_vars, min_term_mag)

Get the order-th order Taylor-expansion terms of this operation that have magnitude above min_term_mag.

This function constructs the terms at the given order which have a magnitude (given by the absolute value of their coefficient) that is greater than or equal to min_term_mag. It calls taylor_order_terms() internally, so that all the terms at order order are typically cached for future calls.

The coefficients of these terms are typically polynomials of the operation’s parameters, where the polynomial’s variable indices index the global parameters of the operation’s parent (usually a Model), not the operation’s local parameter array (i.e. that returned from to_vector).

Parameters
orderint

The order of terms to get (and filter).

max_polynomial_varsint, optional

maximum number of variables the created polynomials can have.

min_term_magfloat

the minimum term magnitude.

Returns
list

A list of Rank1Term objects.

frobeniusdist_squared(other_op, transform=None, inv_transform=None)

Return the squared frobenius difference between this operation and other_op

Optionally transforms this operation first using matrices transform and inv_transform. Specifically, this operation gets transfomed as: O => inv_transform * O * transform before comparison with other_op.

Parameters
other_opDenseOperator

The other operation.

transformnumpy.ndarray, optional

Transformation matrix.

inv_transformnumpy.ndarray, optional

Inverse of transform.

Returns

float

frobeniusdist(other_op, transform=None, inv_transform=None)

Return the frobenius distance between this operation and other_op.

Optionally transforms this operation first using matrices transform and inv_transform. Specifically, this operation gets transfomed as: O => inv_transform * O * transform before comparison with other_op.

Parameters
other_opDenseOperator

The other operation.

transformnumpy.ndarray, optional

Transformation matrix.

inv_transformnumpy.ndarray, optional

Inverse of transform.

Returns

float

residuals(other_op, transform=None, inv_transform=None)

The per-element difference between this DenseOperator and other_op.

Optionally, tansforming this operation first as O => inv_transform * O * transform.

Parameters
other_opDenseOperator

The operation to compare against.

transformnumpy.ndarray, optional

Transformation matrix.

inv_transformnumpy.ndarray, optional

Inverse of transform.

Returns
numpy.ndarray

A 1D-array of size equal to that of the flattened operation matrix.

jtracedist(other_op, transform=None, inv_transform=None)

Return the Jamiolkowski trace distance between this operation and other_op.

Optionally, tansforming this operation first as O => inv_transform * O * transform.

Parameters
other_opDenseOperator

The operation to compare against.

transformnumpy.ndarray, optional

Transformation matrix.

inv_transformnumpy.ndarray, optional

Inverse of transform.

Returns

float

diamonddist(other_op, transform=None, inv_transform=None)

Return the diamond distance between this operation and other_op.

Optionally, tansforming this operation first as O => inv_transform * O * transform.

Parameters
other_opDenseOperator

The operation to compare against.

transformnumpy.ndarray, optional

Transformation matrix.

inv_transformnumpy.ndarray, optional

Inverse of transform.

Returns

float

transform_inplace(s)

Update operation matrix O with inv(s) * O * s.

Generally, the transform function updates the parameters of the operation such that the resulting operation matrix is altered as described above. If such an update cannot be done (because the operation parameters do not allow for it), ValueError is raised.

In this particular case any transform of the appropriate dimension is possible, since all operation matrix elements are parameters.

Parameters
sGaugeGroupElement

A gauge group element which specifies the “s” matrix (and it’s inverse) used in the above similarity transform.

Returns

None

spam_transform_inplace(s, typ)

Update operation matrix O with inv(s) * O OR O * s, depending on the value of typ.

This functions as transform_inplace(…) but is used when this operation is used as a part of a SPAM vector. When typ == “prep”, the spam vector is assumed to be rho = dot(self, <spamvec>), which transforms as rho -> inv(s) * rho, so self -> inv(s) * self. When typ == “effect”, e.dag = dot(e.dag, self) (note that self is NOT self.dag here), and e.dag -> e.dag * s so that self -> self * s.

Parameters
sGaugeGroupElement

A gauge group element which specifies the “s” matrix (and it’s inverse) used in the above similarity transform.

typ{ ‘prep’, ‘effect’ }

Which type of SPAM vector is being transformed (see above).

Returns

None

depolarize(amount)

Depolarize this operation by the given amount.

Generally, the depolarize function updates the parameters of the operation such that the resulting operation matrix is depolarized. If such an update cannot be done (because the operation parameters do not allow for it), ValueError is raised.

Parameters
amountfloat or tuple

The amount to depolarize by. If a tuple, it must have length equal to one less than the dimension of the operation. In standard bases, depolarization corresponds to multiplying the operation matrix by a diagonal matrix whose first diagonal element (corresponding to the identity) equals 1.0 and whose subsequent elements (corresponding to non-identity basis elements) equal 1.0 - amount[i] (or just 1.0 - amount if amount is a float).

Returns

None

rotate(amount, mx_basis='gm')

Rotate this operation by the given amount.

Generally, the rotate function updates the parameters of the operation such that the resulting operation matrix is rotated. If such an update cannot be done (because the operation parameters do not allow for it), ValueError is raised.

Parameters
amounttuple of floats, optional

Specifies the rotation “coefficients” along each of the non-identity Pauli-product axes. The operation’s matrix G is composed with a rotation operation R (so G -> dot(R, G) ) where R is the unitary superoperator corresponding to the unitary operator U = exp( sum_k( i * rotate[k] / 2.0 * Pauli_k ) ). Here Pauli_k ranges over all of the non-identity un-normalized Pauli operators.

mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object

The source and destination basis, respectively. Allowed values are Matrix-unit (std), Gell-Mann (gm), Pauli-product (pp), and Qutrit (qt) (or a custom basis object).

Returns

None

deriv_wrt_params(wrt_filter=None)

The element-wise derivative this operation.

Constructs a matrix whose columns are the vectorized derivatives of the flattened operation matrix with respect to a single operation parameter. Thus, each column is of length op_dim^2 and there is one column per operation parameter. An empty 2D array in the StaticArbitraryOp case (num_params == 0).

Parameters
wrt_filterlist or numpy.ndarray

List of parameter indices to take derivative with respect to. (None means to use all the this operation’s parameters.)

Returns
numpy array

Array of derivatives with shape (dimension^2, num_params)

has_nonzero_hessian()

Whether this operation has a non-zero Hessian with respect to its parameters.

(i.e. whether it only depends linearly on its parameters or not)

Returns

bool

hessian_wrt_params(wrt_filter1=None, wrt_filter2=None)

Construct the Hessian of this operation with respect to its parameters.

This function returns a tensor whose first axis corresponds to the flattened operation matrix and whose 2nd and 3rd axes correspond to the parameters that are differentiated with respect to.

Parameters
wrt_filter1list or numpy.ndarray

List of parameter indices to take 1st derivatives with respect to. (None means to use all the this operation’s parameters.)

wrt_filter2list or numpy.ndarray

List of parameter indices to take 2nd derivatives with respect to. (None means to use all the this operation’s parameters.)

Returns
numpy array

Hessian with shape (dimension^2, num_params1, num_params2)

static convert_to_matrix(m)

Static method that converts a matrix-like object to a 2D numpy array.

Parameters
marray_like

matrix-like object

Returns

numpy array

pygsti.modelmembers.operations.linearop.finite_difference_deriv_wrt_params(operation, wrt_filter, eps=1e-07)

Computes a finite-difference Jacobian for a LinearOperator object.

The returned value is a matrix whose columns are the vectorized derivatives of the flattened operation matrix with respect to a single operation parameter, matching the format expected from the operation’s deriv_wrt_params method.

Parameters

operationLinearOperator

The operation object to compute a Jacobian for.

wrt_filterlist or numpy.ndarray

List of parameter indices to filter the result by (as though derivative is only taken with respect to these parameters).

epsfloat, optional

The finite difference step to use.

Returns

numpy.ndarray

An M by N matrix where M is the number of operation elements and N is the number of operation parameters.

pygsti.modelmembers.operations.linearop.finite_difference_hessian_wrt_params(operation, wrt_filter1, wrt_filter2, eps=0.0001)

Computes a finite-difference Hessian for a LinearOperator object.

The returned value is a tensor whose leading dimension corresponds to the elements of the flattened operation matrix, and whose remaining two dimensions correspond to derivatives with respect to the operations parameters (potentially filtered). This matches the format expected from the operation’s hessian_wrt_params method.

Parameters

operationLinearOperator

The operation object to compute a Hessian for.

wrt_filter1list or numpy.ndarray

List of parameter indices to filter the result by (as though the 1st derivative is only taken with respect to these parameters).

wrt_filter2list or numpy.ndarray

List of parameter indices to filter the result by (as though the 2nd derivative is only taken with respect to these parameters).

epsfloat, optional

The finite difference step to use.

Returns

numpy.ndarray

An M by N1 by N2 tensor where M is the number of operation elements and N1 and N2 are numbers of operation parameters.