From 01c7fd9ab306b41714d563989c7b8911877f7edf Mon Sep 17 00:00:00 2001 From: Riley Murray Date: Mon, 21 Oct 2024 11:05:09 -0400 Subject: [PATCH] migrate internal code from CustomLMOptimizer to SimplerLMOptimizer and custom_leastsq to simplish_leastsq --- pygsti/algorithms/core.py | 9 ++++----- pygsti/algorithms/gaugeopt.py | 4 ++-- pygsti/optimize/simplerlm.py | 15 ++++++++++++++- pygsti/protocols/gst.py | 13 +++++++------ .../{test_customlm.py => test_simplerlm.py} | 16 ++++++++-------- test/unit/protocols/test_gst.py | 4 ++-- 6 files changed, 37 insertions(+), 24 deletions(-) rename test/unit/optimize/{test_customlm.py => test_simplerlm.py} (65%) diff --git a/pygsti/algorithms/core.py b/pygsti/algorithms/core.py index dd0a21ef7..6e0628264 100644 --- a/pygsti/algorithms/core.py +++ b/pygsti/algorithms/core.py @@ -31,8 +31,7 @@ from pygsti.modelmembers import states as _state from pygsti.circuits.circuitlist import CircuitList as _CircuitList from pygsti.baseobjs.resourceallocation import ResourceAllocation as _ResourceAllocation -from pygsti.optimize.customlm import CustomLMOptimizer as _CustomLMOptimizer -from pygsti.optimize.customlm import Optimizer as _Optimizer +from pygsti.optimize.simplerlm import Optimizer as _Optimizer, SimplerLMOptimizer as _SimplerLMOptimizer from pygsti import forwardsims as _fwdsims from pygsti import layouts as _layouts @@ -619,7 +618,7 @@ def run_gst_fit_simple(dataset, start_model, circuits, optimizer, objective_func model : Model the best-fit model. """ - optimizer = optimizer if isinstance(optimizer, _Optimizer) else _CustomLMOptimizer.cast(optimizer) + optimizer = optimizer if isinstance(optimizer, _Optimizer) else _SimplerLMOptimizer.cast(optimizer) objective_function_builder = _objfns.ObjectiveFunctionBuilder.cast(objective_function_builder) array_types = optimizer.array_types + \ objective_function_builder.compute_array_types(optimizer.called_objective_methods, start_model.sim) @@ -666,7 +665,7 @@ def run_gst_fit(mdc_store, optimizer, objective_function_builder, verbosity=0): objfn_store : MDCObjectiveFunction the objective function and store containing the best-fit model evaluated at the best-fit point. """ - optimizer = optimizer if isinstance(optimizer, _Optimizer) else _CustomLMOptimizer.cast(optimizer) + optimizer = optimizer if isinstance(optimizer, _Optimizer) else _SimplerLMOptimizer.cast(optimizer) comm = mdc_store.resource_alloc.comm profiler = mdc_store.resource_alloc.profiler printer = VerbosityPrinter.create_printer(verbosity, comm) @@ -843,7 +842,7 @@ def iterative_gst_generator(dataset, start_model, circuit_lists, (an "evaluated" model-dataset-circuits store). """ resource_alloc = _ResourceAllocation.cast(resource_alloc) - optimizer = optimizer if isinstance(optimizer, _Optimizer) else _CustomLMOptimizer.cast(optimizer) + optimizer = optimizer if isinstance(optimizer, _Optimizer) else _SimplerLMOptimizer.cast(optimizer) comm = resource_alloc.comm profiler = resource_alloc.profiler printer = VerbosityPrinter.create_printer(verbosity, comm) diff --git a/pygsti/algorithms/gaugeopt.py b/pygsti/algorithms/gaugeopt.py index fcd52d267..6b341062a 100644 --- a/pygsti/algorithms/gaugeopt.py +++ b/pygsti/algorithms/gaugeopt.py @@ -290,7 +290,7 @@ def gaugeopt_custom(model, objective_fn, gauge_group=None, gaugeGroupEl = gauge_group.compute_element(x0) # re-used element for evals def _call_objective_fn(gauge_group_el_vec, oob_check=False): - # Note: oob_check can be True if oob_check_interval>=1 is given to the custom_leastsq below + # Note: oob_check can be True if oob_check_interval>=1 is given to the simplish_leastsq below gaugeGroupEl.from_vector(gauge_group_el_vec) return objective_fn(gaugeGroupEl, oob_check) @@ -309,7 +309,7 @@ def _call_jacobian_fn(gauge_group_el_vec): assert(_call_jacobian_fn is not None), "Cannot use 'ls' method unless jacobian is available" ralloc = _baseobjs.ResourceAllocation(comm) # FUTURE: plumb up a resource alloc object? test_f = _call_objective_fn(x0) - solnX, converged, msg, _, _, _, _, _ = _opt.custom_leastsq( + solnX, converged, msg, _, _, _, _, _ = _opt.simplish_leastsq( _call_objective_fn, _call_jacobian_fn, x0, f_norm2_tol=tol, jac_norm_tol=tol, rel_ftol=tol, rel_xtol=tol, max_iter=maxiter, resource_alloc=ralloc, diff --git a/pygsti/optimize/simplerlm.py b/pygsti/optimize/simplerlm.py index ed966c57f..11a3cd8b3 100644 --- a/pygsti/optimize/simplerlm.py +++ b/pygsti/optimize/simplerlm.py @@ -183,6 +183,19 @@ class SimplerLMOptimizer(Optimizer): "per-circuit quantities" computed by the objective function's `.percircuit()` and `.lsvec_percircuit()` methods (`'percircuit'` mode). """ + + @classmethod + def cast(cls, obj): + if isinstance(obj, cls): + return obj + if obj: + try: + return cls(**obj) + except: + from pygsti.optimize.customlm import CustomLMOptimizer + return CustomLMOptimizer(**obj) + return cls() + def __init__(self, maxiter=100, maxfev=100, tol=1e-6, fditer=0, first_fditer=0, uphill_step_threshold=0.0, init_munu="auto", oob_check_interval=0, oob_action="reject", oob_check_mode=0, serial_solve_proc_threshold=100, lsvec_mode="normal"): @@ -624,7 +637,7 @@ def simplish_leastsq( if _np.allclose(x, best_x): rawJTJ_scratch[:, :] = JTJ[:, :] # use pre-allocated memory rawJTJ_scratch[idiag] = undamped_JTJ_diag # no damping; the "raw" JTJ - best_x_state = best_x_state[0:5] + (rawJTJ_scratch,) # update mu,nu,JTJ of initial "best state" + best_x_state = best_x_state[0:4] + (rawJTJ_scratch,) # update mu,nu,JTJ of initial "best state" #determing increment using adaptive damping while True: # inner loop diff --git a/pygsti/protocols/gst.py b/pygsti/protocols/gst.py index 9255943d3..a0695dcba 100644 --- a/pygsti/protocols/gst.py +++ b/pygsti/protocols/gst.py @@ -1249,15 +1249,16 @@ def __init__(self, initial_model=None, gaugeopt_suite='stdgaugeopt', if isinstance(optimizer, _opt.Optimizer): self.optimizer = optimizer - if isinstance(optimizer, _opt.CustomLMOptimizer) and optimizer.first_fditer is None: - #special behavior: can set optimizer's first_fditer to `None` to mean "fill with default" + if hasattr(optimizer,'first_fditer') and optimizer.first_fditer is None: + # special behavior: can set optimizer's first_fditer to `None` to mean "fill with default" self.optimizer = _copy.deepcopy(optimizer) # don't mess with caller's optimizer self.optimizer.first_fditer = default_first_fditer else: - if optimizer is None: optimizer = {} + if optimizer is None: + optimizer = {} if 'first_fditer' not in optimizer: # then add default first_fditer value optimizer['first_fditer'] = default_first_fditer - self.optimizer = _opt.CustomLMOptimizer.cast(optimizer) + self.optimizer = _opt.SimplerLMOptimizer.cast(optimizer) self.objfn_builders = GSTObjFnBuilders.cast(objfn_builders) @@ -1751,12 +1752,12 @@ def __init__(self, modes=('full TP','CPTPLND','Target'), gaugeopt_suite='stdgaug self.target_model = target_model self.gaugeopt_suite = GSTGaugeOptSuite.cast(gaugeopt_suite) self.objfn_builders = GSTObjFnBuilders.cast(objfn_builders) if (objfn_builders is not None) else None - self.optimizer = _opt.CustomLMOptimizer.cast(optimizer) + self.optimizer = _opt.SimplerLMOptimizer.cast(optimizer) self.badfit_options = GSTBadFitOptions.cast(badfit_options) self.verbosity = verbosity if not isinstance(optimizer, _opt.Optimizer) and isinstance(optimizer, dict) \ - and 'first_fditer' not in optimizer: # then a dict was cast into a CustomLMOptimizer above. + and 'first_fditer' not in optimizer: # then a dict was cast into an Optimizer above. # by default, set special "first_fditer=auto" behavior (see logic in GateSetTomography.__init__) self.optimizer.first_fditer = None diff --git a/test/unit/optimize/test_customlm.py b/test/unit/optimize/test_simplerlm.py similarity index 65% rename from test/unit/optimize/test_customlm.py rename to test/unit/optimize/test_simplerlm.py index 699196694..caafed4ed 100644 --- a/test/unit/optimize/test_customlm.py +++ b/test/unit/optimize/test_simplerlm.py @@ -1,7 +1,7 @@ import numpy as np from pygsti.optimize import arraysinterface as _ari -from pygsti.optimize import customlm as lm +from pygsti.optimize import simplerlm as lm from ..util import BaseCase @@ -19,25 +19,25 @@ def gjac(x): return np.array([2 * x[0]],'d') -class CustomLMTester(BaseCase): - def test_custom_leastsq_infinite_objective_fn_norm_at_x0(self): +class LMTester(BaseCase): + def test_simplish_leastsq_infinite_objective_fn_norm_at_x0(self): x0 = np.ones(3, 'd') ari = _ari.UndistributedArraysInterface(2, 3) - xf, converged, msg, *_ = lm.custom_leastsq(f, jac, x0, arrays_interface=ari) + xf, converged, msg, *_ = lm.simplish_leastsq(f, jac, x0, arrays_interface=ari) self.assertEqual(msg, "Infinite norm of objective function at initial point!") - def test_custom_leastsq_max_iterations_exceeded(self): + def test_simplish_leastsq_max_iterations_exceeded(self): x0 = np.ones(3, 'd') ari = _ari.UndistributedArraysInterface(2, 3) - xf, converged, msg, *_ = lm.custom_leastsq(f, jac, x0, max_iter=0, arrays_interface=ari) + xf, converged, msg, *_ = lm.simplish_leastsq(f, jac, x0, max_iter=0, arrays_interface=ari) self.assertEqual(msg, "Maximum iterations (0) exceeded") - def test_custom_leastsq_x_limits(self): + def test_simplish_leastsq_x_limits(self): #perform optimization of g(x), which has minimum at 0, using limits so x must be 10 > x > 1 # and check that the optimal x is near 1.0: x0 = np.array([6.0],'d') xlimits = np.array([[1.0, 10.0]], 'd') ari = _ari.UndistributedArraysInterface(1, 1) - xf, converged, msg, *_ = lm.custom_leastsq(g, gjac, x0, max_iter=100, arrays_interface=ari, + xf, converged, msg, *_ = lm.simplish_leastsq(g, gjac, x0, max_iter=100, arrays_interface=ari, x_limits=xlimits) self.assertAlmostEqual(xf[0], 1.0) diff --git a/test/unit/protocols/test_gst.py b/test/unit/protocols/test_gst.py index 8a0b03b82..2b712f868 100644 --- a/test/unit/protocols/test_gst.py +++ b/test/unit/protocols/test_gst.py @@ -5,7 +5,7 @@ from pygsti.objectivefns.objectivefns import PoissonPicDeltaLogLFunction from pygsti.models.gaugegroup import TrivialGaugeGroup from pygsti.objectivefns import FreqWeightedChi2Function -from pygsti.optimize.customlm import CustomLMOptimizer +from pygsti.optimize.simplerlm import SimplerLMOptimizer from pygsti.protocols import gst from pygsti.protocols.estimate import Estimate from pygsti.protocols.protocol import ProtocolData, Protocol @@ -65,7 +65,7 @@ def test_gaugeopt_suite_raises_on_bad_suite(self): def test_add_badfit_estimates(self): builder = PoissonPicDeltaLogLFunction.builder() - opt = CustomLMOptimizer() + opt = SimplerLMOptimizer() badfit_opts = gst.GSTBadFitOptions(threshold=-10, actions=("robust", "Robust", "robust+", "Robust+", "wildcard", "do nothing")) res = self.results.copy()