:py:mod:`pygsti.optimize.simplerlm` =================================== .. py:module:: pygsti.optimize.simplerlm .. autoapi-nested-parse:: Custom implementation of the Levenberg-Marquardt Algorithm (but simpler than customlm.py) Module Contents --------------- Classes ~~~~~~~ .. autoapisummary:: pygsti.optimize.simplerlm.OptimizerResult pygsti.optimize.simplerlm.Optimizer pygsti.optimize.simplerlm.SimplerLMOptimizer Functions ~~~~~~~~~ .. autoapisummary:: pygsti.optimize.simplerlm.damp_coeff_update pygsti.optimize.simplerlm.jac_guarded pygsti.optimize.simplerlm.simplish_leastsq .. py:class:: OptimizerResult(objective_func, opt_x, opt_f=None, opt_jtj=None, opt_unpenalized_f=None, chi2_k_distributed_qty=None, optimizer_specific_qtys=None) Bases: :py:obj:`object` The result from an optimization. Parameters ---------- objective_func : ObjectiveFunction The objective function that was optimized. opt_x : numpy.ndarray The optimal argument (x) value. Often a vector of parameters. opt_f : numpy.ndarray the optimal objective function (f) value. Often this is the least-squares vector of objective function values. opt_jtj : numpy.ndarray, optional the optimial `dot(transpose(J),J)` value, where `J` is the Jacobian matrix. This may be useful for computing approximate error bars. opt_unpenalized_f : numpy.ndarray, optional the optimal objective function (f) value with any penalty terms removed. chi2_k_distributed_qty : float, optional a value that is supposed to be chi2_k distributed. optimizer_specific_qtys : dict, optional a dictionary of additional optimization parameters. .. py:attribute:: objective_func .. py:attribute:: x .. py:attribute:: f :value: 'None' .. py:attribute:: jtj :value: 'None' .. py:attribute:: f_no_penalties :value: 'None' .. py:attribute:: optimizer_specific_qtys :value: 'None' .. py:attribute:: chi2_k_distributed_qty :value: 'None' .. py:class:: Optimizer Bases: :py:obj:`pygsti.baseobjs.nicelyserializable.NicelySerializable` An optimizer. Optimizes an objective function. .. py:method:: cast(obj) :classmethod: Cast `obj` to a :class:`Optimizer`. If `obj` is already an `Optimizer` it is just returned, otherwise this function tries to create a new object using `obj` as a dictionary of constructor arguments. Parameters ---------- obj : Optimizer or dict The object to cast. Returns ------- Optimizer .. py:class:: SimplerLMOptimizer(maxiter=100, maxfev=100, tol=1e-06, fditer=0, first_fditer=0, init_munu='auto', oob_check_interval=0, oob_action='reject', oob_check_mode=0, serial_solve_proc_threshold=100, lsvec_mode='normal') Bases: :py:obj:`Optimizer` A Levenberg-Marquardt optimizer customized for GST-like problems. Parameters ---------- maxiter : int, optional The maximum number of (outer) interations. maxfev : int, optional The maximum function evaluations. tol : float or dict, optional The tolerance, specified as a single float or as a dict with keys `{'relx', 'relf', 'jac', 'maxdx'}`. A single float sets the `'relf'` and `'jac'` elemments and leaves the others at their default values. fditer : int optional Internally compute the Jacobian using a finite-difference method for the first `fditer` iterations. This is useful when the initial point lies at a special or singular point where the analytic Jacobian is misleading. first_fditer : int, optional Number of finite-difference iterations applied to the first stage of the optimization (only). Unused. init_munu : tuple, optional If not None, a (mu, nu) tuple of 2 floats giving the initial values for mu and nu. oob_check_interval : int, optional Every `oob_check_interval` outer iterations, the objective function (`obj_fn`) is called with a second argument 'oob_check', set to True. In this case, `obj_fn` can raise a ValueError exception to indicate that it is Out Of Bounds. If `oob_check_interval` is 0 then this check is never performed; if 1 then it is always performed. oob_action : {"reject","stop"} What to do when the objective function indicates (by raising a ValueError as described above). `"reject"` means the step is rejected but the optimization proceeds; `"stop"` means the optimization stops and returns as converged at the last known-in-bounds point. oob_check_mode : int, optional An advanced option, expert use only. If 0 then the optimization is halted as soon as an *attempt* is made to evaluate the function out of bounds. If 1 then the optimization is halted only when a would-be *accepted* step is out of bounds. serial_solve_proc_threshold : int, optional When there are fewer than this many processors, the optimizer will solve linear systems serially, using SciPy on a single processor, rather than using a parallelized Gaussian Elimination (with partial pivoting) algorithm coded in Python. Since SciPy's implementation is more efficient, it's not worth using the parallel version until there are many processors to spread the work among. lsvec_mode : {'normal', 'percircuit'} Whether the terms used in the least-squares optimization are the "elements" as computed by the objective function's `.terms()` and `.lsvec()` methods (`'normal'` mode) or the "per-circuit quantities" computed by the objective function's `.percircuit()` and `.lsvec_percircuit()` methods (`'percircuit'` mode). .. py:attribute:: maxiter :value: '100' .. py:attribute:: maxfev :value: '100' .. py:attribute:: tol :value: '1e-06' .. py:attribute:: fditer :value: '0' .. py:attribute:: first_fditer :value: '0' .. py:attribute:: init_munu :value: "'auto'" .. py:attribute:: oob_check_interval :value: '0' .. py:attribute:: oob_action :value: "'reject'" .. py:attribute:: oob_check_mode :value: '0' .. py:attribute:: array_types :value: "('p', 'p', 'p', 'e', 'ep')" .. py:attribute:: called_objective_methods :value: "('lsvec', 'dlsvec')" .. py:attribute:: serial_solve_proc_threshold :value: '100' .. py:attribute:: lsvec_mode :value: "'normal'" .. py:method:: cast(obj) :classmethod: Cast `obj` to a :class:`Optimizer`. If `obj` is already an `Optimizer` it is just returned, otherwise this function tries to create a new object using `obj` as a dictionary of constructor arguments. Parameters ---------- obj : Optimizer or dict The object to cast. Returns ------- Optimizer .. py:method:: run(objective: pygsti.objectivefns.objectivefns.TimeIndependentMDCObjectiveFunction, profiler, printer) Perform the optimization. Parameters ---------- objective : ObjectiveFunction The objective function to optimize. profiler : Profiler A profiler to track resource usage. printer : VerbosityPrinter printer to use for sending output to stdout. .. py:function:: damp_coeff_update(mu, nu, half_max_nu, reject_msg, printer) .. py:function:: jac_guarded(k: int, num_fd_iters: int, obj_fn: Callable, jac_fn: Callable, f, ari, global_x, fdJac_work) .. py:function:: simplish_leastsq(obj_fn, jac_fn, x0, f_norm2_tol=1e-06, jac_norm_tol=1e-06, rel_ftol=1e-06, rel_xtol=1e-06, max_iter=100, num_fd_iters=0, max_dx_scale=1.0, init_munu='auto', oob_check_interval=0, oob_action='reject', oob_check_mode=0, resource_alloc=None, arrays_interface=None, serial_solve_proc_threshold=100, x_limits=None, verbosity=0, profiler=None) An implementation of the Levenberg-Marquardt least-squares optimization algorithm customized for use within pyGSTi. This general purpose routine mimic to a large extent the interface used by `scipy.optimize.leastsq`, though it implements a newer (and more robust) version of the algorithm. Parameters ---------- obj_fn : function The objective function. Must accept and return 1D numpy ndarrays of length N and M respectively. Same form as scipy.optimize.leastsq. jac_fn : function The jacobian function (not optional!). Accepts a 1D array of length N and returns an array of shape (M,N). x0 : numpy.ndarray Initial evaluation point. f_norm2_tol : float, optional Tolerace for `F^2` where `F = `norm( sum(obj_fn(x)**2) )` is the least-squares residual. If `F**2 < f_norm2_tol`, then mark converged. jac_norm_tol : float, optional Tolerance for jacobian norm, namely if `infn(dot(J.T,f)) < jac_norm_tol` then mark converged, where `infn` is the infinity-norm and `f = obj_fn(x)`. rel_ftol : float, optional Tolerance on the relative reduction in `F^2`, that is, if `d(F^2)/F^2 < rel_ftol` then mark converged. rel_xtol : float, optional Tolerance on the relative value of `|x|`, so that if `d(|x|)/|x| < rel_xtol` then mark converged. max_iter : int, optional The maximum number of (outer) interations. num_fd_iters : int optional Internally compute the Jacobian using a finite-difference method for the first `num_fd_iters` iterations. This is useful when `x0` lies at a special or singular point where the analytic Jacobian is misleading. max_dx_scale : float, optional If not None, impose a limit on the magnitude of the step, so that `|dx|^2 < max_dx_scale^2 * len(dx)` (so elements of `dx` should be, roughly, less than `max_dx_scale`). init_munu : tuple, optional If not None, a (mu, nu) tuple of 2 floats giving the initial values for mu and nu. oob_check_interval : int, optional Every `oob_check_interval` outer iterations, the objective function (`obj_fn`) is called with a second argument 'oob_check', set to True. In this case, `obj_fn` can raise a ValueError exception to indicate that it is Out Of Bounds. If `oob_check_interval` is 0 then this check is never performed; if 1 then it is always performed. oob_action : {"reject","stop"} What to do when the objective function indicates (by raising a ValueError as described above). `"reject"` means the step is rejected but the optimization proceeds; `"stop"` means the optimization stops and returns as converged at the last known-in-bounds point. oob_check_mode : int, optional An advanced option, expert use only. If 0 then the optimization is halted as soon as an *attempt* is made to evaluate the function out of bounds. If 1 then the optimization is halted only when a would-be *accepted* step is out of bounds. resource_alloc : ResourceAllocation, optional When not None, an resource allocation object used for distributing the computation across multiple processors. arrays_interface : ArraysInterface An object that provides an interface for creating and manipulating data arrays. serial_solve_proc_threshold : int optional When there are fewer than this many processors, the optimizer will solve linear systems serially, using SciPy on a single processor, rather than using a parallelized Gaussian Elimination (with partial pivoting) algorithm coded in Python. Since SciPy's implementation is more efficient, it's not worth using the parallel version until there are many processors to spread the work among. x_limits : numpy.ndarray, optional A (num_params, 2)-shaped array, holding on each row the (min, max) values for the corresponding parameter (element of the "x" vector). If `None`, then no limits are imposed. verbosity : int, optional Amount of detail to print to stdout. profiler : Profiler, optional A profiler object used for to track timing and memory usage. Returns ------- x : numpy.ndarray The optimal solution. converged : bool Whether the solution converged. msg : str A message indicating why the solution converged (or didn't).