diff --git a/docs/adr/0004-cl-min-for-multi-point-asks.md b/docs/adr/0004-cl-min-for-multi-point-asks.md new file mode 100644 index 0000000..a64c245 --- /dev/null +++ b/docs/adr/0004-cl-min-for-multi-point-asks.md @@ -0,0 +1,90 @@ +# 0004. Use constant-liar (cl_min) for all multi-point asks + +- **Status:** accepted +- **Date:** 2026-06-07 +- **Deciders:** Jakob Langdal + +## Context + +`POST /optimizer` with `extras.experimentSuggestionCount > 1` asks +ProcessOptimizer for a batch of next experiments via `optimizer.ask(n_points=N)`. +Until now we only passed an explicit strategy on the *constrained* branch +(`strategy="cl_min"`, added because Steinerberger sampling raises with +constraints — see issue #83); the *unconstrained* branch fell through to +ProcessOptimizer's default, `strategy="stbr_fill"`. + +For `N > 1` on a model that is already fitted (`#data ≥ initialPoints`), +`stbr_fill` routes points 2…N through `stbr_scipy()` — a Steinerberger +space-filling solver that runs 20 SciPy `minimize` restarts over the +**one-hot-transformed** space, re-transforming every observation on each +objective evaluation. The cost is dominated by the one-hot expansion of +categorical variables. On a real experiment with 5 continuous + 5 categorical +variables (31 transformed dimensions) a single count-2 request consumed ~30 +minutes of CPU — effectively a hang synchronously, or a `WORKER_TIMEOUT` failure +behind the worker. Low-dimensional bundled examples (Catapult: 4 dims, Brownie: +7) masked it because the same path is cheap there. + +The first suggestion is unaffected (it is always `optimizer._ask()`), and the +LHS-based initialization phase is unaffected; only the 2nd…Nth points of a +fitted-model batch took the slow path. + +## Decision + +`cl_min` (constant liar) is the **default** strategy for multi-point asks in +`optimizer._compute_next_experiments`. For `n_points == 1` the strategy is +irrelevant (ProcessOptimizer returns `_ask()` directly), so single-suggestion +behaviour is byte-for-byte unchanged. Constant-liar optimizes the acquisition +function through the GP — which handles categorical dimensions natively — +instead of solving a high-dimensional Steinerberger problem, so it stays fast +(sub-second to a few seconds) regardless of categorical cardinality. + +We still **opt into `stbr_fill`** (Steinerberger space-filling, the nicer +exploration spread for extra batch points) when it is affordable: the request is +unconstrained, `n_points > 1`, and a deterministic cost estimate is within a +wall-clock budget (`STBR_TIME_BUDGET_SECONDS`, default 10s). The estimate +(`_estimate_stbr_seconds`) is calibrated from benchmarks: cost is dominated +*super-linearly* by categorical one-hot dimensions, mildly by continuous/discrete +dimensions, and ~linearly by `n_points - 1`. It is a pure function of the space +and batch size — **not** a wall-clock measurement — so the chosen strategy, and +therefore the suggestions, stay reproducible across machines and load. Anything +over budget (e.g. the motivating 31-dimension experiment) falls back to `cl_min`. + +## Consequences + +- **Performance:** the count-2 hang is gone (~30 min → ~3 s on the motivating + experiment); cost now scales gently with batch size rather than exploding with + categorical dimensionality. +- **Behaviour change (single objective, count > 1, fitted model, no + constraint):** the extra suggestions change character — from space-filling + *exploration* (Steinerberger) to acquisition-driven *batch BO* with the + constant-liar penalty. This is the conventional batch-BO behaviour and matches + what the constrained path always did, but it is a real change for anyone who + relied on the exploration spread. Count 1 is unchanged. +- **Unblocks constrained batches:** constrained multi-point asks via `cl_min` + are fast and respect the constraint (verified: every suggested point satisfies + the `SumEquals`). The frontend previously capped constrained experiments to a + single suggestion because constrained multi-point was unsupported (it errored + on the Steinerberger default); that cap can now be lifted. +- Steinerberger space-filling is retained for cheap experiments via the + server-side budget gate, so low-dimensional experiments keep the nicer + exploration spread while categorical-heavy ones stay fast on `cl_min`. The + cost is a calibrated heuristic (`_estimate_stbr_seconds`) that must be + re-checked if the ProcessOptimizer/SciPy stack or reference hardware changes + materially. + +## Alternatives considered + +- **Always `cl_min`, drop Steinerberger entirely.** Simplest, but loses the + exploration spread for the small experiments where it is both nice and cheap. + Rejected in favour of the budget gate. +- **A wall-clock time budget** (run-then-timeout, or probe-one-restart and + extrapolate). Adapts to hardware, but makes the chosen strategy — and thus the + suggestions — depend on CPU speed and load, breaking reproducibility for a + seeded optimizer. Rejected in favour of a deterministic structural estimate. +- **Gate on raw `transformed_n_dims`.** Rejected: continuous dimensions are + cheap (20 continuous dims ≈ 3s) while categorical one-hot dimensions are not + (10 one-hot ≈ 17s), so a single transformed-dimension threshold mis-predicts. + The estimate weights categorical load specifically. +- **Expose `strategy` as a request/`extras` option** (default `cl_min`). Pushes + a performance-critical, easily-misused knob to clients for a capability almost + no caller needs; rejected in favour of a safe server-side default. diff --git a/optimizerapi/optimizer.py b/optimizerapi/optimizer.py index 088054b..fbdac0a 100644 --- a/optimizerapi/optimizer.py +++ b/optimizerapi/optimizer.py @@ -17,7 +17,7 @@ import matplotlib.pyplot as plt import numpy from ProcessOptimizer import Optimizer, expected_minimum -from ProcessOptimizer.space import Real +from ProcessOptimizer.space import Categorical, Real from ProcessOptimizer.space.constraints import SumEquals from typing import TYPE_CHECKING, Any, cast @@ -87,9 +87,58 @@ def _parse_extras(extras: "Extras", logger: "logging.Logger") -> _ParsedExtras: ) +# Wall-clock budget for the Steinerberger (stbr_fill) space-filling path. We opt +# into stbr_fill only when its estimated runtime fits this budget; otherwise we +# use the always-cheap constant-liar (cl_min). Configurable for slower/faster +# deployments. See ADR 0004. +_STBR_TIME_BUDGET_SECONDS = float(os.getenv("STBR_TIME_BUDGET_SECONDS", "10")) + + +def _estimate_stbr_seconds(space: Any, n_points: int) -> float: + """Estimate the wall-clock cost of the Steinerberger path for this batch. + + For ``n_points > 1`` on a fitted model, ProcessOptimizer's ``stbr_fill`` + strategy computes each extra point with ``stbr_scipy()`` — 20 SciPy + minimisations over the one-hot-transformed space. Benchmarking (CPython + venv, the reference deployment — NOT Pyodide) shows the cost is dominated + *super-linearly* by the categorical one-hot dimensions, only mildly by + continuous/discrete dimensions, and ~linearly by the number of extra points + (``n_points - 1``). The coefficients below are tuned to upper-bound the + measured runtimes near the budget knee, so the estimate is conservative + (it prefers the safe cl_min fallback when in doubt). It is a pure function + of the space and batch size, which keeps the chosen strategy — and thus the + suggestions — independent of machine speed/load (important for + reproducibility). + """ + onehot = sum( + len(d.categories) for d in space.dimensions if isinstance(d, Categorical) + ) + n_other = sum(1 for d in space.dimensions if not isinstance(d, Categorical)) + per_point = 0.5 + 0.15 * n_other + 0.2 * onehot**2 + return max(0, n_points - 1) * per_point + + +def _choose_ask_strategy( + n_points: int, has_constraints: bool, stbr_estimate_seconds: float +) -> str: + """Pick the ``ask`` strategy. + + Constant-liar (``cl_min``) is the safe default: fast and categorical-safe. + We opt into Steinerberger (``stbr_fill``) space-filling — which gives a nicer + exploration spread for the extra batch points — only when the request is + unconstrained, actually a batch, and estimated to finish within budget. + (``stbr_fill`` raises with constraints, and the strategy is irrelevant for a + single point.) + """ + if n_points <= 1 or has_constraints: + return "cl_min" + if stbr_estimate_seconds <= _STBR_TIME_BUDGET_SECONDS: + return "stbr_fill" + return "cl_min" + + def _compute_next_experiments( optimizer: Any, - cfg: "OptimizerConfig", n_points: int, ) -> list[list[str | float]]: """Ask the optimizer for the next N experiments, normalising the shape. @@ -97,11 +146,17 @@ def _compute_next_experiments( ``optimizer.ask`` can return either a single experiment (flat list) or a list of experiments. We always return a list of lists. """ - constraints = cfg.get("constraints", []) - if constraints: - next_exp = optimizer.ask(n_points=n_points, strategy="cl_min") - else: - next_exp = optimizer.ask(n_points=n_points) + has_constraints = optimizer.get_constraints() is not None + strategy = "cl_min" + if n_points > 1 and not has_constraints: + est = _estimate_stbr_seconds(optimizer.space, n_points) + strategy = _choose_ask_strategy(n_points, has_constraints, est) + logging.getLogger(__name__).info( + "multi-point ask n_points=%d -> strategy=%s " + "(stbr est %.1fs, budget %.0fs)", + n_points, strategy, est, _STBR_TIME_BUDGET_SECONDS, + ) + next_exp = optimizer.ask(n_points=n_points, strategy=strategy) if next_exp and not any(isinstance(x, list) for x in next_exp): next_exp = [next_exp] return cast(list[list[str | float]], round_to_length_scales(next_exp, optimizer.space)) @@ -314,7 +369,7 @@ def process_result( parsed = _parse_extras(extras, logging.getLogger(__name__)) result_details["next"] = _compute_next_experiments( - optimizer, cfg, parsed.experiment_suggestion_count + optimizer, parsed.experiment_suggestion_count ) if len(data) >= cfg["initialPoints"]: diff --git a/tests/test_optimizer.py b/tests/test_optimizer.py index 75a2425..61bb062 100644 --- a/tests/test_optimizer.py +++ b/tests/test_optimizer.py @@ -6,7 +6,14 @@ import copy import collections.abc import json +import numpy as np +from ProcessOptimizer import Optimizer from optimizerapi import optimizer_handler as optimizer +from optimizerapi.optimizer import ( + _choose_ask_strategy, + _compute_next_experiments, + _estimate_stbr_seconds, +) from optimizerapi.securepickle import get_crypto, pickleToString, unpickleFromString sampleData = [ @@ -336,12 +343,114 @@ def test_when_using_constraints_strategy_cl_min_should_be_used(mock): instance.ask.assert_called_once_with(n_points=3, strategy="cl_min") -@patch("optimizerapi.optimizer.Optimizer") -def test_when_not_using_constraints_standard_strategy_should_be_used(mock): - instance = mock.return_value - request = brownie_without_constraints - optimizer.run_optimizer(body=request) - instance.ask.assert_called_once_with(n_points=3) +def test_choose_ask_strategy_prefers_cl_min_unless_stbr_affordable(): + # single point: strategy is irrelevant, use the safe default + assert _choose_ask_strategy(1, False, 0.0) == "cl_min" + # constraints: stbr_fill raises with constraints + assert _choose_ask_strategy(3, True, 0.1) == "cl_min" + # unconstrained batch, cheap enough -> opt into Steinerberger + assert _choose_ask_strategy(2, False, 5.0) == "stbr_fill" + # unconstrained batch, over budget -> fall back to cl_min + assert _choose_ask_strategy(2, False, 1000.0) == "cl_min" + + +def test_estimate_stbr_seconds_dominated_by_categorical_load(): + cont = Optimizer([(0.0, 5.0)] * 4, "GP", n_objectives=1).space + cat_heavy = Optimizer( + [(0.0, 5.0)] * 5 + [tuple("abcde")] * 5, "GP", n_objectives=1 + ).space + cont_20 = Optimizer([(0.0, 5.0)] * 20, "GP", n_objectives=1).space + # a small continuous space is affordable; a categorical-heavy one is not + assert _estimate_stbr_seconds(cont, 2) <= 10 + assert _estimate_stbr_seconds(cat_heavy, 2) > 10 + # categorical one-hot load dwarfs an equivalent count of continuous dims + assert _estimate_stbr_seconds(cat_heavy, 2) > _estimate_stbr_seconds(cont_20, 2) + + +def test_cheap_unconstrained_batch_uses_steinerberger(): + # 4 continuous dims, fitted model, count 2: cheap -> stbr_fill, returns 2 pts. + space = [(0.0, 100.0)] * 4 + opt = Optimizer( + space, + "GP", + n_initial_points=4, + acq_func="EI", + acq_func_kwargs={"kappa": 1.96, "xi": 0.01}, + n_objectives=1, + ) + rng = np.random.RandomState(3) + opt.tell( + rng.uniform(0, 100, size=(8, 4)).tolist(), + rng.uniform(0, 1, size=8).tolist(), + ) + assert ( + _choose_ask_strategy(2, False, _estimate_stbr_seconds(opt.space, 2)) + == "stbr_fill" + ) + nxt = _compute_next_experiments(opt, 2) + assert len(nxt) == 2 + assert all(len(point) == len(space) for point in nxt) + + +def test_multi_suggestion_without_constraints_terminates_quickly(): + # Regression for the stbr_fill hang: a fitted model (data >= initialPoints) + # with a categorical-heavy space and experimentSuggestionCount > 1 and no + # constraints. Under the old default strategy this ran for tens of minutes; + # with cl_min it returns in about a second. Guard with a watchdog so a + # regression fails fast instead of hanging the suite. + import signal + + space = [ + {"type": "continuous", "name": "a", "from": 0, "to": 5}, + {"type": "continuous", "name": "b", "from": 0, "to": 5}, + {"type": "continuous", "name": "c", "from": 0, "to": 5}, + {"type": "category", "name": "T", "categories": ["95", "105", "115", "125"]}, + {"type": "category", "name": "pH", "categories": ["4", "5", "6", "7"]}, + {"type": "category", "name": "t", "categories": ["15", "25", "35", "45"]}, + ] + # initialPoints == 4, supply 5 data points so the model is fitted and + # _n_initial_points < 1 (the condition that selects the stbr_scipy branch). + data = [ + {"xi": [1.0, 2.0, 3.0, "95", "4", "15"], "yi": [1.0]}, + {"xi": [2.0, 3.0, 1.0, "105", "5", "25"], "yi": [0.5]}, + {"xi": [3.0, 1.0, 2.0, "115", "6", "35"], "yi": [0.8]}, + {"xi": [4.0, 2.0, 1.0, "125", "7", "45"], "yi": [0.2]}, + {"xi": [0.5, 4.0, 2.0, "95", "5", "35"], "yi": [0.6]}, + ] + request = { + "extras": { + "experimentSuggestionCount": 2, + "graphs": ["single"], + "graphFormat": "json", + "includeModel": "false", + }, + "data": data, + "optimizerConfig": { + "baseEstimator": "GP", + "acqFunc": "EI", + "initialPoints": 4, + "kappa": 1.96, + "xi": 0.01, + "space": space, + "constraints": [], + }, + } + + def _watchdog(signum, frame): + raise AssertionError( + "multi-suggestion run did not terminate within 60s — the slow " + "stbr_fill/stbr_scipy path has regressed" + ) + + signal.signal(signal.SIGALRM, _watchdog) + signal.alarm(60) + try: + result = optimizer.run_optimizer(body=request) + finally: + signal.alarm(0) + + assert len(result["result"]["next"]) == 2 + assert all(len(x) == len(space) for x in result["result"]["next"]) def test_selectedPoint_single_objective_json():