pygsti.modelmembers.operations.lindbladerrorgen
The LindbladErrorgen class and supporting functionality.
Module Contents
Classes
An Lindblad-form error generator. |
|
An object encapsulating a particular way of parameterizing a LindbladErrorgen |
Attributes
- pygsti.modelmembers.operations.lindbladerrorgen.IMAG_TOL = '1e-07'
- class pygsti.modelmembers.operations.lindbladerrorgen.LindbladErrorgen(lindblad_coefficient_blocks, lindblad_basis='auto', mx_basis='pp', evotype='default', state_space=None)
Bases:
pygsti.modelmembers.operations.linearop.LinearOperator
An Lindblad-form error generator.
This error generator consisting of terms that, with appropriate constraints ensurse that the resulting (after exponentiation) operation/layer operation is CPTP. These terms can be divided into “Hamiltonian”-type terms, which map rho -> i[H,rho] and “non-Hamiltonian”/”other”-type terms, which map rho -> A rho B + 0.5*(ABrho + rhoAB).
Initialize a LindbladErrorgen object.
Parameters
- lindblad_coefficient_blockslist of LindbladCoefficientBlock
A list of Lindblad coefficient blocks that define the error generator.
- lindblad_basis{‘auto’, ‘PP’, ‘std’, ‘gm’, ‘qt’} or Basis object, optional
The basis used for Lindblad terms. If ‘auto’, the basis is inferred from the coefficient blocks. Default is ‘auto’.
- mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object, optional
The basis for this error generator’s linear mapping. Default is ‘pp’.
- evotype{“default”, “densitymx”, “svterm”, “cterm”}, optional
The evolution type of the error generator being constructed. Default is “default”.
- state_spaceStateSpace, optional
The state space for the error generator. Default is None.
Raises
- ValueError
If the provided evotype does not support any of the required representations for a LindbladErrorgen.
- property total_term_magnitude
Get the total (sum) of the magnitudes of all this operator’s terms.
The magnitude of a term is the absolute value of its coefficient, so this function returns the number you’d get from summing up the absolute-coefficients of all the Taylor terms (at all orders!) you get from expanding this operator in a Taylor series.
Returns
float
- property total_term_magnitude_deriv
The derivative of the sum of all this operator’s terms.
Computes the derivative of the total (sum) of the magnitudes of all this operator’s terms with respect to the operators (local) parameters.
Returns
- numpy array
An array of length self.num_params
- property num_params
Get the number of independent parameters which specify this operation.
Returns
- int
the number of independent parameters.
- coefficient_blocks
- matrix_basis
- paramvals
- lindblad_term_superops_and_1norms
- combined_lindblad_term_superops
- classmethod from_operation_matrix_and_blocks(op_matrix, lindblad_coefficient_blocks, lindblad_basis='auto', mx_basis='pp', truncate=True, evotype='default', state_space=None)
Create a Lindblad-parameterized error generator from an operation matrix and coefficient blocks.
Parameters
- op_matrixnumpy array or SciPy sparse matrix
A square 2D array that gives the raw operation matrix, assumed to be in the mx_basis basis. The shape of this array sets the dimension of the operation.
- lindblad_coefficient_blockslist
A list of Lindblad coefficient blocks to set from the error generator projections.
- lindblad_basis{‘auto’, ‘PP’, ‘std’, ‘gm’, ‘qt’}, optional
The basis used for Lindblad terms. Default is ‘auto’.
- mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object, optional
The basis for this error generator’s linear mapping. Default is ‘pp’.
- truncatebool, optional
Whether to truncate the projections onto the Lindblad terms in order to meet constraints. Default is True. (e.g. to preserve CPTP) when necessary. If False, then an error is thrown when the Lindblad terms don’t conform to the constrains.
- evotype{“default”, “densitymx”, “svterm”, “cterm”}, optional
The evolution type of the error generator being constructed. Default is “default”.
- state_spaceStateSpace, optional
The state space for the error generator. Default is None.
Returns
LindbladErrorgen
- classmethod from_operation_matrix(op_matrix, parameterization='CPTP', lindblad_basis='PP', mx_basis='pp', truncate=True, evotype='default', state_space=None)
Creates a Lindblad-parameterized error generator from an operation.
Here “operation” means the exponentiated error generator, so this method essentially takes the matrix log of op_matrix and constructs an error generator from this using
from_error_generator()
.Parameters
- op_matrixnumpy array or SciPy sparse matrix
A square 2D array that gives the raw operation matrix, assumed to be in the mx_basis basis. The shape of this array sets the dimension of the operation.
- parameterizationstr, optional (default ‘CPTP’)
Describes how the Lindblad coefficients/projections relate to the error generator’s parameter values. Default is “CPTP”. Supported strings are those castable to LindbladParameterization. See LindbladParameterization for supported options.
- lindblad_basis{‘PP’, ‘std’, ‘gm’, ‘qt’}, optional
The basis used for Lindblad terms. Default is ‘PP’.
- mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object, optional
The basis for this error generator’s linear mapping. Default is ‘pp’. Allowed values are Matrix-unit (std), Gell-Mann (gm), Pauli-product (pp), and Qutrit (qt) (or a custom basis object).
- truncatebool, optional
Whether to truncate the projections onto the Lindblad terms in order to meet constraints. Default is True.
- evotype{“default”, “densitymx”, “svterm”, “cterm”}, optional
The evolution type of the error generator being constructed. Default is “default”.
- state_spaceStateSpace, optional
The state space for the error generator. Default is None.
- classmethod from_error_generator(errgen_or_dim, parameterization='CPTP', lindblad_basis='PP', mx_basis='pp', truncate=True, evotype='default', state_space=None)
Create a Lindblad-parameterized error generator from an error generator matrix or dimension.
Parameters
- errgen_or_dimnumpy array, SciPy sparse matrix, or int
A square 2D array that gives the full error generator or an integer specifying the dimension of a zero error generator.
- parameterizationstr, optional (default ‘CPTP’)
Describes how the Lindblad coefficients/projections relate to the error generator’s parameter values. Default is “CPTP”. Supported strings are those castable to LindbladParameterization. See LindbladParameterization for supported options.
- lindblad_basis{‘PP’, ‘std’, ‘gm’, ‘qt’}, optional
The basis used for Lindblad terms. Default is ‘PP’.
- mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object, optional
The basis for this error generator’s linear mapping. Default is ‘pp’. Allowed values are Matrix-unit (std), Gell-Mann (gm), Pauli-product (pp), and Qutrit (qt) (or a custom basis object).
- truncatebool, optional
Whether to truncate the projections onto the Lindblad terms in order to meet constraints. Default is True.
- evotype{“default”, “densitymx”, “svterm”, “cterm”}, optional
The evolution type of the error generator being constructed. Default is “default”.
- state_spaceStateSpace, optional
The state space for the error generator. Default is None.
Returns
LindbladErrorgen
- classmethod from_error_generator_and_blocks(errgen_or_dim, lindblad_coefficient_blocks, lindblad_basis='PP', mx_basis='pp', truncate=True, evotype='default', state_space=None)
Create a Lindblad-parameterized error generator from an error generator matrix or dimension and coefficient blocks.
Parameters
- errgen_or_dimnumpy array, SciPy sparse matrix, or int
A square 2D array that gives the full error generator or an integer specifying the dimension of a zero error generator.
- lindblad_coefficient_blockslist
A list of Lindblad coefficient blocks to set from the error generator projections.
- lindblad_basis{‘PP’, ‘std’, ‘gm’, ‘qt’}, optional
The basis used for Lindblad terms. Default is ‘PP’.
- mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object, optional
The basis for this error generator’s linear mapping. Default is ‘pp’.
- truncatebool, optional
Whether to truncate the projections onto the Lindblad terms in order to meet constraints. Default is True.
- evotype{“default”, “densitymx”, “svterm”, “cterm”}, optional
The evolution type of the error generator being constructed. Default is “default”.
- state_spaceStateSpace, optional
The state space for the error generator. Default is None.
Returns
LindbladErrorgen
- classmethod from_elementary_errorgens(elementary_errorgens, parameterization='auto', elementary_errorgen_basis='PP', mx_basis='pp', truncate=True, evotype='default', state_space=None)
Create a Lindblad-parameterized error generator from elementary error generators.
Parameters
- elementary_errorgensdict
A dictionary of elementary error generators. Keys are labels specifying the type and basis elements of the elementary error generators, and values are the corresponding coefficients. Keys are (termType, basisLabel1, <basisLabel2>) tuples, where termType is “H” (Hamiltonian), “S” (Stochastic), “C” (Correlation) or “A” (Active). Hamiltonian and Stochastic terms always have a single basis label (so key is a 2-tuple) whereas C and A tuples have 2 basis labels to specify off-diagonal non-Hamiltonian Lindblad terms. Basis labels are pauli strings. Values are coefficients.
- parameterizationstr, optional (default ‘CPTP’)
Describes how the Lindblad coefficients/projections relate to the error generator’s parameter values. Default is “CPTP”. Supported strings are those castable to LindbladParameterization. See LindbladParameterization for supported options.
- elementary_errorgen_basis{‘PP’, ‘std’, ‘gm’, ‘qt’}, optional
The basis used for the elementary error generators. Default is ‘PP’.
- mx_basis{‘std’, ‘gm’, ‘pp’, ‘qt’} or Basis object, optional
The basis for this error generator’s linear mapping. Default is ‘pp’.
- truncatebool, optional
Whether to truncate the projections onto the Lindblad terms in order to meet constraints. Default is True.
- evotype{“default”, “densitymx”, “svterm”, “cterm”}, optional
The evolution type of the error generator being constructed. Default is “default”.
- state_spaceStateSpace, optional
The state space for the error generator. Default is None.
Returns
LindbladErrorgen
- to_dense(on_space='minimal')
Return this error generator 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
- to_sparse(on_space='minimal')
Return the error generator as a sparse matrix.
Returns
scipy.sparse.csr_matrix
- 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()
.
- to_vector()
Extract a vector of the underlying operation parameters from this operation.
Returns
- numpy array
a 1D numpy array with length == num_params().
- from_vector(v, close=False, dirty_value=True)
Initialize the operation using a vector of parameters.
Parameters
- vnumpy array
The 1D vector of operation parameters. Length must == num_params()
- closebool, optional
Whether v is close to this operation’s current set of parameters. Under some circumstances, when this is true this call can be completed more quickly.
- dirty_valuebool, optional
The value to set this object’s “dirty flag” to before exiting this call. This is passed as an argument so it can be updated recursively. Leave this set to True unless you know what you’re doing.
Returns
None
- coefficients(return_basis=False, logscale_nonham=False)
Constructs a dictionary of the Lindblad-error-generator coefficients of this error generator.
Note that these are not necessarily the parameter values, as these coefficients are generally functions of the parameters (so as to keep the coefficients positive, for instance).
Parameters
- return_basisbool
Whether to also return a
Basis
containing the elements with which the error generator terms were constructed.- logscale_nonhambool, optional
Whether or not the non-hamiltonian error generator coefficients should be scaled so that the returned dict contains: (1 - exp(-d^2 * coeff)) / d^2 instead of coeff. This essentially converts the coefficient into a rate that is the contribution this term would have within a depolarizing channel where all stochastic generators had this same coefficient. This is the value returned by
error_rates()
.
Returns
- Ltermdictdict
Keys are (termType, basisLabel1, <basisLabel2>) tuples, where termType is “H” (Hamiltonian), “S” (Stochastic), “C” (Correlation) or “A” (Active). Hamiltonian and Stochastic terms always have a single basis label (so key is a 2-tuple) whereas C and A tuples have 2 basis labels to specify off-diagonal non-Hamiltonian Lindblad terms. Basis labels are pauli strings. Values are coefficients.
- basisBasis
A Basis mapping the basis labels used in the keys of Ltermdict to basis matrices.
- coefficient_labels()
The elementary error-generator labels corresponding to the elements of
coefficients_array()
.Returns
- tuple
A tuple of (<type>, <basisEl1> [,<basisEl2]) elements identifying the elementary error generators of this gate.
- coefficients_array()
The weighted coefficients of this error generator in terms of “standard” error generators.
Constructs a 1D array of all the coefficients returned by
coefficients()
, weighted so that different error generators can be weighted differently when a errorgen_penalty_factor is used in an objective function.Returns
- numpy.ndarray
A 1D array of length equal to the number of coefficients in the linear combination of standard error generators that is this error generator.
- coefficients_array_deriv_wrt_params()
The jacobian of
coefficients_array()
with respect to this error generator’s parameters.Returns
- numpy.ndarray
A 2D array of shape (num_coeffs, num_params) where num_coeffs is the number of coefficients in the linear combination of standard error generators that is this error generator, and num_params is this error generator’s number of parameters.
- error_rates()
Constructs a dictionary of the error rates associated with this error generator.
The error rates pertain to the channel formed by exponentiating this object.
The “error rate” for an individual Hamiltonian error is the angle about the “axis” (generalized in the multi-qubit case) corresponding to a particular basis element, i.e. theta in the unitary channel U = exp(i * theta/2 * BasisElement).
The “error rate” for an individual Stochastic error is the contribution that basis element’s term would have to the error rate of a depolarization channel. For example, if the rate corresponding to the term (‘S’,’X’) is 0.01 this means that the coefficient of the rho -> X*rho*X-rho error generator is set such that if this coefficient were used for all 3 (X,Y, and Z) terms the resulting depolarizing channel would have error rate 3*0.01 = 0.03.
Note that because error generator terms do not necessarily commute with one another, the sum of the returned error rates is not necessarily the error rate of the overall channel.
Returns
- Ltermdictdict
Keys are (termType, basisLabel1, <basisLabel2>) tuples, where termType is “H” (Hamiltonian), “S” (Stochastic), “C” (Correlation) or “A” (Active). Hamiltonian and Stochastic terms always have a single basis label (so key is a 2-tuple) whereas C and A tuples have 2 basis labels to specify off-diagonal non-Hamiltonian Lindblad terms. Basis labels are pauli strings. Values are coefficients. Values are real error rates except for the 2-basis-label case.
- set_coefficients(elementary_errorgens, action='update', logscale_nonham=False, truncate=True)
Sets the coefficients of elementary error generator terms in this error generator.
The dictionary lindblad_term_dict has tuple-keys describing the type of term and the basis elements used to construct it, e.g. (‘H’,’X’).
Parameters
- lindblad_term_dictdict
Keys are (termType, basisLabel1, <basisLabel2>) tuples, where termType is “H” (Hamiltonian), “S” (Stochastic), “C” (Correlation) or “A” (Active). Hamiltonian and Stochastic terms always have a single basis label (so key is a 2-tuple) whereas C and A tuples have 2 basis labels to specify off-diagonal non-Hamiltonian Lindblad terms. Basis labels are pauli strings.
- action{“update”,”add”,”reset”}
How the values in lindblad_term_dict should be combined with existing error-generator coefficients.
- logscale_nonhambool, optional
Whether or not the values in lindblad_term_dict for non-hamiltonian error generators should be interpreted as error rates (of an “equivalent” depolarizing channel, see
errorgen_coefficients()
) instead of raw coefficients. If True, then the non-hamiltonian coefficients are set to -log(1 - d^2*rate)/d^2, where rate is the corresponding value given in lindblad_term_dict. This is what is performed by the functionset_error_rates()
.- truncatebool, optional
Whether to truncate the projections onto the Lindblad terms in order to meet constraints (e.g. to preserve CPTP) when necessary. If False, then an error is thrown when the given coefficients cannot be parameterized as specified.
Returns
None
- set_error_rates(elementary_errorgens, action='update')
Sets the coeffcients of elementary error generator terms in this error generator.
Coefficients are set so that the contributions of the resulting channel’s error rate are given by the values in lindblad_term_dict. See
error_rates()
for more details.Parameters
- lindblad_term_dictdict
Keys are (termType, basisLabel1, <basisLabel2>) tuples, where termType is “H” (Hamiltonian), “S” (Stochastic), “C” (Correlation) or “A” (Active). Hamiltonian and Stochastic terms always have a single basis label (so key is a 2-tuple) whereas C and A tuples have 2 basis labels to specify off-diagonal non-Hamiltonian Lindblad terms. Basis labels are pauli strings.
- action{“update”,”add”,”reset”}
How the values in lindblad_term_dict should be combined with existing error rates.
Returns
None
- coefficient_weights(weights)
Get the non-default coefficient weights.
This method returns a dictionary of coefficient weights that are not equal to the default value of 1.0.
Parameters
- weightsdict
A dictionary where keys are coefficient labels and values are the corresponding weights.
Returns
- dict
A dictionary where keys are coefficient labels and values are the corresponding weights that are not equal to 1.0.
- set_coefficient_weights(weights)
Set the coefficient weights.
This method sets the weights for the coefficients of the error generator. If the coefficient weights array is not initialized, it initializes it to an array of ones.
Parameters
- weightsdict
A dictionary where keys are coefficient labels and values are the corresponding weights to set.
- transform_inplace(s)
Update error generator E with inv(s) * E * 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.
Parameters
- sGaugeGroupElement
A gauge group element which specifies the “s” matrix (and it’s inverse) used in the above similarity transform.
Returns
None
- deriv_wrt_params(wrt_filter=None)
The element-wise derivative this operation.
Construct a matrix whose columns are the vectorized derivatives of the flattened error generator matrix with respect to a single operator parameter. Thus, each column is of length op_dim^2 and there is one column per operation parameter.
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, shape == (dimension^2, num_params)
- hessian_wrt_params(wrt_filter1=None, wrt_filter2=None)
Construct the Hessian of this error generator 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)
- onenorm_upperbound()
Returns an upper bound on the 1-norm for this error generator (viewed as a matrix).
Returns
float
- to_memoized_dict(mmg_memo)
Create a serializable dict with references to other objects in the memo.
Parameters
- mmg_memo: dict
Memo dict from a ModelMemberGraph, i.e. keys are object ids and values are ModelMemberGraphNodes (which contain the serialize_id). This is NOT the same as other memos in ModelMember (e.g. copy, allocate_gpindices, etc.).
Returns
- mm_dict: dict
A dict representation of this ModelMember ready for serialization This must have at least the following fields: module, class, submembers, params, state_space, evotype Additional fields may be added by derived classes.
- class pygsti.modelmembers.operations.lindbladerrorgen.LindbladParameterization(block_types, param_modes, abbrev=None, meta=None)
Bases:
pygsti.baseobjs.nicelyserializable.NicelySerializable
An object encapsulating a particular way of parameterizing a LindbladErrorgen
A LindbladParameterization is a high-level parameterization-type (e.g. “H+S”) that contains two “modes” - one describing the number (and structure) of the non-Hamiltonian Lindblad coefficients (nonham_mode’) and one describing how the Lindblad coefficients are converted to/from parameters (`param_mode).
Parameters
- nonham_modestr
The “non-Hamiltonian mode” describes which non-Hamiltonian Lindblad coefficients are stored in a LindbladOp, and is one of “diagonal” (only the diagonal elements of the full coefficient matrix as a 1D array), “diag_affine” (a 2-by-d array of the diagonal coefficients on top of the affine projections), or “all” (the entire coefficient matrix).
- param_modestr
The “parameter mode” describes how the Lindblad coefficients/projections are converted into parameter values. This can be: “unconstrained” (coefficients are independent unconstrained parameters), “cptp” (independent parameters but constrained so map is CPTP), “depol” (all non-Ham. diagonal coeffs are the same, positive value), or “reldepol” (same as “depol” but no positivity constraint).
- ham_params_allowed, nonham_params_allowedbool
Whether or not Hamiltonian and non-Hamiltonian parameters are allowed.
- block_types
- param_modes
- abbrev = 'None'
- meta = 'None'
- classmethod minimal_from_elementary_errorgens(errs)
Helper function giving minimal Lindblad parameterization needed for specified errors.
Parameters
- errsdict
Error dictionary with keys as (termType, basisLabel) tuples, where termType can be “H” (Hamiltonian), “S” (Stochastic), or “A” (Active), and basisLabel is a string of I, X, Y, or Z, or to describe a Pauli basis element appropriate for the gate (i.e. having the same number of letters as there are qubits in the gate). For example, you could specify a 0.01-radian Z-rotation error and 0.05 rate of Pauli- stochastic X errors on a 1-qubit gate by using the error dictionary: {(‘H’,’Z’): 0.01, (‘S’,’X’): 0.05}.
Returns
- parameterizationstr
Parameterization string for constructing Lindblad error generators.