pygsti.algorithms.gaugeopt

GST gauge optimization algorithms

Module Contents

Functions

gaugeopt_to_target(model, target_model, item_weights=None, cptp_penalty_factor=0, spam_penalty_factor=0, gates_metric='frobenius', spam_metric='frobenius', gauge_group=None, method='auto', maxiter=100000, maxfev=None, tol=1e-08, oob_check_interval=0, return_all=False, comm=None, verbosity=0, check_jac=False)

Optimize the gauge degrees of freedom of a model to that of a target.

gaugeopt_custom(model, objective_fn, gauge_group=None, method='L-BFGS-B', maxiter=100000, maxfev=None, tol=1e-08, oob_check_interval=0, return_all=False, jacobian_fn=None, comm=None, verbosity=0)

Optimize the gauge of a model using a custom objective function.

_create_objective_fn(model, target_model, item_weights=None, cptp_penalty_factor=0, spam_penalty_factor=0, gates_metric='frobenius', spam_metric='frobenius', method=None, comm=None, check_jac=False)

Creates the objective function and jacobian (if available)

_cptp_penalty_size(mdl)

Helper function - same as that in core.py.

_spam_penalty_size(mdl)

Helper function - same as that in core.py.

_cptp_penalty(mdl, prefactor, op_basis)

Helper function - CPTP penalty: (sum of tracenorms of gates),

_spam_penalty(mdl, prefactor, op_basis)

Helper function - CPTP penalty: (sum of tracenorms of gates),

_cptp_penalty_jac_fill(cp_penalty_vec_grad_to_fill, mdl_pre, mdl_post, gauge_group_el, prefactor, op_basis, wrt_filter)

Helper function - jacobian of CPTP penalty (sum of tracenorms of gates)

_spam_penalty_jac_fill(spam_penalty_vec_grad_to_fill, mdl_pre, mdl_post, gauge_group_el, prefactor, op_basis, wrt_filter)

Helper function - jacobian of CPTP penalty (sum of tracenorms of gates)

pygsti.algorithms.gaugeopt.gaugeopt_to_target(model, target_model, item_weights=None, cptp_penalty_factor=0, spam_penalty_factor=0, gates_metric='frobenius', spam_metric='frobenius', gauge_group=None, method='auto', maxiter=100000, maxfev=None, tol=1e-08, oob_check_interval=0, return_all=False, comm=None, verbosity=0, check_jac=False)

Optimize the gauge degrees of freedom of a model to that of a target.

Parameters
  • model (Model) – The model to gauge-optimize

  • target_model (Model) – The model to optimize to. The metric used for comparing models is given by gates_metric and spam_metric.

  • item_weights (dict, optional) – Dictionary of weighting factors for gates and spam operators. Keys can be gate, state preparation, or POVM effect, as well as the special values “spam” or “gates” which apply the given weighting to all spam operators or gates respectively. Values are floating point numbers. Values given for specific gates or spam operators take precedence over “gates” and “spam” values. The precise use of these weights depends on the model metric(s) being used.

  • cptp_penalty_factor (float, optional) – If greater than zero, the objective function also contains CPTP penalty terms which penalize non-CPTP-ness of the gates being optimized. This factor multiplies these CPTP penalty terms.

  • spam_penalty_factor (float, optional) – If greater than zero, the objective function also contains SPAM penalty terms which penalize non-positive-ness of the state preps being optimized. This factor multiplies these SPAM penalty terms.

  • gates_metric ({"frobenius", "fidelity", "tracedist"}, optional) – The metric used to compare gates within models. “frobenius” computes the normalized sqrt(sum-of-squared-differences), with weights multiplying the squared differences (see Model.frobeniusdist()). “fidelity” and “tracedist” sum the individual infidelities or trace distances of each gate, weighted by the weights.

  • spam_metric ({"frobenius", "fidelity", "tracedist"}, optional) – The metric used to compare spam vectors within models. “frobenius” computes the normalized sqrt(sum-of-squared-differences), with weights multiplying the squared differences (see Model.frobeniusdist()). “fidelity” and “tracedist” sum the individual infidelities or trace distances of each “SPAM gate”, weighted by the weights.

  • gauge_group (GaugeGroup, optional) – The gauge group which defines which gauge trasformations are optimized over. If None, then the model’s default gauge group is used.

  • method (string, optional) –

    The method used to optimize the objective function. Can be any method known by scipy.optimize.minimize such as ‘BFGS’, ‘Nelder-Mead’, ‘CG’, ‘L-BFGS-B’, or additionally:

    • ’auto’ – ‘ls’ when allowed, otherwise ‘L-BFGS-B’

    • ’ls’ – custom least-squares optimizer.

    • ’custom’ – custom CG that often works better than ‘CG’

    • ’supersimplex’ – repeated application of ‘Nelder-Mead’ to converge it

    • ’basinhopping’ – scipy.optimize.basinhopping using L-BFGS-B as a local optimizer

    • ’swarm’ – particle swarm global optimization algorithm

    • ’evolve’ – evolutionary global optimization algorithm using DEAP

    • ’brute’ – Experimental: scipy.optimize.brute using 4 points along each dimensions

  • maxiter (int, optional) – Maximum number of iterations for the gauge optimization.

  • maxfev (int, optional) – Maximum number of function evaluations for the gauge optimization. Defaults to maxiter.

  • tol (float, optional) – The tolerance for the gauge optimization.

  • oob_check_interval (int, optional) – If greater than zero, gauge transformations are allowed to fail (by raising any exception) to indicate an out-of-bounds condition that the gauge optimizer will avoid. If zero, then any gauge-transform failures just terminate the optimization.

  • return_all (bool, optional) – When True, return best “goodness” value and gauge matrix in addition to the gauge optimized model.

  • comm (mpi4py.MPI.Comm, optional) – When not None, an MPI communicator for distributing the computation across multiple processors.

  • verbosity (int, optional) – How much detail to send to stdout.

  • check_jac (bool) – When True, check least squares analytic jacobian against finite differences

Returns

  • model if return_all == False

  • (goodnessMin, gaugeMx, model) if return_all == True – where goodnessMin is the minimum value of the goodness function (the best ‘goodness’) found, gaugeMx is the gauge matrix used to transform the model, and model is the final gauge-transformed model.

pygsti.algorithms.gaugeopt.gaugeopt_custom(model, objective_fn, gauge_group=None, method='L-BFGS-B', maxiter=100000, maxfev=None, tol=1e-08, oob_check_interval=0, return_all=False, jacobian_fn=None, comm=None, verbosity=0)

Optimize the gauge of a model using a custom objective function.

Parameters
  • model (Model) – The model to gauge-optimize

  • objective_fn (function) – The function to be minimized. The function must take a single Model argument and return a float.

  • gauge_group (GaugeGroup, optional) – The gauge group which defines which gauge trasformations are optimized over. If None, then the model’s default gauge group is used.

  • method (string, optional) –

    The method used to optimize the objective function. Can be any method known by scipy.optimize.minimize such as ‘BFGS’, ‘Nelder-Mead’, ‘CG’, ‘L-BFGS-B’, or additionally:

    • ’custom’ – custom CG that often works better than ‘CG’

    • ’supersimplex’ – repeated application of ‘Nelder-Mead’ to converge it

    • ’basinhopping’ – scipy.optimize.basinhopping using L-BFGS-B as a local optimizer

    • ’swarm’ – particle swarm global optimization algorithm

    • ’evolve’ – evolutionary global optimization algorithm using DEAP

    • ’brute’ – Experimental: scipy.optimize.brute using 4 points along each dimensions

  • maxiter (int, optional) – Maximum number of iterations for the gauge optimization.

  • maxfev (int, optional) – Maximum number of function evaluations for the gauge optimization. Defaults to maxiter.

  • tol (float, optional) – The tolerance for the gauge optimization.

  • oob_check_interval (int, optional) – If greater than zero, gauge transformations are allowed to fail (by raising any exception) to indicate an out-of-bounds condition that the gauge optimizer will avoid. If zero, then any gauge-transform failures just terminate the optimization.

  • return_all (bool, optional) – When True, return best “goodness” value and gauge matrix in addition to the gauge optimized model.

  • jacobian_fn (function, optional) – The jacobian of objective_fn. The function must take three parameters, 1) the un-transformed Model, 2) the transformed Model, and 3) the GaugeGroupElement representing the transformation that brings the first argument into the second.

  • comm (mpi4py.MPI.Comm, optional) – When not None, an MPI communicator for distributing the computation across multiple processors.

  • verbosity (int, optional) – How much detail to send to stdout.

Returns

  • model if return_all == False

  • (goodnessMin, gaugeMx, model) if return_all == True – where goodnessMin is the minimum value of the goodness function (the best ‘goodness’) found, gaugeMx is the gauge matrix used to transform the model, and model is the final gauge-transformed model.

pygsti.algorithms.gaugeopt._create_objective_fn(model, target_model, item_weights=None, cptp_penalty_factor=0, spam_penalty_factor=0, gates_metric='frobenius', spam_metric='frobenius', method=None, comm=None, check_jac=False)

Creates the objective function and jacobian (if available) for gaugeopt_to_target

pygsti.algorithms.gaugeopt._cptp_penalty_size(mdl)

Helper function - same as that in core.py.

pygsti.algorithms.gaugeopt._spam_penalty_size(mdl)

Helper function - same as that in core.py.

pygsti.algorithms.gaugeopt._cptp_penalty(mdl, prefactor, op_basis)

Helper function - CPTP penalty: (sum of tracenorms of gates), which in least squares optimization means returning an array of the sqrt(tracenorm) of each gate. This function is the same as that in core.py.

Returns

numpy array – a (real) 1D array of length len(mdl.operations).

pygsti.algorithms.gaugeopt._spam_penalty(mdl, prefactor, op_basis)

Helper function - CPTP penalty: (sum of tracenorms of gates), which in least squares optimization means returning an array of the sqrt(tracenorm) of each gate. This function is the same as that in core.py.

Returns

numpy array – a (real) 1D array of length _spam_penalty_size(mdl)

pygsti.algorithms.gaugeopt._cptp_penalty_jac_fill(cp_penalty_vec_grad_to_fill, mdl_pre, mdl_post, gauge_group_el, prefactor, op_basis, wrt_filter)

Helper function - jacobian of CPTP penalty (sum of tracenorms of gates) Returns a (real) array of shape (len(mdl.operations), gauge_group_el.num_params).

pygsti.algorithms.gaugeopt._spam_penalty_jac_fill(spam_penalty_vec_grad_to_fill, mdl_pre, mdl_post, gauge_group_el, prefactor, op_basis, wrt_filter)

Helper function - jacobian of CPTP penalty (sum of tracenorms of gates) Returns a (real) array of shape (_spam_penalty_size(mdl), gauge_group_el.num_params).